Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

add reactive_power_sigma attribute #409

Merged
merged 44 commits into from
Nov 3, 2023
Merged
Show file tree
Hide file tree
Changes from 37 commits
Commits
Show all changes
44 commits
Select commit Hold shift + click to select a range
a2c6a1d
add reactive_power_sigma attribute
mgovers Oct 19, 2023
fe83b9c
update docs
mgovers Oct 19, 2023
0cbab46
update docs
mgovers Oct 19, 2023
040e626
fix tests
mgovers Oct 19, 2023
ece11aa
update docs
mgovers Oct 19, 2023
2546126
minor
mgovers Oct 19, 2023
30c5d39
minor
mgovers Oct 19, 2023
89a19a4
update documentation
mgovers Oct 19, 2023
c7976c2
change reactive_power_sigma to p_sigma and q_sigma
mgovers Oct 19, 2023
6a6ca95
simplify docs
mgovers Oct 20, 2023
323103e
separate voltage and power sensor calc params
mgovers Oct 20, 2023
0ac2296
uniform complex random variable instead of complete calc param
mgovers Oct 25, 2023
04a68de
complex p, real p and q variances
mgovers Oct 27, 2023
74d03f0
partial fix build
mgovers Oct 27, 2023
e9bcebc
update all logic with asymmetric + split p and q variances
mgovers Oct 30, 2023
d740d24
fix tests pt1
mgovers Oct 30, 2023
d8460b0
Merge branch 'main' into feature/reactive-power-sigma
mgovers Oct 30, 2023
089bf02
Merge branch 'feature/reactive-power-sigma' of https://github.com/Pow…
mgovers Oct 30, 2023
6ebd35c
fix measurements tests
mgovers Oct 30, 2023
66d1f80
fix all tests
mgovers Oct 30, 2023
67c1062
test
mgovers Oct 31, 2023
c221f15
Merge branch 'main' into feature/reactive-power-sigma
mgovers Oct 31, 2023
48b09fc
Merge branch 'feature/reactive-power-sigma' of https://github.com/Pow…
mgovers Oct 31, 2023
38fe56a
fix all existing tests, add measured values tests
mgovers Nov 1, 2023
83a0ebb
Merge remote-tracking branch 'origin/main' into feature/reactive-powe…
mgovers Nov 1, 2023
09db890
Merge remote-tracking branch 'origin/feature/refactor-test-power-sens…
mgovers Nov 1, 2023
ce12a94
const correctness
mgovers Nov 1, 2023
3855384
make clang tidy happy
mgovers Nov 1, 2023
822f52f
also for voltage sensor
mgovers Nov 1, 2023
c206cd1
fix protected
mgovers Nov 1, 2023
3bf259d
fix more clang-tidy
mgovers Nov 1, 2023
5a92661
fix clang build
mgovers Nov 1, 2023
454b976
fix lifetime issue in eigen division
mgovers Nov 1, 2023
bf602c1
Merge branch 'main' into feature/reactive-power-sigma
TonyXiang8787 Nov 1, 2023
faafd04
tests for variance for different p and q sigma
mgovers Nov 1, 2023
8a5563d
Merge branch 'feature/reactive-power-sigma' of https://github.com/Pow…
mgovers Nov 1, 2023
8024b3e
fix sonar cloud
mgovers Nov 1, 2023
8c6e14e
resolve comments
mgovers Nov 2, 2023
3eb324e
efficient diagonal tensor multiplication
mgovers Nov 3, 2023
28a7352
pseudo abstract implementation for sensor
mgovers Nov 3, 2023
6ff436e
also ensure that it's not constructible from derived classes
mgovers Nov 3, 2023
bbc430b
correct optimization steps
mgovers Nov 3, 2023
fe5a30e
fix clang-tidy
mgovers Nov 3, 2023
8703d81
more clang-tidy
mgovers Nov 3, 2023
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
10 changes: 9 additions & 1 deletion code_generation/data/attribute_classes/input.json
Original file line number Diff line number Diff line change
Expand Up @@ -360,7 +360,7 @@
{
"data_type": "double",
"names": "power_sigma",
"description": "sigma of error margin of power measurement"
"description": "sigma of error margin of apparent power measurement"
}
]
},
Expand All @@ -376,6 +376,14 @@
"q_measured"
],
"description": "measured active/reactive power"
},
{
"data_type": "RealValue<sym>",
"names": [
"p_sigma",
"q_sigma"
],
"description": "sigma of error margin of active/reactive power measurement"
}
]
},
Expand Down
8 changes: 8 additions & 0 deletions code_generation/data/attribute_classes/update.json
Original file line number Diff line number Diff line change
Expand Up @@ -162,6 +162,14 @@
"q_measured"
],
"description": "measured active/reactive power"
},
{
"data_type": "RealValue<sym>",
"names": [
"p_sigma",
"q_sigma"
],
"description": "sigma of error margin of active/reactive power measurement"
}
]
},
Expand Down
5 changes: 4 additions & 1 deletion docs/_static/css/custom.css
Original file line number Diff line number Diff line change
Expand Up @@ -2,4 +2,7 @@

