LTspiceで三相PWM整流器を作る～d-q座標系編～

はじめに

秋分を過ぎ，本格的に季節が巡って参りましたね🍂 前記事よりまた間が空いてしまいましたが，三相PWM整流器の直流電圧制御について，私の理解が非常に浅かったことが分かり，改めて考え直していたところです。本記事では，まず，前記事で宿題となっていたd-q座標系での等価回路を導出します。本題である直流電圧制御系とその安定性については，長くなりそうなので次記事に譲ることとします。

何を導出したいか

ここでは，系統と変換器が抵抗とインダクタンスを介して向き合っている三相回路をd-q座標系で表現した等価回路を導出します。図1に三相PWM整流器（単線図）とフェーザ図を示します。

f:id:s-inoue2010:20200929174821p:plain:w420 — 図1: 三相PWM整流器とフェーザ図

スイッチング周波数成分を無視し，電圧形インバータ部分も正弦波と見なして3線図とすれば，図2(a)の等価回路を得ます。図2(b), (c)はd-q座標系での等価回路（d軸・q軸等価回路）です。これはどのように導出するのでしょうか？ d軸とq軸の干渉項である−ωLi_qやωLi_dは一体どこから出てくるのか…？

f:id:s-inoue2010:20200929175041p:plain:w600 — 図2: 三相等価回路とd-q座標系での等価回路

なお，本記事では三相等価回路である図2(a)から出発し，座標変換にて図2(b), (c)が得られることを示しますが， α-β変換（Clarke変換）やd-q変換（Park変換）そのものの概念や効用については記載していません…💦 アニメーションを使った説明が秀逸なMathWorksさんの記事へのリンクを貼っておきます。

jp.mathworks.com

三相座標系からd-q座標系へ

前記事にも記載しましたが，図2(a)では次式が回路方程式として成り立ちます。

$\displaystyle \begin{bmatrix} v_{u} \\ v_{v} \\ v_{w} \end{bmatrix} = R \begin{bmatrix} i_{u} \\ i_{v} \\ i_{w} \end{bmatrix} + L \frac{d}{dt} \begin{bmatrix} i_{u} \\ i_{v} \\ i_{w} \end{bmatrix} + \begin{bmatrix} v_{su} \\ v_{sv} \\ v_{sw} \end{bmatrix} \tag{1}$

(1)式は実は図1(b)のフェーザ図と同じことを言っています。少々天下りですが，零相成分を含めた逆α-β変換と逆d-q変換*1を施す行列をそれぞれ[C_αβ]，[C_dq]とします（「逆」変換で行列を定義しておくと，以下の数式展開が自然に進みます）。

$\displaystyle [\boldsymbol{C}_{\alpha \beta} ] = \sqrt{\frac{2}{3}} \begin{bmatrix} \frac{1}{\sqrt{2}} & 1 & 0 \\ \frac{1}{\sqrt{2}} & -\frac{1}{2} & \frac{\sqrt{3}}{2} \\ \frac{1}{\sqrt{2}} & -\frac{1}{2} & -\frac{\sqrt{3}}{2} \end{bmatrix} \tag{2}$

$\displaystyle [\boldsymbol{C}_{dq} ] = \begin{bmatrix} 1 & 0 & 0 \\ 0 & \cos \theta & -\sin \theta \\ 0 & \sin \theta & \cos \theta \end{bmatrix} \tag{3}$

したがって，(2)式と(3)式の積は

$\displaystyle [\boldsymbol{C}] = [\boldsymbol{C}_{\alpha \beta} ] [\boldsymbol{C}_{dq} ] = \sqrt{\frac{2}{3}} \begin{bmatrix} \frac{1}{\sqrt{2}} & \cos \theta & -\sin \theta \\ \frac{1}{\sqrt{2}} & \cos \left (\theta - \frac{2 \pi}{3} \right ) & -\sin \left (\theta - \frac{2 \pi}{3} \right ) \\ \frac{1}{\sqrt{2}} & \cos \left (\theta + \frac{2 \pi}{3} \right ) & -\sin \left (\theta + \frac{2 \pi}{3} \right ) \end{bmatrix} \tag{4}$

と書けます。[C]を用いると，三相座標系とd-q座標系の変換を簡潔に書くことができ，

$\displaystyle \begin{bmatrix} v_{u} \\ v_{v} \\ v_{w} \end{bmatrix} = [\boldsymbol{C}] \begin{bmatrix} v_{0} \\ v_{d} \\ v_{q} \end{bmatrix} \tag{5}$

$\displaystyle \begin{bmatrix} v_{su} \\ v_{sv} \\ v_{sw} \end{bmatrix} = [\boldsymbol{C}] \begin{bmatrix} v_{s0} \\ v_{sd} \\ v_{sq} \end{bmatrix} \tag{6}$

