(m3b9)=

# Time-domain Simulation of Electric Power Systems

*Author: [Ruizhi Yu](https://github.com/rzyu45)*

Synchronous machine is the core of electric power systems. The classic second-order models of the machines are

```{math}
\left\{
    \begin{aligned}
        &\dot{\omega} = \frac{P_\text{m}-P_\text{e}-D(\omega-1)}{T_j}\\
        &\dot{\delta} = (\omega-\omega_\text{COI})\omega_\text{B}\\
        &E_\text{d}'=\sin\delta\qty(U_\text{x}+r_\text{a}I_\text{x}-X_\text{q}'I_\text{y})-\cos\delta\qty(U_\text{y}+r_\text{a}I_\text{y}+X_\text{q}'I_\text{x})\\
        &E_\text{q}'=\cos\delta\qty(U_\text{x}+r_\text{a}I_\text{x}-X_\text{d}'I_\text{y})+\sin\delta\qty(U_\text{y}+r_\text{a}I_\text{y}+X_\text{d}'I_\text{x})\\
        &I_{\text{x},i}-\sum \qty(G_{ij}U_{\text{x},j}-B_{ij}U_{\text{y},j})=0\\
        &I_{\text{y},i}-\sum \qty(G_{ij}U_{\text{y},j}+B_{ij}U_{\text{x},j})=0\\
        &P_\text{e}=U_\text{x}I_\text{x}+U_\text{y}I_\text{y}+\qty(I_\text{x}^2+I_\text{y}^2)r_\text{a}
    \end{aligned}
\right.
```
where $\omega$ is the rotor speed; $\delta$ is the rotor angle; $P_\text{m}$ and $P_\text{e}$ are mechanical and electric power respectively; $T_j$ is the inertial constant; $D$ is the damping coefficient; $\omega_\text{COI}$ is the center-of-inertial rotor speed; $\omega_\text{B}$ is the base synchronous frequency; $E_\text{d}'$ and $E_\text{q}'$ are the d-axis and q-axis transient internal voltages; $U_\text{x}$ and $U_\text{y}$ are respectively the real and imaginary parts of voltages; $I_\text{x}$ and $I_\text{y}$ are respectively the real and imaginary parts of currents; $r_\text{a}$, $X_\text{d}'$ and $X_\text{q}'$ are generator parameters; $G$ and $B$ are respectively the conductance and susceptance matrices.


Below is the diagram of an electric power system with three machines and nine buses.

![omega](fig/m3b9.png)

We want to know, through simulation, if the three machines can keep synchronous in the scenario where the three-phase short-cut fault occurs on bus 6 and is cleared 30ms later. This is a typical security analysis usage of time-domain simulation.

The three-phase short-cut fault can be mathematically interpreted as the surge of node conductance of bus 6.

In Solverz, this can be realized by setting $G_{66}$ as a `TimeSeriesParam`:

```python
m.G66 = TimeSeriesParam('G66',
                        v_series=[G[6, 6], 10000, 10000, G[6, 6], G[6, 6]],
                        time_series=[0, 0.002, 0.03, 0.032, 10])
```
During the fault, the value of $G_{66}$ is set to be 10000 and is restored after having cleared the fault. The time series specify the fault-clearing time.


The Solverz implementation of the above time-domain simulation is as follows wherein the `test_m3b9.xlsx` file can be found in the [directory of the source file](https://github.com/rzyu45/Solverz-Cookbook/tree/main/docs/source/dae/m3b9).

```{literalinclude} src/plot_m3b9.py
```

We have

```{eval-rst}
.. plot:: dae/m3b9/src/plot_m3b9.py
```

We find the system remains stable after the three-phase fault disturbance.