MATLAB/Simulinkでフェランチ効果をシミュレーションする

# はじめに
以前に[潮流解析](https://qiita.com/BLITZ1127/items/f931c5e98f8fe9ffb821)を実施しました。そこでは、Power GuiのLoad Flow Analyzer を活用しました。本記事では、潮流解析にもう少し踏み込み、フェランチ効果を確認してみたいと思います。

# フェランチ効果とは
フェランチ効果とは、送配電線路において、軽負荷時や無負荷時（主に夜間）などに受電端電圧が送電端電圧よりも高くなる現象です。通常の送配電線路では電圧降下が発生するために、送電端電圧よりも受電端電圧の方が低くなっています。しかし、フェランチ効果では、長距離送電線ケーブルのような長い系統のような場合において電線路によっては受電端電圧よりも送電端電圧の方が高くなることがあります。その様子について、検討していきたいと思います。

# 受電端電圧の解析解の導出
まず、送配電線路を等価回路モデルに置き替えます。その様子を図1に示します。前回の記事では、送電線の抵抗を無視していましたが、今回は厳密に計算するために記載することにします。この時の受電端の特性については、前回同様の考えで計算することができます。

![ferranti_effect_circuit.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3035813/2b0edbc5-3286-c644-d396-0d50ae846a52.png)
図1　送配電線路の等価回路モデル


まず、送電線の電圧降下に注目します。電圧降下は以下の(1)式で表されます。ここで、$\dot{A}$のようなマークは複素数表示を表すこととします。
　$\dot{V_s} －\dot{V_r} =\dot{I}×(R+jX) （1）$
ただし、相差角を考慮すると、(2), (3)式の関係があります。
　$\dot{V_s} =V_s \cos{\delta} +j V_s \sin{\delta} 　(2)$
　$\dot{V_r}  =V_r　(3)$
したがって、(1)式は(4)式のように展開できます。
　$\dot{I}  =\frac{\dot{V_s}－\dot{V_r}}{R+jX}=\frac{V_s \cos{\delta} +j V_s \sin{\delta} -V_r}{R+jX}
=\frac{RV_s \cos{\delta} -RV_r+XV_s \sin{\delta}}{R^2+X^2} + j\frac{-XV_s \cos{\delta} +XV_r+RV_s \sin{\delta}}{R^2+X^2}　(4)$
次に、受電端を考えます。受電端での皮相電力は以下の(5)式で表されます。ここで$\bar{A}$  のようなマークは複素数表示の複素共役を表すこととします。
　$S =\bar{\dot{I}} ×\dot{V_r} 　(5)$
先ほどの(4)式から(6)式のように展開できます。
　$S=\bar{\dot{I}} ×\dot{V_r} 
=\frac{RV_s \cos{\delta} -RV_r+XV_s \sin{\delta}}{R^2+X^2}V_r + j\frac{XV_s \cos{\delta} -XV_r-RV_s \sin{\delta}}{R^2+X^2}V_r
　(6)$
したがって有効電力$P$と無効電力$Q$は以下の(7),(8)式で表されます。
　$P = Re{(S)} =\frac{RV_s \cos{\delta} -RV_r+XV_s \sin{\delta}}{R^2+X^2}V_r (7)$
　$Q =Im{(S)} = \frac{XV_s \cos{\delta} -XV_r-RV_s \sin{\delta}}{R^2+X^2}V_r(8)$

せっかくMATLABを使用しているので、$V_r$および$\delta$について、$V_s,P,Q,X,R$に関する式として解いてしまいます。symbolic math toolboxを使用します。具体的には以下のようなコードを書きます。
```
eq1 = P== (R*Vs*cos(delta)-R*Vr+X*Vs*sin(delta))*Vr/(R^2+X^2);
eq2 = Q== (X*Vs*cos(delta)-X*Vr-R*Vs*sin(delta))*Vr/(R^2+X^2);
eqns = [eq1 eq2];
sol_param = solve(eqns,Vr,delta);
```
これで、sol_paramに解が導入されます。
具体的な値を代入し、結果を得るには以下のようにします。
```
sol_param_sub = subs(sol_param,Vs,1.0);
sol_param_sub = subs(sol_param_sub,P,0.5);
sol_param_sub = subs(sol_param_sub,Q,0);
sol_param_sub = subs(sol_param_sub,X,0.5);
sol_param_sub = subs(sol_param_sub,R,0.1);
Vrs = vpa(sol_param_sub.Vr);
deltas = vpa(sol_param_sub.delta);
```
このようにすることで、$V_s=1.0$pu, $P=0.5$pu, $Q=0$pu, $X=0.5$pu, $R=0.1$puのときの$V_r$と$\delta$を得ることができます。

# 解析解のプロット
　先ほどの解析解を用いて、プロットします。変数がいくつかあるので、$P=0.5$pu, $V_s＝1.0$pu, $R=0.1$pu, $X=0.5$puと設定します。この時に、$Q$を変動させたときの$V_r$および$δ$の様子をプロットします。-1.0~1.0puまでを$Q$の変動範囲とします。なお、$Q$は無効電力を消費する向きを正としています。
　その結果を図2に示します。$Q$のレンジはおよそ+0.3puまでとなっています。これより大きくすると、解が存在しなくなるためです。$Q$を進み方向（負の方向）に大きくすると電圧が上昇しているのが分かるかと思います。グラフから、$Q=-0.16$pu付近で$V_s=V_r=1.0$puとなっており、それより$Q$が小さい場合には送電端よりも受電端の電圧が上昇していることがわかると思います。これがフェランチ効果です。また、$δ$に関しては、無効電力を消費するほど大きくなっていることが分かります。

![Plot.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3035813/a4255ab9-c6cc-7f57-0414-d82a74fe9882.png)

図2　無効電力と受電端電圧、相差角の関係

# MATLAB/Simulinkによるフェランチ効果の解析
　それでは、Power GuiのLoad Flow Analyzer を活用し、解析解と一致しているか確認します。回路は図3のとおりとします。また、解析条件は表1のとおりとします。条件1は基本の運転点を、条件2は負荷が運転しており無効電力を消費している状況を、条件3は進相コンデンサなどが投入され、無効電力が負荷から供給されている状況を示しています。

![スクリーンショット_20221226_134656.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3035813/b8896cdd-02e2-db15-840c-1f88f7b0ad53.png)
図3　シミュレーション回路図

表1　回路パラメータ

||単位法|非単位法|
|:-:|:-:|:-:|
|電源電圧|1pu|66kV|
|周波数||60Hz|
|送電線抵抗 R|0.1pu|4.356Ω|
|送電線リアクタンス X|0.5pu|21.78Ω(57.775mH)|

表2　解析条件

||有効電力 P|無効電力 Q|
|:-:|:-:|:-:|
|条件1|50MW(0.5pu)|0Mvar(0.0pu)|
|条件2|50MW(0.5pu)|-20Mvar(-0.2pu)(遅れ)|
|条件3|50MW(0.5pu)|+20Mvar(0.2pu)(進み)|


# 解析結果の比較
　それでは、Simulinkで計算した値を解析解と比較しましょう。比較結果を表3に示します。Load Flow Analyzerでの計算結果が解析解と一致していることが分かります。符号が違うのは、向きの取り方の問題です。また、相差角について若干の誤差がありますが、これについては詳細原理を把握できていません。
　解析解は汎用性が高いことがメリットです。一方、Load Flow AnalyzerはGUIで直感的にわかりやすいので活用しやすく、特定の系統においての潮流を計算するために優れたツールであるといえます。

表3　解析結果

||受電端電圧$V_r$|相差角$\delta$|
|:-:|:-:|:-:|
|条件1　解析解|0.9060pu|16.0182deg|
|条件1　シミュレーション解|0.9060pu|-16.0184deg|
|条件2　解析解|0.7530pu|17.7847deg|
|条件2　シミュレーション解|0.7530pu|-17.7854deg|
|条件3　解析解|1.0132pu|15.4553deg|
|条件3　シミュレーション解|1.0132pu|-15.4543deg|


# まとめ
本記事では、フェランチ効果について解説しました。フェランチ効果は受電端の電圧が送電端の電圧よりも高くなる現象であり、その原因は無効電力にあることが分かりました。また、フェランチ効果（電圧上昇）がどの程度発生するか、解析解を算出するとともに、Load Flow Analyzerでも計算できることを確認しました。