Skip to content

Commit ed825ab

Browse files
committed
Fix bad parameter initialization and polish docs
1 parent 78aec28 commit ed825ab

5 files changed

Lines changed: 497 additions & 390 deletions

File tree

‎GridKit/Model/PhasorDynamics/Stabilizer/IEEEST/Ieeest.hpp‎

Lines changed: 1 addition & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -138,34 +138,18 @@ namespace GridKit
138138
RealT Vcu_{0};
139139
RealT Tdelay_{0};
140140

141-
RealT a0_{1};
142141
RealT a1_{0};
143142
RealT a2_{0};
144143
RealT a3_{0};
145144
RealT a4_{0};
146145

147-
// Precomputed masks and safe inverse coefficients for branch-free degenerate paths.
148-
RealT use_notch_{0};
149-
RealT bypass_notch_{1};
150-
RealT use_4th_order_{0};
151-
RealT use_3rd_order_{0};
152-
RealT use_2nd_order_{0};
153-
RealT safe_inv_a4_{0};
154-
RealT safe_inv_a3_{0};
155-
RealT safe_inv_a2_{0};
156-
RealT use_T2_block_{1};
157-
RealT bypass_T2_block_{0};
158-
RealT use_T4_block_{1};
159-
RealT bypass_T4_block_{0};
160-
RealT use_T6_block_{1};
161-
RealT bypass_T6_block_{0};
162-
163146
ComponentSignals<ScalarT, IdxT, IeeestInternalVariables, IeeestExternalVariables> signals_;
164147

165148
std::unique_ptr<MonitorT> monitor_;
166149

167150
void initializeParameters(const ModelDataT& data);
168151
void initializeMonitor();
152+
void setDerivedParameters();
169153

170154
std::vector<ScalarT> ws_;
171155
std::vector<IdxT> ws_indices_;

‎GridKit/Model/PhasorDynamics/Stabilizer/IEEEST/IeeestImpl.hpp‎

Lines changed: 62 additions & 43 deletions
Original file line numberDiff line numberDiff line change
@@ -43,6 +43,15 @@ namespace GridKit
4343
{
4444
}
4545

46+
template <typename scalar_type, typename index_type>
47+
void Ieeest<scalar_type, index_type>::setDerivedParameters()
48+
{
49+
a1_ = A1_ + A3_;
50+
a2_ = A2_ + A4_ + A1_ * A3_;
51+
a3_ = A1_ * A4_ + A2_ * A3_;
52+
a4_ = A2_ * A4_;
53+
}
54+
4655
template <typename scalar_type, typename index_type>
4756
void Ieeest<scalar_type, index_type>::initializeParameters(const ModelDataT& data)
4857
{
@@ -119,31 +128,7 @@ namespace GridKit
119128
Tdelay_ = std::get<RealT>(data.parameters.at(ModelDataT::Parameters::Tdelay));
120129
}
121130

122-
a0_ = 1;
123-
a1_ = A1_ + A3_;
124-
a2_ = A2_ + A4_ + A1_ * A3_;
125-
a3_ = A1_ * A4_ + A2_ * A3_;
126-
a4_ = A2_ * A4_;
127-
128-
// Precompute masks and safe inverse coefficients so the residual stays branch-free.
129-
use_notch_ = static_cast<RealT>(a2_ != 0.0 || a3_ != 0.0 || a4_ != 0.0);
130-
bypass_notch_ = 1.0 - use_notch_;
131-
132-
use_4th_order_ = static_cast<RealT>(a4_ != 0.0);
133-
use_3rd_order_ = static_cast<RealT>(a4_ == 0.0 && a3_ != 0.0);
134-
use_2nd_order_ = static_cast<RealT>(a4_ == 0.0 && a3_ == 0.0 && a2_ != 0.0);
135-
safe_inv_a4_ = use_4th_order_ / (a4_ + (1.0 - use_4th_order_));
136-
safe_inv_a3_ = use_3rd_order_ / (a3_ + (1.0 - use_3rd_order_));
137-
safe_inv_a2_ = use_2nd_order_ / (a2_ + (1.0 - use_2nd_order_));
138-
139-
use_T2_block_ = static_cast<RealT>(T2_ != 0.0);
140-
bypass_T2_block_ = 1.0 - use_T2_block_;
141-
142-
use_T4_block_ = static_cast<RealT>(T4_ != 0.0);
143-
bypass_T4_block_ = 1.0 - use_T4_block_;
144-
145-
use_T6_block_ = static_cast<RealT>(T6_ != 0.0);
146-
bypass_T6_block_ = 1.0 - use_T6_block_;
131+
setDerivedParameters();
147132
}
148133