.wy-nav-content {
max-width: none;
}
}
.wy-table-responsive table td, .wy-table-responsive table th {
white-space: inherit;
}
TonyXiang8787 marked this conversation as resolved.
Show resolved Hide resolved
8 changes: 6 additions & 2 deletions docs/advanced_documentation/c-api.md
Original file line number Diff line number Diff line change
Expand Up @@ -89,8 +89,10 @@ In this way, the library is handling the memory (de-)allocation.
e.g., `aligned_alloc` and `free`.
You need to first call `PGM_meta_*` functions to retrieve the size and alignment of a component.

NOTE: Do not mix these two methods in creation and destruction.
```{warning}
Do not mix these two methods in creation and destruction.
You cannot use `PGM_destroy_buffer` to release a buffer created in your own code, or vice versa.
```

### Set and Get Attribute

Expand Down Expand Up @@ -127,8 +129,10 @@ The newly added attribute will have rubbish value in the memory.

Therefore, it is your choice of trade-off: maximum performance or backwards compatibility.

NOTE: you do not need to call `PGM_buffer_set_nan` on output buffers,
```{note}
You do not need to call `PGM_buffer_set_nan` on output buffers,
because the buffer will be overwritten in the calculation core with the real output data.
```

## Dataset views

Expand Down
53 changes: 34 additions & 19 deletions docs/user_manual/calculations.md
Original file line number Diff line number Diff line change
Expand Up @@ -21,10 +21,12 @@ Power flow is a "what-if" based grid calculation that will calculate the node vo
Some typical use-cases are network planning and contingency analysis.

Input:

- Network data: topology + component attributes
- Assumed load/generation profile

Output:

- Node voltage magnitude and angle
- Power flow through branches

Expand All @@ -35,10 +37,12 @@ network data and measurements. Measurements meaning power flow or voltage values
either measured, estimated or forecasted.

Input:

- Network data: topology + component attributes
- Power flow / voltage measurements with uncertainty

Output:

- Node voltage magnitude and angle
- Power flow through branches
- Deviation between measurement values and estimated state
Expand All @@ -63,30 +67,34 @@ $$
$$

The number of measurements can be found by the sum of the following:

- number of nodes with a voltage sensor with magnitude only
- two times the number of nodes with a voltage sensor with magnitude and angle
- two times the number of nodes without appliances connected
- two times the number of nodes where all connected appliances are measured by a power sensor
- two times the number of branches with a power sensor

Note: enough measurements doesn't necessarily mean that the system is observable. The location of the measurements is also
```{note}
Having enough measurements does not necessarily mean that the system is observable. The location of the measurements is also
of importance. Also, there should be at least one voltage measurement. The [iterative linear](#iterative-linear)
state estimation algorithm assumes voltage angles to be zero when not given. This might result in the calculation succeeding, but giving
a faulty outcome instead of raising a singular matrix error.
a faulty outcome instead of raising a singular matrix error.
```

#### Short Circuit Calculations


Short circuit calculation is carried out to analyze the worst case scenario when a fault has occurred.
The currents flowing through branches and node voltages are calculated.
Some typical use-cases are selection or design of components like conductors or breakers and power system protection, e.g. relay co-ordination.

