2.1 Mathematical model

The electric network is modelled as coupled synchronous machines combined with their respective frequency controllers (secondary) and by considering the voltage dynamics, as follows:

(2.1) \begin{equation}\left\{ \begin{array}{ll} \dot{\boldsymbol{\theta}}_i=\omega_i\\ \dot{\omega}_i=-\alpha_{i}\dot{\boldsymbol{\theta}}_i+P_{i}^{*}-\sum\limits_{j=1}^N E^{\prime}_{q,i}B_{ij}E^{\prime}_{q,j}\sin(\theta_{i}-\theta_{j})+ u_{i}\\ T^{\prime}_{{d,i}_{0}}\dot{\overline{E}}^{\prime}_{q,i}=E_{f,i}-E^{\prime}_{q,i}+ (X_{d,i}-X^{\prime}_{d,i})\sum\limits_{j=1}^N B_{ij}E^{\prime}_{q,j}\cos(\theta_{i}-\theta_{j})\\ \tau_{g_i}\dot{u}_i= -u_{i} -\gamma_{i}\theta_{i} - \beta_{i}\omega_i, \end{array}\right. \end{equation}

where $i\in {1,...,N}$ denotes the node index in the network, $\theta$ the voltage phase angle, $\omega$ the frequency, E the voltage amplitude and u the control. $\alpha$ is a damping constant, $P^*$ the power consumed/generated at a node, B the susceptibility matrix, while T and $\tau_{g}$ are time constants, $E_{f}$ is the rotor’s field voltage, and $X_d$ and $X^{\prime}_d$ are voltage dynamic parameters [Reference Schmietendorf, Peinke, Friedrich and Kamps18]. The nodes (islanded power grids or synchronous machines), which compose the network, are assumed to be all to all connected. For simplicity, we assume that the frequency controller acts instantaneously ( $\tau_{g}$ =0), such that it is described as the proportional derivative control given by:

$u_{i} =-\gamma_{i}\theta_{i} - \beta_{i}\omega_i$ . The derivative term $\beta_{i}\dot{\boldsymbol{\theta}}_{i}$ can be absorbed into the damping term $\alpha_{i}\dot{\boldsymbol{\theta}}_{i}$ of the swing equation. One can then also rewrite the set of equation (2.1) by equation (2.2) substituting $(X_{d,i}-X^{\prime}_{d,i})$ by $X_{i}$ , $T^{\prime}_{{d,i}_{0}}$ and $E^{\prime}_{q,i}$ by respectively $T_{{d,i}}$ and $E_{i}$ without loss of generality

(2.2) \begin{equation}\left\{ \begin{array}{ll} \dot{\boldsymbol{\theta}}_i=\omega_i\\[2pt] \dot{\omega}_i=-\alpha_{i}\dot{\boldsymbol{\theta}}_i-\gamma_{i}\theta_{i}+P_{i}^{*}-\sum\limits_{j=1}^N E_{i}B_{ij}E_{j}\sin(\theta_{i}-\theta_{j})\\[5pt] T_{{d,i}}\dot{\overline{E}}_{i}=E_{f,i}-E_{i}+ X_{i}\sum\limits_{j=1}^N B_{ij}E_{j}\cos(\theta_{i}-\theta_{j}).\\ \end{array}\right. \end{equation}

2.2 Linear stability analysis

The stable operation of the electric network requires that the frequency, voltage and difference of phases between connected nodes is constant, meaning that the network is operating in a synchronous regime. Such system is said to be linearly stable if, subjected to a small perturbation, it regains its stable operation, in the case of a power grid its synchronous state. Konstantin Sharafutdinov et al. have extensively studied the stability of such electric network model considering the voltage dynamics as described in equation (2.2), but without secondary control. Thus, based on their work and mainly on the necessary and sufficient conditions of the uncontrolled system to be linearly stable, the effects of the secondary control on the stability of the network will be investigated.

Thus, to analyse the stability of the system with respect to small perturbation, we linearise equation (2.2) around a steady state $(\theta_{i}^{*},\omega_{i}^{*},E_{i}^{*})$ . We denote small perturbation around the steady state as $\theta_{i}=\theta_{i}^{*}+\delta\theta_{i}$ , $\omega_{i}=\omega_{i}^{*}+\delta\omega_{i}$ , $E_{i}=E_{i}^{*}+\delta E_{i}$ . The linearisation of the equation (2.2) leads to the equation (2.3), where $\textbf{X}_{\textbf{1}}$ , $\textbf{X}_{\textbf{2}}$ and $\textbf{X}_{\textbf{3}}$ are n-dimensional vectors of $\delta\theta_{i}$ , $\delta\omega_{i}$ and $\delta E_{i}$ .

(2.3) \begin{align}\left\{ \begin{array}{l} \dot{\textbf{X}}_{\textbf{1}}= \textbf{X}_{\textbf{2}},\\[3pt] \dot{\textbf{X}}_{\textbf{2}}=-(\textbf{P}+{\boldsymbol\Gamma})\textbf{X}_{\textbf{1}}-\textbf{A}\textbf{X}_{\textbf{2}} - {\boldsymbol{\Lambda}}\textbf{X}_{\textbf{3}} \\[3pt] \dot{\textbf{X}}_{\textbf{3}}=\textbf{T}^{-1}{\boldsymbol{\chi}}{\boldsymbol{\Lambda}}\textbf{X}_{\textbf{1}} + \textbf{T}^{-1}({\boldsymbol{\chi}}\textbf{C}-\unicode{x1D7D9})\textbf{X}_{\textbf{3}}, \end{array}\right. \end{align}

where $\textbf{T}^{-1}$ , ${\boldsymbol\Gamma}$ , $\textbf{A}$ and ${\boldsymbol{\chi}}$ are diagonal matrices with elements $T_{ii}$ = $\frac{1}{T_{i}}$ , ${\boldsymbol\Gamma}_{\textbf{ii}}$ = $\gamma_i$ , $A_{ii}$ = $\alpha_i$ and $\chi_{ii}$ = $X_{i,}$ respectively, representing the relaxation time of the transient voltage dynamics matrix, the control, the damping and the transient reactance matrix of the synchronous machine. Matrices $\textbf{P}$ , ${\boldsymbol{\Lambda}}$ , $\textbf{C}$ , ${\boldsymbol\Gamma}$ $\in$ $\unicode{x211D}^{n \times n}$ whose the elements are respectively defined as follows:

(2.4) \begin{equation}P_{ij}=\left\{ \begin{array}{l@{\quad}l} -{E_{i}}^{*}B_{ij}{E_{j}}^{*}\cos({\theta_{i}}^{*}-{\theta_{j}}^{*}), & i \ne j,\\[5pt] \sum\limits_{l \ne i}^N {E_{i}}^{*}B_{il}{E_{l}}^{*}\cos({\theta_{l}}^{*}-{\theta_{i}}^{*}), & i=j, \end{array}\right.\end{equation}

(2.5) \begin{equation}{\Lambda}_{ij}=\left\{ \begin{array}{l@{\quad}l} {E_{j}}^{*}B_{ij}\sin({\theta_{i}}^{*}-{\theta_{j}}^{*}), & i \ne j,\\[5pt] -\sum\limits_{l= i}^N {E_{l}}^{*}B_{il}\sin({\theta_{l}}^{*}-{\theta_{i}}^{*}), & i=j, \end{array}\right.\end{equation}

(2.6) \begin{equation}{C}_{ij}=B_{ij}\cos({\theta_{i}}^{*}-{\theta_{j}}^{*}).\end{equation}

This later set of equations can be rewritten into the following compact form

(2.7) \begin{equation}\frac{\text{d}}{\text{d}t}\left( {\begin{array}{c@{\quad}c@{\quad}c} \textbf{X}_{\textbf{1}}\\[4pt] \textbf{X}_{\textbf{2}}\\[4pt] \textbf{X}_{\textbf{3}}\end{array} } \right)= \underbrace{\left[ {\begin{array}{ccc} \textbf{0} & \unicode{x1D7D9} & \textbf{0} \\[4pt] -(\textbf{P}+{\boldsymbol\Gamma}) & -\textbf{A} & - {\boldsymbol{\Lambda}}\\[4pt] \textbf{T}^{-1}{\boldsymbol{\chi}}{\boldsymbol{\Lambda}} & \textbf{0} &\textbf{T}^{-1}({\boldsymbol{\chi}}\textbf{C}-\unicode{x1D7D9})\end{array} } \right]}_{:=\,\textbf{J}}\left[ {\begin{array}{cc} \textbf{X}_{\textbf{1}}\\[4pt] \textbf{X}_{\textbf{2}}\\ \textbf{X}_{\textbf{3}}\end{array} } \right].\end{equation}

The evaluation of the Jacobian matrix $\textbf{J}$ around an existing fixed point $(\theta_i^{*}, \omega_i^{*}, E_i^{*})$ determines the stability of the system around this fixed point. That is done by computing the eigenvalues $\mu_{i}$ associated with this fixed point. The system is then linearly stable if the maximal real part of the eigenvalues is negative. In this case, perturbed, the system will regain its fixed point. But, if at least one eigenvalue has a positive real part, then, the system is said to be linearly unstable. This previous analysis is possible once we can determine the steady state of the system and therefore compute the eigenvalues.

In most of the cases, the fixed points can be difficult to determine, hence complicating the determination of eigenvalues. In such cases, the determination of condition of existence of the fixed point is necessary and the stability conditions of the system can be derived through some mathematical formulations.

Sharafutdinov et al. in [Reference Sharafutdinov, Rydin Gorjão, Matthiae, Faulwasser and Witthaut19] have provided for the uncontrolled system two sufficient and necessary stability conditions grouped into two propositions (proposition I and II in [Reference Sharafutdinov, Rydin Gorjão, Matthiae, Faulwasser and Witthaut19]), which both have to hold to have a linearly stable fixed point. In fact, it is worth noting that the Jacobian matrix $\textbf{J}$ has one eigenvector (1,0,0) with eigenvalue $\mu_{1}= 0$ , which corresponds to a global phase shift of the synchronous machine [Reference Manik, Witthaut, Schäfer, Matthiae, Sorge, Rohden, Katifori and Timme13, Reference Sharafutdinov, Rydin Gorjão, Matthiae, Faulwasser and Witthaut19]. This particular case, which does not have any physical implications for the stability, is excluded from the stability analysis. Thereby, we reduce the space of possible states to $S_{\bot}$ , which is the space perpendicular to the solution $\theta_i^{*}$ +c(1,0,0)/c $\in \unicode{x211D}$ and defined by

(2.8) \begin{equation} S_{\bot}=\left\{ (\textbf{X}_{\textbf{1}},\textbf{X}_{\textbf{2}},\textbf{X}_{\textbf{3}})\in \unicode{x211D}^{3n}/\textbf{1}\textbf{X}_{\textbf{1}}= 0\right\}.\end{equation}

The Proposition I in [Reference Sharafutdinov, Rydin Gorjão, Matthiae, Faulwasser and Witthaut19] states that a given steady state $(\theta_i^{*}, \omega_i^{*}, E_i^{*})$ is linearly stable if and only if

• • given the space $S_{\bot}^{1}=\left\{ (\textbf{X}_{\textbf{1}},\textbf{X}_{\textbf{2}}) \in \unicode{x211D}^{2N}/\textbf{1}\textbf{X}_{\textbf{1}}= 0\right\}$ , the matrix $\textbf{P}+{\boldsymbol\Gamma}$ is positive definite on $S_{\bot}^{1}$ , furthermore,

• • the matrix $\textbf{C}-{\boldsymbol{\chi}}^{-1}+ {\boldsymbol{\Lambda}}\textbf{P}^{+}{\boldsymbol{\Lambda}}^{T} $ is negative definite.

where T denotes the transpose of a matrix, which is identical to its inverse if it is an orthogonal matrix. In addition, the matrix $\textbf{P}^{+}$ represents the Moore-Penrose pseudo-inverse of $\textbf{P}$ . From this proposition I, the only difference appearing in comparison to the study made in Sharafutdinov et al. in [Reference Sharafutdinov, Rydin Gorjão, Matthiae, Faulwasser and Witthaut19] is the matrix ${\boldsymbol\Gamma}$ in the sum $\textbf{P}+{\boldsymbol\Gamma}$ . Thus, from the proposition I, it is clear that the secondary control represented by the matrix ${\boldsymbol\Gamma}$ can improve the stability of the network with the necessary condition to be positive definite. In our case, where the nodes are assumed to have the same system parameters such are $\alpha_{i}=\alpha$ , $\gamma_{i}=\gamma$ , $T_{d,i}=T_{d}$ and $E_{f,i}=E_{f}$ . The matrix ${\boldsymbol\Gamma}$ is a diagonal matrix, which elements $\gamma$ are always positive. Thus, all the eigenvalues of the symmetric matrix ${\boldsymbol\Gamma}$ are positive; hence, the matrix ${\boldsymbol\Gamma}$ is positive definite (See Appendix A).

2.3 Bulk dynamics

In order to analyse the stability of the network, we simplify the analysis by focusing on the ensemble dynamics, that is we consider the bulk or average dynamics of the network. As in the previous stability analysis, we wish to evaluate the impact of the secondary control through its parameter $\gamma$ on the dynamics of the system. Thus, we first consider a system with a constant voltage, which is a well known and studied case. Next, we include the voltage dynamics and analyse the effects of the secondary control on the voltage.

2.3.1 Constant voltage

Let’ s consider equation (2.2) and assume the voltage to be constant. The resulting equation is the simple Kuramoto model with secondary control, discussed for example in [Reference Tchuisseu, Gomila, Colet, Witthaut, Timme and Schäfer22]:

(2.9) \begin{equation}\left\{ \begin{array}{ll} \dot{\boldsymbol{\theta}}_i=\omega_i,\\ \dot{\omega}_i=-\alpha_{i}\dot{\boldsymbol{\theta}}_i-\gamma_{i}\theta_{i}+P_{i}^{*}-\sum\limits_{j=1}^N E_{i}B_{ij}E_{j}\sin(\theta_{i}-\theta_{j}).\\ \end{array}\right. \end{equation}

Now, we take the average of the variables of equation (2.9) over the number of nodes n that form the network. We further assume that the node parameters $\alpha_i=\alpha$ , $\gamma_i=\gamma$ are identical for all the nodes. Thus, we obtain the following differential equations:

(2.10) \begin{equation}\left\{ \begin{array}{ll} \dot{\overline{\theta}}=\overline{\omega},\\[5pt] \dot{\overline{\omega}}=-\alpha\dot{\overline{\theta}}-\gamma\overline{\theta}+\frac{1}{N}\sum_{i=1}^N {P_{i}^{*}}.\\ \end{array}\right. \end{equation}

Equation (2.10) can easily be solved, and the average value of the principal variable is obtained as follows:

\begin{equation*} \overline{\theta}(t)=C_{1}\exp{r_1t}+ C_{2}\exp{r_2t}+\frac{\sum_{i=1}^N {P_{i}^{*}}}{N\gamma},\end{equation*}

where $C_1$ and $C_2$ are constants which are determined by the initial conditions and $r_1$ and $r_2$ are expressed as follows:

(2.11) \begin{equation}\left\{ \begin{array}{ll} r_1=-\dfrac{\alpha}{2}+\dfrac{\sqrt{\alpha^2 - 4\gamma}}{2},\\[12pt] r_2=-\dfrac{\alpha}{2}-\dfrac{\sqrt{\alpha^2 - 4\gamma}}{2}.\\ \end{array}\right.\end{equation}

Thus, it is clear that for time tending to infinity and for any value of the secondary control parameter $\gamma$ ( $\gamma>0$ ), the average value of the angle of rotation of each node in the network is constant and the average value of the frequency is zero:

\begin{equation*} \overline{\theta}(t)=\frac{\sum_{i=1}^N {P_{i}^{*}}}{N\gamma}, \quad \bar{\omega}=0.\end{equation*}

As we can see, the average of the angle $\overline{\theta}$ for the controlled system depends on the control parameter $\gamma$ , the balance of the network as well as the size of the network. For a balanced network, that is with $\sum_{i=1}^N {P_{i}^{*}}=0$ , the mean value of the angles is zero. For an imbalanced and fixed network, the mean angle is inversely proportional to the secondary control parameter; thus for a large control, the average of the angles tends to zero.

For an uncontrolled electric network ( $\gamma=0$ ), the average of the angles is time varying and given by the following expression:

\begin{equation*} \overline{\theta}(t)=D_{1}+ D_{2}\exp{(-\alpha t)}+\frac{\sum_{i=1}^N {P_{i}^{*}}}{N\alpha}t,\end{equation*}

where $D_1=-D_2=-{\sum_{i=1}^N{{P_{i}^{*}}}/{N\alpha^2}}$ are constants determined by the initial conditions. Thus, the mean value of the frequency of the network as a function of time t is given by the equation (2.12):

(2.12) \begin{equation} \overline{\omega}(t)= -\alpha{D_{2}}\exp{(-\alpha t)}+\frac{\sum_{i=1}^N {P_{i}^{*}}}{N\alpha}, \end{equation}

this basically corresponds to the frequency deviation for an uncontrolled and imbalanced network.

This imbalance can be easily absorbed in large networks (large N) and those with large primary control (given by the $\alpha$ parameter).

We have so far shown that, when the voltage is constant, any imbalance in the network induces the deviation of the mean frequency, which is constant and non-zero for a network without control. This implies that the average value of the angle is linearly increasing with the time. For the same network, but in presence of secondary control in all the nodes, one notices that the average frequency in the network is always zero after a long period of time and the average angular in the network is a constant value proportional to the disturbance and inversely proportional to the size of the network and the secondary control parameter. We have therefore seen the effects of the secondary control in such an electric network in which the voltage dynamics is not considered. But how does the voltage dynamics change these results?

2.3.2 Dynamical voltage

In this part, we investigate the effects of the secondary control on the dynamics of the voltage in an electric network. Thus, we consider the full equation (2.2), summing up all three equations and dividing them by the total number of nodes. Again, we assume that the parameters of the nodes are constant and identical for all the nodes, obtaining the following equations:

(2.13) \begin{equation}\left\{ \begin{array}{lll} \dot{\overline{\theta}}=\overline{\omega}\\[4pt] \dot{\overline{\omega}}=-\alpha\dot{\overline{\theta}}-\gamma\overline{\theta}+\frac{1}{N}\sum_{i=1}^N {P_{i}^{*}}\\[6pt] T_{d}\dot{\overline{E}}=E_f -\overline{E} + \frac{X}{N}\sum_{i=1}^N\sum\limits_{j=1}^N B_{ij}E_{j}\cos\left(\theta_{i}-\theta_{j}\right) \end{array}\right.\end{equation}

As shown in equation (2.13), the equations describing the dynamics of the mean voltage (the third equation) contain the difference of phases of connected nodes as argument of a cosine function. The presence of this term makes the determination of the mean voltage fixed point very difficult. Nevertheless, one can find the range of variation of $\overline{E}$ considering the following relations:

\begin{equation*} \begin{cases} B_{ij}=B_{ji}=\left\{ \begin{array}{l@{\quad}l} B_{0} & \mbox{if } i=j \\ B_{1} & \mbox{else} \end{array}\right.\\[14pt] -1 \leq \cos{(\theta_{i}-\theta_{j})}\leq 1. \end{cases}\end{equation*}

Thus, the term with cosine in equation (2.13) is approximated by the following relation:

(2.14) \begin{equation}\begin{cases}-X(B_{0}+(N-1)B_{1})\overline{E}\leq \frac{X}{N}\sum_{i=1}^N\sum\limits_{j=1}^N B_{ij}E_{j}\cos(\theta_{i}-\theta_{j}), \\[6pt] \frac{X}{N}\sum_{i=1}^N\sum\limits_{j=1}^N B_{ij}E_{j}\cos(\theta_{i}-\theta_{j})\leq X(B_{0}+(n-1)B_{1})\overline{E}.\end{cases}\end{equation}

The range of variation of the mean value can then be determined, solving the following inequalities:

(2.15) \begin{equation} \begin{cases}T_{d}\dot{\overline{E}}+\overline{E}-E_f\leq X[(B_{0}+(N-1)B_{1})]\ \overline{E},\\[6pt] T_{d}\dot{\overline{E}}+\overline{E}-E_f \geq X[(B_{0}-(N-1)B_{1})]\ \overline{E}. \end{cases} \end{equation}

Thus, the average of the voltage is bounded as follows:

(2.16) \begin{equation} \begin{cases} \overline{E\left(t\right)} \leq E_{0}\exp{\left(\dfrac{-\left(1-\sigma_1\right)t}{T_{d}}\right)}+\dfrac{E_{f}}{T_{d}\left(1-\sigma_1\right)}\left(\ 1-\exp{\left(\dfrac{-\left(1-\sigma_1\right)t}{T_{d}}\right)}\right)\,\\[16pt] \overline{E\left(t\right)} \geq E_{0}\exp{\left(\dfrac{-\left(1+\sigma_2\right)t}{T_{d}}\right)}+\dfrac{E_{f}}{T_{d}\left(1+\sigma_2\right)}\left(\ 1-\exp{\left(\dfrac{-\left(1+\sigma_2\right)t}{T_{d}}\right)}\right)\ , \end{cases} \end{equation}

where $E_{0}$ is the initial value of the mean voltage: $\overline{E}(t=0)=E_{0}$ , and $\sigma_2$ and $\sigma_1$ are given by:

\begin{equation*} \begin{cases} \sigma_1 =X[(B_{0}+(N-1)B_{1})],\ \\[5pt] \sigma_2 =-X[(B_{0}+(N-1)B_{1})].\ \end{cases}\end{equation*}

Evaluating all previous inequalities, we observe that the dynamics of the mean voltage is strictly related to the size of the network. Thus, according to equation (2.16), the mean voltage will be bounded between real values if and only if:

(2.17) \begin{equation} \begin{cases} 1-\sigma_1 =1-X[(B_{0}+(N-1)B_{1})]\geq 0, \\[5pt] 1-\sigma_2 =1+X[(B_{0}-(N-1)B_{1})] \geq 0.\\ \end{cases}\end{equation}

Thus, the stability of the electric voltage is achieved if the network size and the network coupling are chosen in the range given by the following inequality:

(2.18) \begin{equation} 1 -X\frac{1-B_{0}}{B_{1}}\leq N\leq 1+X\frac{1-B_{0}}{B_{1}}\end{equation}

In the present study, the values of the system’s parameters are set as: $X=1$ , $B_{0}=-0.8$ and $B_{1}=1$ , following literature values [Reference Auer, Kleis, Schultz, Kurths and Hellmann1, Reference Schmietendorf, Peinke, Friedrich and Kamps18]. Thus, the mean voltage of the network of size N will be bounded if and only if:

(2.19) \begin{equation} {N= 2}\end{equation}

Indeed, let’s consider an electric network with the same parameters as the studied case and composed of n nodes, where each node represents a small region or a city. The result above states that an extension of the current network composed of $N=2$ nodes by connecting it with additional nodes leads to the instability of the new extended network in terms of voltage.

Thus, from this bulk analysis, it appears that the secondary frequency control has explicitly no effects on the stabilisation of the mean voltage, but does have an effect when it comes to stabilising the mean frequency of the network. This implies therefore the necessity of having another form of control for the mean voltage in larger electric networks.

3.1 2-node system

The two node system considered here consists of a generator node of power $P_{1}>0$ and a consumer node of power $P_{2}<0$ . The two nodes represent simply two interconnected power grids (or synchronous machines) with equal control parameters $\gamma_{1}=\gamma_{2}=\gamma$ and are described by the equation (2.2) and parameters, such as $T_d$ , inspired by literature values [Reference Auer, Kleis, Schultz, Kurths and Hellmann1, Reference Schmietendorf, Peinke, Friedrich and Kamps18]. Initially, without disturbance, the two systems evolve till reaching a stable final state, in which the variables of each node tend to a steady state as shown in Figure 1 without control (first column in black) and with control (second column in red). One observes that in both cases, the network is synchronised, meaning that the frequency at each node is equal to $\omega_{syn}= 50$ Hz and the difference of phases between the connected node tends to a constant value. The final states reached by the voltages in the controlled case are sightly greater than the ones of the uncontrolled case. This slight increase results from the effects of the secondary control parameter at each node, which affects the dynamics of the voltage through the angles. In fact, at a constant voltage, the secondary control is basically the governor which activates the online or offline substations, thereby increasing or reducing the power generated in order to balance the power in the system ( $\omega=0$ ) [Reference Tchuisseu, Gomila, Colet, Witthaut, Timme and Schäfer22–Reference Tchuisseu, Gomila and Colet24]. Thus, once the frequency is brought back to its nominal value, the phases in each node tend to constant values, which basically implies a constant difference of phases between connected nodes, hence constant voltage. Without control, the two nodes of the system evolve until sharing a constant amount of power and reach a power balanced state. Next, we consider a perturbation to the power system: We assume that from time $t=40$ s to $t=42$ s, the power $P_{1}$ at node 1 experiences a gradual increase of power from 0 to $P_{dist}=1$ , resulting in a sudden increase of its corresponding frequency, as shown in Figure 2 (dashed lines) in the controlled (red) and uncontrolled (black) case. This perturbation leads to an instability of the uncontrolled case, visible by the emergence of periodic oscillation in the frequency and the voltages around states far away from their original states, while the corresponding phases are rising in opposite directions. On the other hand, the variables in the controlled system (plot in red in Figure 2) converge to new steady states for the voltages and phases, while the corresponding frequencies return to their synchronous states as before the perturbation. Thus, the secondary control not only stabilised the frequencies but also stabilised the voltage, thereby avoiding the system to drop into an unstable regime. Thus, the secondary frequency control acts as a damping for the voltages. This is a very similar behaviour as the primary frequency control has on the frequency through the damping coefficient $\alpha$ on the frequency of the uncontrolled ( $\gamma=0$ ) network. Hence, we can consider secondary frequency control acts somehow as a primary control for the voltage.

Figure 1.

An unperturbed system converges to a steady state with and without control. [Description of figure] parameters of simulation are the following: $\alpha_{1}=\alpha_{2}=0.2$ , $T_{d,1}=T_{d,2}=2$ , $E_{f,1}=E_{f,2}=1.$ , $E_{1}(0)=E_{2}(0)=1.14$ , $\omega_{1}(0)=\omega_{2}(0)=50$ , $\theta_{1}(0)=\theta_{2}(0)=0$ , $P_{dist}$ =0, $\gamma_{1}=\gamma_{2}= 0$ (first column), $\gamma_{1}=\gamma_{2}= 1$ (second column).

Figure 2.

Secondary control stabilises a 2 node system. [Description of figure] parameters of simulation are the following: $\alpha_{1}=\alpha_{2}=0.2$ , $T_{d,1}=T_{d,2}=2$ , $E_{f,1}=E_{f,2}=1.$ , $E_{1}(0)=E_{2}(0)=1.14$ , $\omega_{1}(0)=\omega_{2}(0)=50$ , $\theta_{1}(0)=\theta_{2}(0)=0$ , $P_{dist}$ =1, $\gamma_{1}=\gamma_{2}= 0$ (first column), $\gamma_{1}=\gamma_{2}= 1$ (second column).

Figure 3.

Uncontrolled networks display frequency instability and large voltage oscillations. Illustration of average of the phase (a), frequency deviation (b) and the voltage (c) as a function of the time for different size of the network, with the corresponding zoom of the voltage dynamics given at the bottom of the figure. Parameter values: $\gamma=0$ , $\alpha=0.2$ , $Pdist=1$ , $X=1$ , $E_{f}=1$ , $E_{i}(t=0)=1.14$ , $T_{d}=1$ , $B_{0}=-0.8$ and $B_{1}=1$ .

3.2 Bulk dynamics

Moving towards the $N>2$ , case, we now present the numerical results complementing the analytical bulk analysis. Thus, the question to be answered is the following: how does the secondary control affect the mean voltage dynamics in networks? Thus, we consider the electric network, whose dynamics is governed by the equation (2.2), and we analyse the dynamics of the mean voltage and frequency for different network sizes. We consider for these purposes that we have an interconnected network constituted of n nodes, which can represent n interconnected isolated power grids or synchronous machines. So far, we assumed that the nodes are all-to-all-coupled, each node having equal power $P_{i}$ in absolute value, and the network without perturbation is power balanced ( $\sum_i {P_{i}=0}$ ). Such all-to-all-coupled networks naturally emerge after Kron reduction of any network topology [Reference Dorfler and Bullo5]. In all the studied scenarios, the system is perturbed from the time $t=40$ s to $t=42$ s by gradually increasing the power at a single node from 0 to $P_{dist}$ .

First, we consider the uncontrolled ( $\gamma=0$ ) case and plot the deviation of the mean phase, frequency and voltage as function of time for different network size, and without secondary control in Figure 3. We observe a clear agreement with the analytical predictions given by the expression of the frequency deviation in equation (2.12), meaning that the mean frequency for a long time tends to a constant. Since we assumed that the primary control parameter/the inertia of each power plant is constant, we observe that the mean of the frequency deviation decreases when the size of the network increases. The corresponding curve of the average voltage is plotted at the third column of the first row. Complementing the analytical results, we observe large transient voltage dynamics, which increase in amplitude with network size. To better observe the dynamics of the mean value of the voltage towards the end of our simulation window, we provide zooms in the lower row of Figure 3 for each size n of the network. For the two nodes system, the average frequency tends to a constant value. The mean voltage conversely is oscillating around a stable value $\overline{E} \approx 0.55$ V. For larger values of n on the other hand, the mean voltage is fluctuating and displays large periodic peaks. These regular peaks persist throughout the running time of our simulation and are observed for N $\textbf{=}$ 10 and N $\textbf{=}$ 20. For N $\textbf{=}$ 50, after the transient, the mean voltage mostly fluctuates around a constant value and only displays one peak, likely pointing to a longer periodicity of these peaks due to the larger network. The small fluctuations observed are probably due to the difference of phases existing between the connected nodes. In this uncontrolled studied case, the system is clearly not a synchronous stable electric network.

Now, we include secondary control ( $\gamma>0$ ) and repeat the same simulations as before. Figure 4 (top) shows the time evolution of the mean frequency and the mean voltage of the considered network in presence of the secondary control. First, we notice that the mean frequency deviation tends to zero for every network size, as we predicted in our linear stability analysis. Hence, we focus on the mean voltage dynamics, for which we could not derive an equality from the equation (2.13) but only constrained its range of fluctuations. Thus, only these numerical simulation can tell us how the voltage evolves with time. Similar to the frequency, after an initial transient phase, the mean voltage tends to a constant value for all the considered sizes of the network.

Figure 4.

Secondary control stabilises the frequency and voltage following initial oscillations. Illustration of average of the phase (a), frequency deviation (b) and the voltage (c) as a function of the time for different size of the network, with the corresponding zoom of the voltage dynamics given at the bottom of the figure. Parameter values: $\gamma=1$ , $\alpha=0.2$ , $Pdist=1$ .

Again, to highlight the voltage dynamics towards the end of our simulation window, we provide a zoom of the voltages in the second row of Figure 4. For the two node system, we observe that the mean voltage tends to steady state, and it is not oscillating as in the corresponding uncontrolled case. Thus, the secondary control clearly stabilises the voltage. In fact, without secondary control, the phase of each node is increasing continuously with the time (see Figure 2), leading to the oscillation of the voltage. The secondary frequency control on other hand stabilises the phases; hence, the cosine function in the expression of each voltage becomes a constant, leading to constant voltage at each node. For larger networks, it may happen that due to the power flow between connected nodes, the phases slightly vary, leading at some times to the fluctuations of the voltage. This justifies some rare peaks observed in the mean voltage. Nevertheless, the magnitude of the voltage and its peaks in the controlled network mostly remains lower than the one in the uncontrolled network. In addition, we observe that the steady state reached in the end decreases when the network size increases.

3.3 Relaxation time

We have shown in the previous analysis how the dynamics of the phase, the frequency, the voltage as well as their corresponding mean values evolve with the secondary parameter $\gamma$ , once the network is perturbed. In fact, once perturbed from its stable state, the variables of the system vary until they reach new steady states, often different from their original steady states without perturbation. Now, we investigate how long it takes the system to relax to its (new) steady state, by computing the return time or relaxation time. According to the bulk dynamics of the network, the return time of the frequency and phase average depends strongly on the damping of the system, and slightly on the secondary control parameter, which impacts mostly the oscillatory regime of these variables. If the return time of the mean frequency and phase can be obtained explicitly, this cannot so easily be achieved for the mean value of the voltage. Hence, we again use numerical computations to obtain the return time for different control parameters.

Figure 5.

Secondary control stabilises the system by reducing the return time. We plot the return time of the mean voltage when the control parameter increases in the two node system. The parameters of simulation are the following: $\alpha_{1}=\alpha_{2}=0.2$ , $T_{d,1}=T_{d,2}=2$ , $E_{f,1}=E_{f,2}=1$ , $P_{dist}=1$ $N=2$ nodes.

In particular, we consider the previously described 2 nodes system, including voltage dynamics. For this system, we aim to determine the evolution of the return time as a function of the secondary control parameter $\gamma$ . The network is perturbed as previously, by gradually decreasing the power at the node 1 from $P_1$ to $P_{1}+P_{dist}$ ( $P_{dist}=-1$ ) from time $t=40$ s to $t=42$ s and all other parameters also remain unchanged. The relaxation time or return time is numerically computed by recording the time taken by the voltage at a node to regain a steady state after perturbation. Thus, we define by $\overline{E}(t)$ the mean voltage at the time t and $\overline{E}(t-T)$ the mean at the time $t=t-T$ , where T is a characteristic time. We define also by $\xi$ the numerical tolerance. The mean voltage is considered stable after perturbation of the network if and only if $|\overline{E}(t)-\overline{E}(t-T)|\leq$ $\xi$ . Note that we focus here on the relaxation and therefore of the voltage, not the frequency, as we are interested in quantifying the impact of frequency control on voltage stability. The so computed return times decrease with increasing secondary control parameter $\gamma$ , see Figure 5. First, we note that for any non-zero value of the secondary control parameter $\gamma$ , the system after perturbation always regains a steady state, that is the secondary frequency control guarantees stability and return. Secondly, we observe that the return time is reduced by increasing the control amplitude, approaching zero for sufficiently large control $\gamma$ . This means that the secondary frequency control is strong enough to immediately compensate the introduced disturbance and the system never leaves its original fixed point.

3.4 Case study: Heterogeneous parameters and network influence

So far we have assumed networks that are all-to-all-coupled and have homogeneous parameters, that is identical loads, generation and control values. We used this simplification to derive analytical results. To demonstrate that our results are in principle also applicable to more general systems, we consider one non-all-to-all-coupled system with heterogeneous parameters, see Figure 6(a). In particular, we simulate the dynamical behaviour of 20 nodes connected to a common bus with power values P randomly drawn between $-0.7$ and 0.7. We assume nodes with large absolute value of P to provide more control power and hence set the control values proportional to the absolute power value: $\alpha_i=0.2|P_i|$ , secondary control $\gamma_i=|P_i|$ , where $|...|$ denotes the absolute value. Analogue to the earlier analysis, we observe that the average frequency can only be restored once secondary control is used, see Figure 6(b). Meanwhile, the voltage dynamics is not controlled and oscillates (Figure 6(c)) since we have more than two nodes, consistent with our analytical results. How exactly non-homogeneous parameters and network topology affect the dynamics and controlability of both frequency and voltage is beyond the scope of this study.

Figure 6.

Exploring heterogeneous parameters and non-all-to-all-coupling. Illustration of the topology (a), the average of the frequency deviation (b) and the voltage (c) as a function of the time for the controlled (black) and controlled (red) network. Parameter values: $P=[-0.4, 0.53, -0.51, 0.56,$ $0.52, 0.48, -0.55, -0.45, 0.491, 0.509, -0.482, -0.518, -0.46,-0.64, 0.42, 0.58, -0.5, 0.5, 0.35, -0.45]$ , primary control $\alpha_i=0.2|P_i|$ , secondary control $\gamma_i=|P_i|$ , while $P_{dist}$ , X and B are chosen as in the other simulations.