$\displaystyle \begin{bmatrix} i_{u} \\ i_{v} \\ i_{w} \end{bmatrix} = [\boldsymbol{C}] \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} \tag{7}$

と表現できます。 (5), (6), (7)式を用いて(1)式を書き直せば，次式を得ます（「逆」変換で行列を定義しておいたメリットがここで活きてきます）。

$\displaystyle [\boldsymbol{C}] \begin{bmatrix} v_{0} \\ v_{d} \\ v_{q} \end{bmatrix} = R [\boldsymbol{C}] \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} + L \frac{d}{dt} \left ( [\boldsymbol{C}] \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} \right ) + [\boldsymbol{C}] \begin{bmatrix} v_{s0} \\ v_{sd} \\ v_{sq} \end{bmatrix} \tag{8}$

(8)式の両辺に[C]の逆行列[C]⁻¹を左から掛けると，

$\displaystyle \begin{bmatrix} v_{0} \\ v_{d} \\ v_{q} \end{bmatrix} = R \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} + L [\boldsymbol{C}]^{-1} \frac{d}{dt} \left ( [\boldsymbol{C}] \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} \right ) + \begin{bmatrix} v_{s0} \\ v_{sd} \\ v_{sq} \end{bmatrix} \tag{9}$

となり，[C]のほとんどが消えますが，右辺第2項にまだ残っています。電流[i₀, i_d, i_q]^Tが時刻tに依存するのはもちろんのこと， [C]もθ = ωtを含んでおり，やはり時刻tに依存します。したがって，(9)式右辺第2項を積の時間微分 *2として

$\displaystyle L [\boldsymbol{C}]^{-1} \left \{ \left ( \frac{d}{dt} [\boldsymbol{C}] \right ) \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} + [\boldsymbol{C}] \left ( \frac{d}{dt} \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} \right ) \right \} \tag{10}$

と書いて，さらに，L[C]⁻¹をカッコ内に入れ込むと，(10)式は

$\displaystyle \left ( L [\boldsymbol{C}]^{-1} \frac{d}{dt} [\boldsymbol{C}] \right ) \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} + L [\boldsymbol{C}]^{-1} [\boldsymbol{C}] \left ( \frac{d}{dt} \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} \right ) \tag{11}$

と整理できます。上式第2項の[C]⁻¹[C]は単位行列になりますから，真に計算すべきは[C]⁻¹ d/dt [C]の部分だけです。まじめに書くと煩雑すぎるので，θ = ωtに注意しながら計算した結果だけ示すと次式となります。

$\displaystyle [\boldsymbol{C}]^{-1} \frac{d}{dt} [\boldsymbol{C}] = \begin{bmatrix} 0 & 0 & 0 \\ 0 & 0 & -\omega \\ 0 & \omega & 0 \end{bmatrix} \tag{12}$

おお…何というシンプルな式…❗ 腕力ご興味のある方は手計算で試してみてもいいかもしれません。私はSymPyを使います…💦 導出するんじゃなかったんかい❗─というツッコミが期待される（？）ので，SymPyでの計算の様子を後述します。

d-q座標系での回路方程式へ

(12)式を(9)式に戻すと，次式を得ます。

$\displaystyle \begin{bmatrix} v_{0} \\ v_{d} \\ v_{q} \end{bmatrix} = R \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} + L \begin{bmatrix} 0 & 0 & 0 \\ 0 & 0 & -\omega \\ 0 & \omega & 0 \end{bmatrix} \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} + L \frac{d}{dt} \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} + \begin{bmatrix} v_{s0} \\ v_{sd} \\ v_{sq} \end{bmatrix} \tag{13}$

右辺第1項のRをRと単位行列の積と考えて，第1項と第2項をまとめると，

$\displaystyle \begin{bmatrix} v_{0} \\ v_{d} \\ v_{q} \end{bmatrix} = \begin{bmatrix} R & 0 & 0 \\ 0 & R & -\omega L\\ 0 & \omega L & R \end{bmatrix} \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} + L \frac{d}{dt} \begin{bmatrix} i_{0} \\ i_{d} \\ i_{q} \end{bmatrix} + \begin{bmatrix} v_{s0} \\ v_{sd} \\ v_{sq} \end{bmatrix} \tag{14}$

と書けます。ここで，三相3線では零相電流i₀は零ですので，(14)式の1行目，1列目を省略して，

$\displaystyle \begin{bmatrix} v_{d} \\ v_{q} \end{bmatrix} = \begin{bmatrix} R & -\omega L\\ \omega L & R \end{bmatrix} \begin{bmatrix} i_{d} \\ i_{q} \end{bmatrix} + L \frac{d}{dt} \begin{bmatrix} i_{d} \\ i_{q} \end{bmatrix} + \begin{bmatrix} v_{sd} \\ v_{sq} \end{bmatrix} \tag{15}$