149134
template <typename scalar_type, typename index_type>
@@ -182,6 +167,8 @@ namespace GridKit
182167
&y_[11], &(this->getVariableIndex(11)));
183168
}
184169

170+
tagDifferentiable();
171+
185172
return 0;
186173
}
187174

@@ -204,7 +191,7 @@ namespace GridKit
204191
ret += 1;
205192
}
206193

207-
if (a4_ == 0 && a3_ == 0 && a2_ == 0 && a1_ != 0)
194+
if (a2_ == ZERO<RealT> && a3_ == ZERO<RealT> && a4_ == ZERO<RealT> && a1_ != ZERO<RealT>)
208195
{
209196
Log::error() << "Ieeest: a2, a3, and a4 are all zero - no valid notch filter\n";
210197
ret += 1;
@@ -222,19 +209,53 @@ namespace GridKit
222209
yp_[static_cast<size_t>(i)] = 0.0;
223210
}
224211

212+
ScalarT u{0.0};
213+
if (signals_.template isAttached<IeeestExternalVariables::U>())
214+
{
215+
u = signals_.template readExternalVariable<IeeestExternalVariables::U>();
216+
ws_[0] = u;
217+
ws_indices_[0] = signals_.template readExternalVariableIndex<IeeestExternalVariables::U>();
218+
}
219+
220+
const ScalarT zero{0.0};
221+
const ScalarT x1 = u;
222+
const ScalarT x2 = zero;
223+
const ScalarT x3 = zero;
224+
const ScalarT x4 = zero;
225+
const ScalarT v4 = x1 + A5_ * x2 + A6_ * x3;
226+
const ScalarT x5 = v4;
227+
const ScalarT v5 = v4;
228+
const ScalarT x6 = v5;
229+
const ScalarT v6 = v5;
230+
const ScalarT x7 = v6;
231+
const ScalarT v7 = (T6_ == ZERO<RealT>) ? Ks_ * v6 : zero;
232+
233+
y_[0] = x1;
234+
y_[1] = x2;
235+
y_[2] = x3;
236+
y_[3] = x4;
237+
y_[4] = x5;
238+
y_[5] = x6;
239+
y_[6] = x7;
240+
y_[7] = v4;
241+
y_[8] = v5;
242+
y_[9] = v6;
243+
y_[10] = v7;
244+
y_[11] = Math::clamp(v7, Lsmin_, Lsmax_);
245+
225246
return 0;
226247
}
227248