Input:

- Network data: topology + component attributes
- Fault type and impedance.
- In the API call: choose between `minimum` and `maximum` voltage scaling to calculate the minimum or maximum short circuit currents (according to IEC 60909).

Output:

- Node voltage magnitude and angle
- Current flowing through branches and fault.

Expand All @@ -104,11 +112,11 @@ The table below can be used to pick the right algorithm. Below the table a more
| [Linear](calculations.md#linear) | &#10004; | | `CalculationMethod.linear` |
| [Linear current](calculations.md#linear-current) | &#10004; | | `CalculationMethod.linear_current` |

```note
```{note}
By default, the Newton-Raphson method is used.
```

```note
```{note}
When all the load/generation types are of constant impedance, power-grid-model uses the [Linear](#linear) method regardless of the input provided by the user.
This is because this method will then be accurate and fastest.
```
Expand Down Expand Up @@ -228,14 +236,15 @@ This algorithm is a Jacobi-like method for powerflow analysis.
It has linear convergence as opposed to quadratic convergence in the Newton-Raphson method. This means that the number of iterations will be greater. Newton-Raphson will also be more robust in achieving convergence in case of greater meshed configurations. However, the iterative current algorithm will be faster most of the time.

The algorithm is as follows:

1. Build $Y_{bus}$ matrix
2. Initialization of $U_N^0$ to $1$ plus the intrinsic phase shift of transformers
3. Calculate injected currents: $I_N^i$ for $i^{th}$ iteration. The injected currents are calculated as per ZIP model of loads and generation using $U_N$.
$
\begin{eqnarray}
I_N = \overline{S_{Z}} \cdot U_{N} + \overline{(\frac{S_{I}}{U_{N}})} \cdot |U_{N}| + \overline{(\frac{S_{P}}{U_N})}
\end{eqnarray}
$
$
\begin{eqnarray}
I_N = \overline{S_{Z}} \cdot U_{N} + \overline{(\frac{S_{I}}{U_{N}})} \cdot |U_{N}| + \overline{(\frac{S_{P}}{U_N})}
\end{eqnarray}
$
4. Solve linear equation: $YU_N^i = I_N^i$
5. Check convergence: If maximum voltage deviation from the previous iteration is greater than the tolerance setting (ie. $u^{(i-1)}_\sigma > u_\epsilon$), then go back to step 3.

Expand All @@ -254,8 +263,9 @@ It will be more accurate when most of the load/generation types are of constant
When all the load/generation types are of constant impedance, power-grid-model uses Linear method regardless of the input provided by the user. This is because this method will then be accurate and fastest.

The algorithm is as follows:

1. Assume injected currents by loads $I_N=0$ since we model loads/generation as impedance.
Admittance of each load/generation is calculated from rated power as $y=-\overline{S}$
Admittance of each load/generation is calculated from rated power as $y=-\overline{S}$
2. Build $Y{bus}$. Add the admittances of loads/generation to the diagonal elements.
3. Solve linear equation: $YU_N = I_N$ for $U_N$

Expand Down Expand Up @@ -355,20 +365,21 @@ $$
\begin{eqnarray}
\underline{S} = \sum_{k=1}^{N_{appliance}} \underline{S}_k
\quad\text{and}\quad
\sigma^2 = \sum_{k=1}^{N_{appliance}} \sigma_k^2
\sigma_P^2 = \sum_{k=1}^{N_{appliance}} \sigma_{P,k}^2
\sigma_Q^2 = \sum_{k=1}^{N_{appliance}} \sigma_{Q,k}^2
\end{eqnarray}
$$

Where $S_k$ and $\sigma_k$ are the measured value and the standard deviation of the individual appliances.
Where $S_k$ and $\sigma_{P,k}$ and $\sigma_{Q,k}$ are the measured value and the standard deviation of the individual appliances.

#### Iterative linear

Linear WLS requires all measurements to be linear. This is only possible is all measurements are phasor unit measurements,
which is not realistic in a distribution grid. Therefore, traditional measurements are linearized before the algorithm is performed:

- Bus voltage: Linear WLS requires a voltage phasor including a phase angle. Given that the phase shift in the distribution grid is very small,
it is assumed that the angle shift is zero plus the intrinsic phase shift of transformers. For a certain bus `i`, the voltage
magnitude measured at that bus is translated into a voltage phasor, where $\theta_i$ is the intrinsic transformer phase shift:
it is assumed that the angle shift is zero plus the intrinsic phase shift of transformers. For a certain bus `i`, the voltage
magnitude measured at that bus is translated into a voltage phasor, where $\theta_i$ is the intrinsic transformer phase shift:

$$
\begin{eqnarray}
Expand All @@ -395,6 +406,8 @@ $$
\end{eqnarray}
$$

The aggregated apparent power flow is considered as a single measurement, with variance $\sigma_S^2 = \sigma_P^2 + \sigma_Q^2$.

The assumption made in the linearization of measurements introduces a system error to the algorithm, because the phase shifts of
bus voltages are ignored in the input measurement data. This error is corrected by applying an iterative approach to the linear WLS
algorithm. In each iteration, the resulted voltage phase angle will be applied as the phase shift of the measured voltage phasor
Expand All @@ -421,8 +434,10 @@ the system error of the phase shift converges to zero. Because the matrix is pre
pre-factorized, the computation cost of each iteration is much smaller than Newton-Raphson
method, where the Jacobian matrix needs to be constructed and factorized each time.

NOTE: Since the algorithm will assume angles to be zero if not given, this might result in not having a
```{warning}
Since the algorithm will assume angles to be zero if not given, this might result in not having a
crash due to an unobservable system, but succeeding with the calculations and giving faulty results.
```

Algorithm call: `CalculationMethod.iterative_linear`

Expand All @@ -441,11 +456,11 @@ The assumptions used for calculations in power-grid-model are aligned to the one
- The pre-fault voltage is considered in the calculation and is calculated based on the grid parameters and topology. (Excl. loads and generation)
- The calculations are assumed to be time-independent. (Voltages are sine throughout with the fault occurring at a zero crossing of the voltage, the complexity of rotating machines and harmonics are neglected, etc.)
- To account for the different operational conditions, a voltage scaling factor of `c` is applied to the voltage source while running short circuit calculation function.
The factor `c` is determined by the nominal voltage of the node that the source is connected to and the API option to calculate the `minimum` or `maximum` short circuit currents.
The table to derive `c` according to IEC 60909 is shown below.
The factor `c` is determined by the nominal voltage of the node that the source is connected to and the API option to calculate the `minimum` or `maximum` short circuit currents.
The table to derive `c` according to IEC 60909 is shown below.

| Algorithm | c_max | c_min |
|----------------|-------|-------|
| -------------- | ----- | ----- |
| `U_nom` <= 1kV | 1.10 | 0.95 |
| `U_nom` > 1kV | 1.10 | 1.00 |

Expand Down
31 changes: 21 additions & 10 deletions docs/user_manual/components.md
Original file line number Diff line number Diff line change
Expand Up @@ -107,7 +107,7 @@ usually permanently connects two joints. In this case, the attribute `from_statu

### Line

* type name: 'line'
* type name: `line`

`line` is a {hoverxreftooltip}`user_manual/components:branch` with specified serial impedance and shunt admittance. A cable is
also modeled as `line`. A `line` can only connect two nodes with the same rated voltage.
Expand Down Expand Up @@ -463,15 +463,14 @@ measuring a `source`, a positive `p_measured` indicates that the active power fl
2. The node injection power sensor gets placed on a node.
In the state estimation result, the power from this injection is distributed equally among the connected appliances at that node.
Because of this distribution, at least one appliance is required to be connected to the node where an injection sensor is placed for it to function.

```

##### Input

| name | data type | unit | description | required | update | valid values |
| ------------------------ | ----------------------------------------------------------------------------- | ---------------- | --------------------------------------------------------------------------------------------------------------- | :--------------------------------: | :------: | :--------------------------------------------------: |
| `measured_terminal_type` | {py:class}`MeasuredTerminalType <power_grid_model.enum.MeasuredTerminalType>` | - | indicate if it measures an `appliance` or a `branch` | &#10004; | &#10060; | the terminal type should match the `measured_object` |
| `power_sigma` | `double` | volt-ampere (VA) | standard deviation of the measurement error. Usually this is the absolute measurement error range divided by 3. | &#10024; only for state estimation | &#10004; | `> 0` |
| name | data type | unit | description | required | update | valid values |
| ------------------------ | ----------------------------------------------------------------------------- | ---------------- | ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------- | :---------------------------------: | :------: | :--------------------------------------------------: |
| `measured_terminal_type` | {py:class}`MeasuredTerminalType <power_grid_model.enum.MeasuredTerminalType>` | - | indicate if it measures an `appliance` or a `branch` | &#10004; | &#10060; | the terminal type should match the `measured_object` |
| `power_sigma` | `double` | volt-ampere (VA) | standard deviation of the measurement error. Usually this is the absolute measurement error range divided by 3. See {hoverxreftooltip}`user_manual/components:Power Sensor Complete Types`. | &#10024; only for state estimation. | &#10004; | `> 0` |

#### Power Sensor Concrete Types

Expand All @@ -485,10 +484,22 @@ the meaning of `RealValueInput` is different, as shown in the table below.

##### Input

| name | data type | unit | description | required | update |
| ------------ | ---------------- | -------------------------- | ----------------------- | :--------------------------------: | :------: |
| `p_measured` | `RealValueInput` | watt (W) | measured active power | &#10024; only for state estimation | &#10004; |
| `q_measured` | `RealValueInput` | volt-ampere-reactive (var) | measured reactive power | &#10024; only for state estimation | &#10004; |
| name | data type | unit | description | required | update |
| ------------ | ---------------- | -------------------------- | ------------------------------------------------------------------------------------------------------------------------------ | :---------------------------------: | :------: |
| `p_measured` | `RealValueInput` | watt (W) | measured active power | &#10024; only for state estimation | &#10004; |
| `q_measured` | `RealValueInput` | volt-ampere-reactive (var) | measured reactive power | &#10024; only for state estimation | &#10004; |
| `p_sigma` | `RealValueInput` | watt (W) | standard deviation of the active power measurement error. Usually this is the absolute measurement error range divided by 3. | &#10060; see the explanation below. | &#10004; | `> 0` |
| `q_sigma` | `RealValueInput` | volt-ampere-reactive (var) | standard deviation of the reactive power measurement error. Usually this is the absolute measurement error range divided by 3. | &#10060; see the explanation below. | &#10004; | `> 0` |

```{note}
1. If both `p_sigma` and `q_sigma` are provided, they represent the standard deviation of the active and reactive power, respectively, and the value of `power_sigma` is ignored. Any infinite component invalidates the entire measurement.

2. If neither `p_sigma` nor `q_sigma` are provided, `power_sigma` represents the standard deviation of the apparent power.

3. Providing only one of `p_sigma` and `q_sigma` results in undefined behaviour.
```

See the documentation on [state estimation calculation methods](calculations.md#state-estimation-algorithms) for details per method on how the variances are taken into account for both the active and reactive power and for the individual phases.

##### Steady state output

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -179,13 +179,15 @@ using AsymVoltageSensorInput = VoltageSensorInput<false>;

struct GenericPowerSensorInput : SensorInput {
MeasuredTerminalType measured_terminal_type; // type of measured terminal
double power_sigma; // sigma of error margin of power measurement
double power_sigma; // sigma of error margin of apparent power measurement
};

template <bool sym>
struct PowerSensorInput : GenericPowerSensorInput {
RealValue<sym> p_measured; // measured active/reactive power
RealValue<sym> q_measured; // measured active/reactive power
RealValue<sym> p_sigma; // sigma of error margin of active/reactive power measurement
RealValue<sym> q_sigma; // sigma of error margin of active/reactive power measurement
};
using SymPowerSensorInput = PowerSensorInput<true>;
using AsymPowerSensorInput = PowerSensorInput<false>;
Expand Down