と書いても良いでしょう。さらに別の表現の仕方として行列をもっとまとめれば*3，

$\displaystyle \begin{bmatrix} v_{d} \\ v_{q} \end{bmatrix} = \begin{bmatrix} R + L \frac{d}{dt} & -\omega L \\ \omega L & R + L \frac{d}{dt} \end{bmatrix} \begin{bmatrix} i_{d} \\ i_{q} \end{bmatrix} + \begin{bmatrix} v_{sd} \\ v_{sq} \end{bmatrix} \tag{16}$

とも書けますね。 (15)式や(16)式は，そう！図2(b), (c)を数式として表現している──逆に言えば図2(b), (c)は(15)式や(16)式を回路図として表現したものとなっています！ 以上より，座標変換にて三相等価回路からd軸，q軸等価回路を導出することができました。

[C]⁻¹ d/dt [C]の導出

腕力に自信がないので，SymPyを使うことにします。早速ですが，下記のようなプログラムを書きました。

import sympy as sp

# シンボルの定義
omega, t = sp.symbols('omega t', real = True, positive = True)

# α-β逆変換行列の定義
C_alpha_beta = sp.sqrt(sp.Rational(2, 3)) * sp.Matrix([
    [1 / sp.sqrt(2), 1, 0],
    [1 / sp.sqrt(2), -sp.Rational(1, 2), sp.sqrt(3) / 2],
    [1 / sp.sqrt(2), -sp.Rational(1, 2), -sp.sqrt(3) / 2]
])
print('α-β逆変換:')
print(C_alpha_beta)

# d-q逆変換行列の定義
C_dq = sp.Matrix([
    [1, 0, 0],
    [0, sp.cos(omega * t), -sp.sin(omega * t)],
    [0, sp.sin(omega * t), sp.cos(omega * t)]
])
print('d-q逆変換:')
print(C_dq)

# 両逆変換の積
C = sp.simplify(sp.expand_trig(C_alpha_beta * C_dq))
print('両逆変換の積C:')
print(C)

# 干渉項を導くC**(-1) * d/dt C
C_coupling = sp.simplify(sp.expand_trig(C**(-1) * sp.diff(C, t)))
print('干渉項:')
print(C_coupling)

これを実行すると，

α-β逆変換:
Matrix([[sqrt(3)/3, sqrt(6)/3, 0], [sqrt(3)/3, -sqrt(6)/6, sqrt(2)/2], [sqrt(3)/3, -sqrt(6)/6, -sqrt(2)/2]])
d-q逆変換:
Matrix([[1, 0, 0], [0, cos(omega*t), -sin(omega*t)], [0, sin(omega*t), cos(omega*t)]])
両逆変換の積C:
Matrix([[sqrt(3)/3, sqrt(6)*cos(omega*t)/3, -sqrt(6)*sin(omega*t)/3], [sqrt(3)/3, -sqrt(6)*cos(omega*t + pi/3)/3, sqrt(6)*sin(omega*t + pi/3)/3], [sqrt(3)/3, -sqrt(6)*sin(omega*t + pi/6)/3, -sqrt(6)*cos(omega*t + pi/6)/3]])
干渉項:
Matrix([[0, 0, 0], [0, 0, -omega], [0, omega, 0]])

が得られます。最後のMatrix([[0, 0, 0], [0, 0, -omega], [0, omega, 0]])の部分がまさに(12)式に相当します。

ちなみに，分数の「1 / 2」はsp.Raitonal(1, 2)と書かないと，実数の0.5に展開されてしまうようです。今回は罠にはまってしばらく悩んでしまいました…。

次記事の予告

次記事では，三相PWM整流器の直流電圧制御系について，私なりに小信号モデルを立て，安定性を検討した結果について書きたいと思っています。

とは言っても，d-q座標系での電流制御についても，まだまだ私として理解・咀嚼できていないことが多いです。例えば，d-q非干渉制御をしていたとしても，i_dとi_qもサンプリングされたものであれば， −ωLi_qやωLi_dを完全に打ち消すことはできないのではないか？また，d-q座標系では，正相と逆相に対して異なる応答をするはずですし，d-q座標系と三相座標系では周波数がωだけシフトしているはずですが，制御理論としてはどのように扱うべきなのか…私の中で整理ができていないことがあります（世の中では公知なのだと思いますが…💦）。これらについても，いずれかの時点でアタックしてみようと思っています😅

*1:それぞれ「逆0-α-β変換」，「逆0-d-q変換」という言い方が本来は正しいですが，私個人としては目くじらを立てることはないかと思っています。

*2:釈迦に説法とは存じますが，高校数学を思い出せば{x(t) y(t)}' = x'(t) y(t) + x(t) y'(t)ですよね。

*3:数学的に若干厳密さを欠いている気がしますが，微分演算子d/dtが単体で行列に入ることを許容しています。なお，前記事ではp = d/dtとおいていました。