228249
template <typename scalar_type, typename index_type>
229250
int Ieeest<scalar_type, index_type>::tagDifferentiable()
230251
{
231-
tag_[0] = true;
232-
tag_[1] = true;
233-
tag_[2] = true;
234-
tag_[3] = true;
235-
tag_[4] = (T2_ != 0.0);
236-
tag_[5] = (T4_ != 0.0);
237-
tag_[6] = (T6_ != 0.0);
252+
tag_[0] = (a2_ != ZERO<RealT> || a3_ != ZERO<RealT> || a4_ != ZERO<RealT>);
253+
tag_[1] = tag_[0];
254+
tag_[2] = (a3_ != ZERO<RealT> || a4_ != ZERO<RealT>);
255+
tag_[3] = (a4_ != ZERO<RealT>);
256+
tag_[4] = (T2_ != ZERO<RealT>);
257+
tag_[5] = (T4_ != ZERO<RealT>);
258+
tag_[6] = (T6_ != ZERO<RealT>);
238259
tag_[7] = false;
239260
tag_[8] = false;
240261
tag_[9] = false;
@@ -294,19 +315,17 @@ namespace GridKit
294315

295316
ScalarT u = ws[0];
296317

297-
f[0] = -x1_dot + use_notch_ * x2;
298-
f[1] = -x2_dot + (use_4th_order_ + use_3rd_order_) * x3
299-
+ use_2nd_order_ * (-a0_ * x1 - a1_ * x2 + u) * safe_inv_a2_;
300-
f[2] = -x3_dot + use_4th_order_ * x4
301-
+ use_3rd_order_ * (-a0_ * x1 - a1_ * x2 - a2_ * x3 + u) * safe_inv_a3_;
302-
f[3] = -x4_dot + use_4th_order_ * (-a0_ * x1 - a1_ * x2 - a2_ * x3 - a3_ * x4 + u) * safe_inv_a4_;
318+
f[0] = -tag_[0] * x1_dot + x2;
319+
f[1] = -tag_[1] * x2_dot + x3;
320+
f[2] = -tag_[2] * x3_dot + x4;
321+
f[3] = -a4_ * x4_dot - x1 - a1_ * x2 - a2_ * x3 - a3_ * x4 + u;
303322
f[4] = -T2_ * x5_dot - x5 + v4;
304323
f[5] = -T4_ * x6_dot - x6 + v5;
305324
f[6] = -T6_ * x7_dot - x7 + v6;
306-
f[7] = -v4 + bypass_notch_ * u + use_notch_ * (x1 + A5_ * x2 + (use_4th_order_ + use_3rd_order_) * A6_ * x3);
307-
f[8] = use_T2_block_ * (-T2_ * (v5 - x5) + T1_ * (v4 - x5)) + bypass_T2_block_ * (v4 - v5);
308-
f[9] = use_T4_block_ * (-T4_ * (v6 - x6) + T3_ * (v5 - x6)) + bypass_T4_block_ * (v5 - v6);
309-
f[10] = use_T6_block_ * (-T6_ * v7 + Ks_ * T5_ * (v6 - x7)) + bypass_T6_block_ * (Ks_ * v6 - v7);
325+
f[7] = -v4 + x1 + A5_ * x2 + A6_ * x3;
326+
f[8] = tag_[4] * (-T2_ * (v5 - x5) + T1_ * (v4 - x5)) + (1 - tag_[4]) * (v4 - v5);
327+
f[9] = tag_[5] * (-T4_ * (v6 - x6) + T3_ * (v5 - x6)) + (1 - tag_[5]) * (v5 - v6);
328+
f[10] = tag_[6] * (-T6_ * v7 + Ks_ * T5_ * (v6 - x7)) + (1 - tag_[6]) * (Ks_ * v6 - v7);
310329
f[11] = -vss + Math::clamp(v7, Lsmin_, Lsmax_);
311330

312331
return 0;

‎GridKit/Model/PhasorDynamics/Stabilizer/IEEEST/README.md‎

Lines changed: 105 additions & 35 deletions
Original file line numberDiff line numberDiff line change
@@ -3,6 +3,12 @@
33
Standard IEEE power system stabilizer: 4th-order notch filter, two lead–lag
44
blocks, washout, and output limiter.
55

6+
Notes:
7+
- $V_{cl}$, $V_{cu}$, and $T_{delay}$ are accepted for input-format
8+
compatibility but are not modeled.
9+
- A zero denominator time constant bypasses its corresponding lead–lag or
10+
washout block.
11+
612
## Block Diagram
713

814
![](../../../../../docs/Figures/stabilizer_ieeest_diagram.png)
@@ -11,40 +17,68 @@ Figure 1: Stabilizer IEEEST model. Figure courtesy of [PowerWorld](https://www.p
1117

1218
## Model Parameters
1319

14-
Symbol | Units | Description | Typical Value
15-
------------|--------|--------------------------------------|--------------
16-
$A_1$ | [s] | Notch denominator coefficient | 1.013
17-
$A_2$ | [s²] | Notch denominator coefficient | 0.013
18-
$A_3$ | [s] | Notch denominator coefficient | 0.0
19-
$A_4$ | [s²] | Notch denominator coefficient | 0.0
20-
$A_5$ | [s] | Notch numerator coefficient | 1.013
21-
$A_6$ | [s²] | Notch numerator coefficient | 0.113
22-
$T_1$ | [s] | Lead–lag 1 numerator time constant | 0.0
23-
$T_2$ | [s] | Lead–lag 1 denominator time constant | 0.02
24-
$T_3$ | [s] | Lead–lag 2 numerator time constant | 0.0
25-
$T_4$ | [s] | Lead–lag 2 denominator time constant | 0.0
26-
$T_5$ | [s] | Washout numerator time constant | 1.65
27-
$T_6$ | [s] | Washout denominator time constant | 1.65
28-
$K_s$ | [p.u.] | Stabilizer gain | 3.0
29-
$L_s^{\min}$ | [p.u.] | Minimum stabilizer output limit | -0.1
30-
$L_s^{\max}$ | [p.u.] | Maximum stabilizer output limit | 0.1
31-
32-
The IEEE 421.5 IEEEST also defines a cutout window ($V_{cl}$, $V_{cu}$) and an
33-
input delay ($T_{delay}$). These parameters are accepted for input-format
34-
compatibility but are not modeled here.
35-
36-
### Derived Parameters
20+
Symbol | Units | JSON | Description | Typical Value | Note
21+
-------------|--------|----------|--------------------------------------|---------------|------
22+
$A_1$ | [s] | `A1` | Notch denominator coefficient | 1.013 |
23+
$A_2$ | [s²] | `A2` | Notch denominator coefficient | 0.013 |
24+
$A_3$ | [s] | `A3` | Notch denominator coefficient | 0.0 |
25+
$A_4$ | [s²] | `A4` | Notch denominator coefficient | 0.0 |
26+
$A_5$ | [s] | `A5` | Notch numerator coefficient | 1.013 |
27+
$A_6$ | [s²] | `A6` | Notch numerator coefficient | 0.113 |
28+
$T_1$ | [s] | `T1` | Lead–lag 1 numerator time constant | 0.0 |
29+
$T_2$ | [s] | `T2` | Lead–lag 1 denominator time constant | 0.02 |
30+
$T_3$ | [s] | `T3` | Lead–lag 2 numerator time constant | 0.0 |
31+
$T_4$ | [s] | `T4` | Lead–lag 2 denominator time constant | 0.0 |
32+
$T_5$ | [s] | `T5` | Washout numerator time constant | 1.65 |
33+
$T_6$ | [s] | `T6` | Washout denominator time constant | 1.65 |
34+
$K_s$ | [p.u.] | `Ks` | Stabilizer gain | 3.0 |
35+
$L_s^{\min}$ | [p.u.] | `Lsmin` | Minimum stabilizer output limit | -0.1 |
36+
$L_s^{\max}$ | [p.u.] | `Lsmax` | Maximum stabilizer output limit | 0.1 |
37+
$V_{cl}$ | [p.u.] | `Vcl` | Lower input cutout threshold | 0.0 | Accepted but not modeled
38+
$V_{cu}$ | [p.u.] | `Vcu` | Upper input cutout threshold | 0.0 | Accepted but not modeled
39+
$T_\text{delay}$ | [s] | `Tdelay` | Input delay | 0.0 | Accepted but not modeled
40+
41+
### Parameter Validation
42+
43+
The fixed realization rejects a first-order-only notch denominator:
44+
45+
```math
46+
\begin{aligned}
47+
a_2 \ne 0 \lor a_3 \ne 0 \lor a_4 \ne 0 \lor a_1 = 0
48+
\end{aligned}
49+
```
50+
51+
### Model Derived Parameters
52+
53+
The notch-filter denominator expands to:
3754

3855
```math
3956
\begin{aligned}
40-
a_0 &= 1 \\
4157
a_1 &= A_1 + A_3 \\
4258
a_2 &= A_2 + A_4 + A_1 A_3 \\
4359
a_3 &= A_1 A_4 + A_2 A_3 \\
4460
a_4 &= A_2 A_4
4561
\end{aligned}
4662
```
4763

64+
The binary DAE selectors choose the active notch-filter order:
65+
66+
```math
67+
\begin{aligned}
68+
\delta_1 &= \delta_2 =
69+
\begin{cases}
70+
1 & a_2 \ne 0 \lor a_3 \ne 0 \lor a_4 \ne 0 \\
71+
0 & \text{otherwise}
72+
\end{cases}
73+
\\
74+
\delta_3 &=
75+
\begin{cases}
76+
1 & a_3 \ne 0 \lor a_4 \ne 0 \\
77+
0 & \text{otherwise}
78+
\end{cases}
79+
\end{aligned}
80+
```
81+
4882
## Model Variables
4983

5084
### Internal Variables
@@ -58,6 +92,9 @@ $x_5$ | [-] | Lead–lag 1 state
5892
$x_6$ | [-] | Lead–lag 2 state
5993
$x_7$ | [-] | Washout state
6094

95+
For reduced-order notch filters, unused notch states remain in the fixed
96+
component state vector and are pinned by algebraic residuals.
97+
6198
#### Algebraic
6299

63100
Symbol | Units | Description
@@ -82,10 +119,10 @@ $u$ | [p.u.] | Stabilizer input signal
82119

83120
```math
84121
\begin{aligned}
85-
0 &= -\dot{x}_1 + x_2 \\
86-
0 &= -\dot{x}_2 + x_3 \\
87-
0 &= -\dot{x}_3 + x_4 \\
88-
0 &= -\dot{x}_4 - \dfrac{a_0}{a_4}x_1 - \dfrac{a_1}{a_4}x_2 - \dfrac{a_2}{a_4}x_3 - \dfrac{a_3}{a_4}x_4 + \dfrac{1}{a_4}u \\
122+
0 &= -\delta_1\dot{x}_1 + x_2 \\
123+
0 &= -\delta_2\dot{x}_2 + x_3 \\
124+
0 &= -\delta_3\dot{x}_3 + x_4 \\
125+
0 &= -a_4\dot{x}_4 - x_1 - a_1x_2 - a_2x_3 - a_3x_4 + u \\
89126
0 &= -T_2 \dot{x}_5 - x_5 + v_4 \\
90127
0 &= -T_4 \dot{x}_6 - x_6 + v_5 \\
91128
0 &= -T_6 \dot{x}_7 - x_7 + v_6
@@ -96,10 +133,22 @@ $u$ | [p.u.] | Stabilizer input signal
96133

97134
```math
98135
\begin{aligned}
99-
0 &= -v_4 + x_1 + A_5 x_2 + A_6 x_3 \\
100-
0 &= -T_2(v_5 - x_5) + T_1(v_4 - x_5) \\
101-
0 &= -T_4(v_6 - x_6) + T_3(v_5 - x_6) \\
102-
0 &= -T_6 v_7 + K_s T_5(v_6 - x_7) \\
136+
0 &= -v_4 + x_1 + A_5x_2 + A_6x_3 \\
137+
0 &=
138+
\begin{cases}
139+
-T_2(v_5 - x_5) + T_1(v_4 - x_5) & T_2 \ne 0 \\
140+
v_4 - v_5 & T_2 = 0
141+
\end{cases} \\
142+
0 &=
143+
\begin{cases}
144+
-T_4(v_6 - x_6) + T_3(v_5 - x_6) & T_4 \ne 0 \\
145+
v_5 - v_6 & T_4 = 0
146+
\end{cases} \\
147+
0 &=
148+
\begin{cases}
149+
-T_6 v_7 + K_s T_5(v_6 - x_7) & T_6 \ne 0 \\
150+
K_s v_6 - v_7 & T_6 = 0
151+
\end{cases} \\
103152
0 &= -V_{ss} + \text{clamp}(v_7, L_s^{\min}, L_s^{\max})
104153
\end{aligned}
105154
```
@@ -109,6 +158,27 @@ The output limiter uses GridKit's smooth
109158

110159
## Initialization
111160

112-
All states and their derivatives initialize to zero. The stabilizer comes
113-
online at rest and produces signal only in response to deviations in the input
114-
$u$.
161+
States and derivatives initialize to the steady state implied by the attached
162+
input $u$:
163+
164+
```math
165+
\begin{aligned}
166+
x_1 &= v_4 = x_5 = v_5 = x_6 = v_6 = x_7 = u \\
167+
x_2 &= x_3 = x_4 = 0 \\
168+
v_7 &=
169+
\begin{cases}
170+
K_su & T_6 = 0 \\
171+
0 & \text{otherwise}
172+
\end{cases}
173+
&
174+
V_{ss} &= \text{clamp}(v_7, L_s^{\min}, L_s^{\max})
175+
\end{aligned}
176+
```
177+
178+
All internal derivatives initialize to zero.
179+
180+
## Model Outputs
181+
182+
Output | Units | Description | Note
183+
-------|--------|---------------------------|-----
184+
`vss` | [p.u.] | Limited stabilizer signal | Exported through `output` when assigned

0 commit comments

Comments
 (0)