Make PhasorDynamics models Hessenberg-ready

[lukelowry]

Summary

Track PhasorDynamics models whose state equations are not in explicit derivative form, especially parameter-switched cases where a state can become algebraic.

Rationale

In #439, IEEET1 supports $T_R = 0$ by retagging $V_{ts}$ as algebraic. That works for fully implicit DAE residuals, but can block future Hessenberg-form solvers that need separate differential right-hand-side and algebraic residual methods.

Description

Define a convention for these parametrically piecewise models, then update each affected model so it has a Hessenberg-compatible representation.

Models to cover:

• IEEET1
• TGOV1
• IEEEST
• EXAC1
• EXAC2
• ESAC6A
• ESDC2A
• ESST4B
• EXPIC1
• SCRX
• REGCA
• REECA
• IEEEG1
• GGOV1

Additional information

Any changes would be mostly mechanical tasks; unclear what the final desired algebraic/differential variant implementation looks like at this time.

We also need to determine:

• how to document conditionally algebraic states/equations
• how to reduce the count of duplicate equations if we have multiple variations of each model

[Steven-Roberts]

Here are some thoughts on what Hessenberg functionality would be needed for the time integrators:

• We'll need to access and update differential and algebraic subvectors separately. Perhaps this could be done by adding y_diff, y_p_diff, y_alg, and y_p_alg functions to Evaluator. We'll probably want to ensure differential variables for a model are grouped into a contiguous memory block, and likewise for the algebraic variables.
• The partitioned integrators use IDA for solving components. I think the existing residual and Jacobian functions should suffice for this.
• The partitioned integrators also perform a nonlinear solve using all algebraic constraints. This ought to necessitate a function like Evaluator::getAlgebraicJacobian which computes the Jacobian of the algebraic RHS with respect to algebraic variables. For models which are compositions of other models, e.g., SystemModel, the getAlgebraicJacobian function would produce a block diagonal matrix. At this point, I don't think other blocks of the Jacobian are needed.

[pelesh]

The challenge is that with this modeling approach what equations are differentiable and what algebraic is determined at runtime based on model parameter setting. The model in #460 apparently can produce 31 different combinations of differential and algebraic equations.

[pelesh]

Per @Steven-Roberts comment:

After giving some more thought to the approach of adding a dummy variable to handle a variable that is diff/alg depending on a parameter, I realized this makes the index of the DAE parameter dependent. In the worst case, the DAE would be index-2. IDA can handle index-2 DAEs, although I think the accuracy of the algebraic variables drops by one order. I suspect the partitioned DAE integrators would have more serious order reduction.

I suggest we set minimum values for inertia parameters as the short-term solution.

[lukelowry]

WECC240 case

Comparing (current) $T\geq10^{-3}$ vs $T\geq10^{-5}$ on all models with a lower bound on inertial/time constant $T$

Validation error against reference solution in power world

field	metric	$10^{-3}$	$10^{-5}$	change
Δω	abs max	5.715e-06	4.669e-06	-18.3 %
Δω	rel max	0.3423 %	0.2797 %	-18.3 %
Δω	abs rmse	9.826e-07	8.206e-07	-16.5 %
Δω	rel rmse	0.6462 %	0.5397 %	-16.5 %
V	abs max	1.409e-04	1.416e-04	+0.5 %
V	rel max	0.0124 %	0.0125 %	+0.5 %
V	abs rmse	1.117e-05	1.056e-05	-5.4 %
V	rel rmse	0.0011 %	0.0010 %	-5.4 %
P	abs max	1.693e-02	1.063e-02	-37.2 %
P	rel max	0.0217 %	0.0136 %	-37.2 %
P	abs rmse	8.279e-04	6.418e-04	-22.5 %
P	rel rmse	0.0049 %	0.0038 %	-22.5 %
Q	abs max	3.545e-02	3.574e-02	+0.8 %
Q	rel max	0.1648 %	0.1662 %	+0.8 %
Q	abs rmse	2.003e-03	1.925e-03	-3.9 %
Q	rel rmse	0.0994 %	0.0955 %	-3.9 %

This costs performance, however. The implicit form with algebraic tagging observes the same error as using $T=10^{-5}$ but is dramatically faster. Sometimes twice as fast