本技術ノートは、電力系統の潮流計算における主要なアルゴリズムの数学的基礎、実装方法、収束特性、適用場面を包括的に解説します。理論的厳密性と実践的有用性の両立を目指し、学習者から実務者まで幅広いニーズに対応します。
電力系統の定常状態解析において、各母線 $i$ における電力保存則(Kirchhoffの電流則)から出発する。
Step 1: 複素電力の定義
母線 $i$ への複素電力注入は: $S_i = P_i + jQ_i = V_i \cdot I_i^*$
ここで、$V_i$ は母線 $i$ の複素電圧、$I_i^*$ は注入電流の共役である。
Step 2: アドミタンス行列による電流表現
電流は電圧とアドミタンス行列の関係により: $I_i = \sum_{j=1}^n Y_{ij} V_j$
したがって、複素電力は: $S_i = V_i \left(\sum_{j=1}^n Y_{ij} V_j\right)^* = V_i^* \sum_{j=1}^n Y_{ij}^* V_j^*$
Step 3: 実部・虚部分離
| $V_i = | V_i | e^{j\theta_i}$、$Y_{ij} = G_{ij} + jB_{ij}$ として: |
$P_i = \text{Re}[S_i] = \text{Re}\left[V_i^* \sum_{j=1}^n Y_{ij} V_j\right] = P_{Gi} - P_{Di}$ … (1)
$Q_i = \text{Im}[S_i] = \text{Im}\left[V_i^* \sum_{j=1}^n Y_{ij} V_j\right] = Q_{Gi} - Q_{Di}$ … (2)
| 複素電圧を極座標で表現:$V_i = | V_i | e^{j\theta_i}$ |
導出過程:
| $P_i = \sum_{j=1}^n | V_i | V_j | Y_{ij} | \cos(\theta_{ij} - \theta_i + \theta_j)$ … (3) |
| $= \sum_{j=1}^n | V_i | V_j | [G_{ij}\cos(\theta_i - \theta_j) + B_{ij}\sin(\theta_i - \theta_j)]$ … (4) |
| $Q_i = \sum_{j=1}^n | V_i | V_j | Y_{ij} | \sin(\theta_{ij} - \theta_i + \theta_j)$ … (5) |
| $= \sum_{j=1}^n | V_i | V_j | [G_{ij}\sin(\theta_i - \theta_j) - B_{ij}\cos(\theta_i - \theta_j)]$ … (6) |
| ここで、$\theta_{ij} = \arg(Y_{ij})$、$ | Y_{ij} | = \sqrt{G_{ij}^2 + B_{ij}^2}$ である。 |
送電線要素:
送電線を $R + jX$ で表現すると、送電線アドミタンスは:
$y_{ij} = \frac{1}{R_{ij} + jX_{ij}} = \frac{R_{ij}}{R_{ij}^2 + X_{ij}^2} - j\frac{X_{ij}}{R_{ij}^2 + X_{ij}^2} = g_{ij} + jb_{ij}$ … (7)
母線アドミタンス行列:
$Y_{ii} = \sum_{k \in \mathcal{N}i} y{ik} + y_{sh,i}$ (自己アドミタンス) … (8)
$Y_{ij} = -y_{ij}$ (相互アドミタンス, $i \neq j$) … (9)
ここで、$\mathcal{N}i$ は母線 $i$ に接続する母線の集合、$y{sh,i}$ は母線 $i$ の分路アドミタンスである。
微小変動 $\Delta\mathbf{x}$ に対する電力変化は:
| $[\Delta\mathbf{P}; \Delta\mathbf{Q}] = [\frac{\partial \mathbf{P}}{\partial \boldsymbol{\theta}}, \frac{\partial \mathbf{P}}{\partial | \mathbf{V} | }; \frac{\partial \mathbf{Q}}{\partial \boldsymbol{\theta}}, \frac{\partial \mathbf{Q}}{\partial | \mathbf{V} | }] [\Delta\boldsymbol{\theta}; \Delta | \mathbf{V} | ] = \mathbf{J} \Delta\mathbf{x}$ … (10) |
これが潮流計算の線形化方程式の基礎となる。
| 与えられる量: | V | , θ (通常 θ = 0°) |
| 与えられる量: P, | V |
| 求める量: | V | , θ |
P_i = ∑_{j=1}^n |V_i||V_j||Y_{ij}|cos(θ_{ij} - θ_i + θ_j)
Q_i = ∑_{j=1}^n |V_i||V_j||Y_{ij}|sin(θ_{ij} - θ_i + θ_j)
簡略表記:
P_i = ∑_{j=1}^n |V_i||V_j|[G_{ij}cos(θ_i - θ_j) + B_{ij}sin(θ_i - θ_j)]
Q_i = ∑_{j=1}^n |V_i||V_j|[G_{ij}sin(θ_i - θ_j) - B_{ij}cos(θ_i - θ_j)]
Newton-Raphson法は多変数非線形方程式系 $\mathbf{f}(\mathbf{x}) = \mathbf{0}$ の数値解法である。
Step 1: Taylor展開
点 $\mathbf{x}^{(k)}$ における1次Taylor展開:
$\mathbf{f}(\mathbf{x}^{(k)} + \Delta\mathbf{x}^{(k)}) \approx \mathbf{f}(\mathbf{x}^{(k)}) + \mathbf{J}(\mathbf{x}^{(k)}) \Delta\mathbf{x}^{(k)}$ … (11)
ここで、$\mathbf{J}(\mathbf{x})$ はJacobian行列:
$J_{ij}(\mathbf{x}) = \frac{\partial f_i(\mathbf{x})}{\partial x_j}$ … (12)
Step 2: Newton方向の決定
$\mathbf{f}(\mathbf{x}^{(k)} + \Delta\mathbf{x}^{(k)}) = \mathbf{0}$ を要求すると:
$\mathbf{f}(\mathbf{x}^{(k)}) + \mathbf{J}(\mathbf{x}^{(k)}) \Delta\mathbf{x}^{(k)} = \mathbf{0}$ … (13)
したがって、Newton方向は:
$\Delta\mathbf{x}^{(k)} = -[\mathbf{J}(\mathbf{x}^{(k)})]^{-1} \mathbf{f}(\mathbf{x}^{(k)})$ … (14)
Step 3: 反復更新
$\mathbf{x}^{(k+1)} = \mathbf{x}^{(k)} + \Delta\mathbf{x}^{(k)} = \mathbf{x}^{(k)} - [\mathbf{J}(\mathbf{x}^{(k)})]^{-1} \mathbf{f}(\mathbf{x}^{(k)})$ … (15)
変数ベクトル:
| $\mathbf{x} = [\boldsymbol{\theta}_{PV,PQ}; | \mathbf{V} | _{PQ}] = [\theta_2, \ldots, \theta_n; | V | _1, \ldots, | V | _m]^T$ … (16) |
ここで、$n$ は総母線数、$m$ はPQ母線数である。
ミスマッチベクトル:
$\mathbf{f}(\mathbf{x}) = [\Delta\mathbf{P}; \Delta\mathbf{Q}] = [\mathbf{P}{spec} - \mathbf{P}{calc}(\mathbf{x}); \mathbf{Q}{spec} - \mathbf{Q}{calc}(\mathbf{x})]$ … (17)
Newton更新式:
| $[\Delta\boldsymbol{\theta}^{(k)}; \Delta | \mathbf{V} | ^{(k)}] = -[\mathbf{J}{P\theta}, \mathbf{J}{P | V | }; \mathbf{J}{Q\theta}, \mathbf{J}{Q | V | }]^{-1} [\Delta\mathbf{P}^{(k)}; \Delta\mathbf{Q}^{(k)}]$ … (18) |
Jacobian行列は4つのブロックから構成される:
| $\mathbf{J} = [\frac{\partial \mathbf{P}}{\partial \boldsymbol{\theta}}, \frac{\partial \mathbf{P}}{\partial | \mathbf{V} | }; \frac{\partial \mathbf{Q}}{\partial \boldsymbol{\theta}}, \frac{\partial \mathbf{Q}}{\partial | \mathbf{V} | }] = [\mathbf{J}{P\theta}, \mathbf{J}{P | V | }; \mathbf{J}{Q\theta}, \mathbf{J}{Q | V | }]$ … (19) |
| 式(4)より:$P_i = \sum_{j=1}^n | V_i | V_j | [G_{ij}\cos(\theta_i - \theta_j) + B_{ij}\sin(\theta_i - \theta_j)]$ |
対角要素($i = j$):
| $\frac{\partial P_i}{\partial \theta_i} = \sum_{j=1, j \neq i}^n | V_i | V_j | [-G_{ij}\sin(\theta_i - \theta_j) + B_{ij}\cos(\theta_i - \theta_j)]$ … (20) |
式(6)と比較すると:
| $\frac{\partial P_i}{\partial \theta_i} = -Q_i - | V_i | ^2 B_{ii}$ … (21) |
非対角要素($i \neq j$):
| $\frac{\partial P_i}{\partial \theta_j} = | V_i | V_j | [G_{ij}\sin(\theta_i - \theta_j) - B_{ij}\cos(\theta_i - \theta_j)]$ … (22) |
対角要素:
| $\frac{\partial P_i}{\partial | V_i | } = \frac{1}{ | V_i | } \sum_{j=1}^n | V_j | [G_{ij}\cos(\theta_i - \theta_j) + B_{ij}\sin(\theta_i - \theta_j)]$ … (23) |
式(4)を用いると:
| $\frac{\partial P_i}{\partial | V_i | } = \frac{P_i + | V_i | ^2 G_{ii}}{ | V_i | }$ … (24) |
非対角要素:
| $\frac{\partial P_i}{\partial | V_j | } = | V_i | [G_{ij}\cos(\theta_i - \theta_j) + B_{ij}\sin(\theta_i - \theta_j)]$ … (25) |
同様の手順で無効電力についても導出:
$\mathbf{J}_{Q\theta}$ ブロック:
| $\frac{\partial Q_i}{\partial \theta_i} = P_i - | V_i | ^2 G_{ii}$ … (26) |
| $\frac{\partial Q_i}{\partial \theta_j} = - | V_i | V_j | [G_{ij}\cos(\theta_i - \theta_j) + B_{ij}\sin(\theta_i - \theta_j)]$ … (27) |
| **$\mathbf{J}_{Q | V | }$ ブロック:** |
| $\frac{\partial Q_i}{\partial | V_i | } = \frac{Q_i - | V_i | ^2 B_{ii}}{ | V_i | }$ … (28) |
| $\frac{\partial Q_i}{\partial | V_j | } = | V_i | [G_{ij}\sin(\theta_i - \theta_j) - B_{ij}\cos(\theta_i - \theta_j)]$ … (29) |
Jacobian行列の構築:$O(n^2)$
LU分解:$O(n^3)$
前進・後進代入:$O(n^2)$
総計算量:$O(n^3)$ per iteration
1. 初期化: θ_i^(0) = 0, |V_i|^(0) = 1.0 (PQ母線)
2. 反復 k = 0, 1, 2, ...
a) ミスマッチ計算: f^(k) = [ΔP^(k); ΔQ^(k)]
b) 収束判定: ||f^(k)||_∞ < ε なら終了
c) Jacobian構築: J^(k)
d) 線形方程式求解: J^(k) Δx^(k) = -f^(k)
e) 変数更新: x^(k+1) = x^(k) + Δx^(k)
function newtonRaphsonPowerFlow(Ybus, Pbus, Qbus, V0) {
let V = [...V0];
const tolerance = 1e-8;
const maxIter = 20;
for (let iter = 0; iter < maxIter; iter++) {
// 電力ミスマッチ計算
const Scalc = calculatePowerInjection(V, Ybus);
const deltaP = subtractReal(Pbus, Scalc);
const deltaQ = subtractImag(Qbus, Scalc);
// 収束判定
const maxMismatch = Math.max(
Math.max(...deltaP.map(Math.abs)),
Math.max(...deltaQ.map(Math.abs))
);
if (maxMismatch < tolerance) {
return {converged: true, voltage: V, iterations: iter};
}
// Jacobian行列構築
const J = buildJacobian(V, Ybus);
// 線形方程式求解
const mismatch = [...deltaP, ...deltaQ];
const dx = solveLinearSystem(J, mismatch);
// 電圧更新
V = updateVoltage(V, dx);
}
return {converged: false, voltage: V, iterations: maxIter};
}
function buildJacobian(V, Ybus) {
const n = V.length;
const J = new Array(2*n-1).fill(0).map(() => new Array(2*n-1).fill(0));
for (let i = 1; i < n; i++) { // Slack母線除く
for (let j = 1; j < n; j++) {
if (i === j) {
// 対角要素
J[i-1][j-1] = -calculateQi(i, V, Ybus) -
Ybus[i][i].imag * Math.pow(Math.abs(V[i]), 2);
} else {
// 非対角要素
const Vi = Math.abs(V[i]);
const Vj = Math.abs(V[j]);
const theta_ij = Math.arg(V[i]) - Math.arg(V[j]);
const Y_ij = Ybus[i][j];
J[i-1][j-1] = Vi * Vj * (Y_ij.real * Math.sin(theta_ij) -
Y_ij.imag * Math.cos(theta_ij));
}
}
}
return J;
}
高速分離解法は、実用送電系統の物理的特性に基づく近似により、Newton-Raphson法を簡素化した手法である。
送電系統の以下の物理的特性を仮定:
仮定 A1: $r_{ij} \ll x_{ij}$ (送電線抵抗 $\ll$ リアクタンス)
仮定 A2: $|\theta_i - \theta_j| \ll 1$ (母線間位相角差は小さい)
仮定 A3: $|V_i| \approx 1.0$ (電圧大きさは基準値付近)
仮定 A4: $Q_i/|V_i|^2 \ll B_{ii}$ (無効電力項 $\ll$ 自己サセプタンス)
仮定A1, A2により、式(22), (25)は次のように近似できる:
$G_{ij}\sin(\theta_i - \theta_j) \approx G_{ij}(\theta_i - \theta_j) \approx 0$ ($\because G_{ij} \ll B_{ij}$)
$B_{ij}\cos(\theta_i - \theta_j) \approx B_{ij}$
したがって:
| $\frac{\partial P_i}{\partial \theta_j} \approx - | V_i | V_j | B_{ij}$ … (30) |
同様に、仮定A3により:
| $\frac{\partial P_i}{\partial | V_j | } \approx B_{ij}(\theta_i - \theta_j) \approx 0$ … (31) |
近似により、Jacobian行列は以下のように分離される:
$\mathbf{J} \approx [-\mathbf{B}’, \mathbf{0}; \mathbf{0}, -\mathbf{B}’’]$ … (32)
ここで:
| $\mathbf{B}’’$: Q- | V | カップリング行列 |
B’行列(有効電力-位相角):
$B’{ii} = \sum{j \in \mathcal{N}i} \frac{1}{x{ij}}$ … (33)
$B’{ij} = -\frac{1}{x{ij}}$ ($i \neq j$, 隣接母線) … (34)
B’‘行列(無効電力-電圧大きさ):
$B’‘{ii} = B’{ii} + b_{sh,i}$ … (35)
$B’‘{ij} = B’{ij}$ … (36)
ここで、$b_{sh,i}$ は母線 $i$ の分路サセプタンスである。
近似により、連立方程式は2つの独立な部分問題に分離される:
P-θ 部分問題:
| $\Delta\mathbf{P} = -\mathbf{B}’ \Delta\boldsymbol{\theta} \quad \Rightarrow \quad \Delta\boldsymbol{\theta} = -(\mathbf{B}’)^{-1} \frac{\Delta\mathbf{P}}{ | \mathbf{V} | }$ … (37) |
| **Q- | V | 部分問題:** |
| $\Delta\mathbf{Q} = -\mathbf{B}’’ \Delta | \mathbf{V} | \quad \Rightarrow \quad \Delta | \mathbf{V} | = -(\mathbf{B}’’)^{-1} \frac{\Delta\mathbf{Q}}{ | \mathbf{V} | }$ … (38) |
初期化:
反復計算($k = 0, 1, 2, \ldots$):
Step 1: P-θ サブ問題
| $\Delta\mathbf{P}^{(k)} = \mathbf{P}{spec} - \mathbf{P}{calc}(\boldsymbol{\theta}^{(k)}, | \mathbf{V} | ^{(k)})$ … (39) |
| $\Delta\boldsymbol{\theta}^{(k)} = -(\mathbf{L}_1 \mathbf{U}_1)^{-1} \frac{\Delta\mathbf{P}^{(k)}}{ | \mathbf{V} | ^{(k)}}$ … (40) |
$\boldsymbol{\theta}^{(k+1)} = \boldsymbol{\theta}^{(k)} + \Delta\boldsymbol{\theta}^{(k)}$ … (41)
| Step 2: Q- | V | サブ問題 |
| $\Delta\mathbf{Q}^{(k)} = \mathbf{Q}{spec} - \mathbf{Q}{calc}(\boldsymbol{\theta}^{(k+1)}, | \mathbf{V} | ^{(k)})$ … (42) |
| $\Delta | \mathbf{V} | ^{(k)} = -(\mathbf{L}_2 \mathbf{U}_2)^{-1} \frac{\Delta\mathbf{Q}^{(k)}}{ | \mathbf{V} | ^{(k)}}$ … (43) |
| $ | \mathbf{V} | ^{(k+1)} = | \mathbf{V} | ^{(k)} + \Delta | \mathbf{V} | ^{(k)}$ … (44) |
収束判定:
| $\max( | \Delta\mathbf{P}^{(k)} | _\infty, | \Delta\mathbf{Q}^{(k)} | _\infty) < \varepsilon$ … (45) |
B'_{ii} = ∑_{j≠i} X_{ij}^(-1)
B'_{ij} = -X_{ij}^(-1) (i ≠ j)
B''_{ii} = B'_{ii} + B_{shi}
B''_{ij} = B'_{ij}
1. 初期化・前処理
- B', B''行列構築
- LU分解: [L₁U₁] = B', [L₂U₂] = B''
2. 反復計算 k = 0, 1, 2, ...
a) P-θ サブ問題
- ΔP^(k) = P_spec - P_calc^(k)
- Δθ^(k) = (B')^(-1) (ΔP^(k) / |V|^(k))
- θ^(k+1) = θ^(k) + Δθ^(k)
b) Q-|V| サブ問題
- ΔQ^(k) = Q_spec - Q_calc^(k)
- Δ|V|^(k) = (B'')^(-1) (ΔQ^(k) / |V|^(k))
- |V|^(k+1) = |V|^(k) + Δ|V|^(k)
c) 収束判定
function fastDecoupledPowerFlow(Ybus, Pbus, Qbus, V0) {
// B'およびB''行列構築
const Bp = buildBprimeMatrix(Ybus);
const Bpp = buildBdoubleprimeMatrix(Ybus);
// LU分解(一度のみ)
const LU_Bp = luDecomposition(Bp);
const LU_Bpp = luDecomposition(Bpp);
let V = [...V0];
const tolerance = 1e-8;
const maxIter = 30;
for (let iter = 0; iter < maxIter; iter++) {
// P-θ サブ問題
const deltaP = calculateActiveMismatch(V, Pbus, Ybus);
const deltaPoverV = deltaP.map((dp, i) => dp / Math.abs(V[i]));
const deltaTheta = solveLU(LU_Bp, deltaPoverV);
// 位相角更新
for (let i = 0; i < V.length; i++) {
const mag = Math.abs(V[i]);
const newAngle = Math.arg(V[i]) + deltaTheta[i];
V[i] = mag * Math.exp(new Complex(0, newAngle));
}
// Q-|V| サブ問題
const deltaQ = calculateReactiveMismatch(V, Qbus, Ybus);
const deltaQoverV = deltaQ.map((dq, i) => dq / Math.abs(V[i]));
const deltaV = solveLU(LU_Bpp, deltaQoverV);
// 電圧大きさ更新
for (let i = 0; i < V.length; i++) {
const newMag = Math.abs(V[i]) + deltaV[i];
const angle = Math.arg(V[i]);
V[i] = newMag * Math.exp(new Complex(0, angle));
}
// 収束判定
const maxP = Math.max(...deltaP.map(Math.abs));
const maxQ = Math.max(...deltaQ.map(Math.abs));
if (Math.max(maxP, maxQ) < tolerance) {
return {converged: true, voltage: V, iterations: iter};
}
}
return {converged: false, voltage: V, iterations: maxIter};
}
反復法の最も基本的な形:
V_i^(k+1) = (1/Y_{ii}) [S_i*/V_i* - ∑_{j≠i} Y_{ij}V_j^(k+1/k)]
ここで:
S_i* = P_i - jQ_i: 指定電力の共役V_i*: 電圧の共役j < i: 既に更新済みの電圧を使用j > i: 前回反復の電圧を使用1. 初期化: V_i^(0) = 1.0∠0° (PQ母線)
2. 反復 k = 0, 1, 2, ...
for i = 2 to n (Slack母線以外)
if (母線iがPQ母線) then
V_i^(k+1) = (1/Y_{ii})[S_i*/V_i* - ∑_{j≠i} Y_{ij}V_j^(k+1/k)]
else if (母線iがPV母線) then
V_temp = (1/Y_{ii})[S_i*/V_i* - ∑_{j≠i} Y_{ij}V_j^(k+1/k)]
V_i^(k+1) = |V_i|^spec × (V_temp/|V_temp|)
収束判定
V_i^(k+1) = (1-ω)V_i^(k) + ω·V_i^new
ここで、ωは緩和係数(通常 1.0 < ω < 2.0)
function gaussSeidelPowerFlow(Ybus, busData, V0) {
let V = [...V0];
const tolerance = 1e-8;
const maxIter = 100;
const omega = 1.6; // SOR加速係数
for (let iter = 0; iter < maxIter; iter++) {
const V_old = [...V];
for (let i = 1; i < V.length; i++) { // Slack母線除く
const bus = busData[i];
// 指定電力
const Sspec = new Complex(bus.PG - bus.PD, bus.QG - bus.QD);
// 他母線からの影響項計算
let sum = new Complex(0, 0);
for (let j = 0; j < V.length; j++) {
if (j !== i) {
sum = sum.add(Ybus[i][j].multiply(V[j]));
}
}
// 新しい電圧値計算
const Vnew = Sspec.conjugate().divide(V[i].conjugate())
.subtract(sum)
.divide(Ybus[i][i]);
if (bus.type === 'PQ') {
// PQ母線: そのまま更新
V[i] = V[i].multiply(1 - omega).add(Vnew.multiply(omega));
} else if (bus.type === 'PV') {
// PV母線: 電圧大きさ制約
const Vmag_spec = bus.V;
const Vnew_normalized = Vnew.multiply(Vmag_spec / Vnew.magnitude());
V[i] = V[i].multiply(1 - omega).add(Vnew_normalized.multiply(omega));
}
}
// 収束判定
let maxMismatch = 0;
for (let i = 1; i < V.length; i++) {
const mismatch = calculatePowerMismatch(i, V, busData, Ybus);
maxMismatch = Math.max(maxMismatch, Math.abs(mismatch.P), Math.abs(mismatch.Q));
}
if (maxMismatch < tolerance) {
return {converged: true, voltage: V, iterations: iter};
}
}
return {converged: false, voltage: V, iterations: maxIter};
}
| ** | V_i | = 1.0**: すべての母線電圧大きさ = 1.0pu |
AC潮流方程式:
P_i = ∑_{j=1}^n |V_i||V_j|[G_{ij}cos(θ_i - θ_j) + B_{ij}sin(θ_i - θ_j)]
DC近似:
P_i ≈ ∑_{j=1}^n B_{ij}(θ_i - θ_j) = ∑_{j=1}^n B_{ij}θ_i - ∑_{j=1}^n B_{ij}θ_j
行列形式:
P = Bθ
B_{ii} = ∑_{j≠i} (1/X_{ij})
B_{ij} = -1/X_{ij} (i ≠ j)
1. サセプタンス行列B構築
2. 指定有効電力ベクトルP構築
3. 基準母線(Slack)除去
4. 線形方程式求解: θ = B^(-1) P
5. 送電線潮流計算: P_{ij} = (θ_i - θ_j)/X_{ij}
function dcPowerFlow(branchData, busData) {
const n = busData.length;
// サセプタンス行列構築
const B = buildBMatrix(branchData, n);
// 指定有効電力ベクトル
const P = busData.map(bus => bus.PG - bus.PD);
// Slack母線(通常母線1)を除去
const B_reduced = removeSlackBus(B, 0);
const P_reduced = removeSlackBus(P, 0);
// 線形方程式求解
const theta_reduced = solveLinearSystem(B_reduced, P_reduced);
// 完全な位相角ベクトル復元
const theta = [0, ...theta_reduced]; // Slack母線は0
// 送電線潮流計算
const branchFlows = branchData.map(branch => {
const i = branch.fbus - 1;
const j = branch.tbus - 1;
const flow = (theta[i] - theta[j]) / branch.x;
return {
from: branch.fbus,
to: branch.tbus,
flow: flow
};
});
return {
bus: busData.map((bus, i) => ({
id: bus.id,
angle: theta[i] * 180 / Math.PI, // 度に変換
voltage: 1.0 // DC仮定
})),
branch: branchFlows
};
}
function buildBMatrix(branchData, nbus) {
const B = new Array(nbus).fill(0).map(() => new Array(nbus).fill(0));
branchData.forEach(branch => {
const i = branch.fbus - 1;
const j = branch.tbus - 1;
const b = 1.0 / branch.x;
B[i][i] += b;
B[j][j] += b;
B[i][j] -= b;
B[j][i] -= b;
});
return B;
}
| 電圧情報: | V | , Q の情報なし |
Newton法の収束性改善:
x^(k+1) = x^(k) - [J^T J + μI]^(-1) J^T f(x^(k))
μ: 制御パラメータ(大→Gauss法、小→Newton法)電圧安定性解析のための手法:
F(x,λ) = 0
x: 状態変数(電圧)λ: 負荷パラメータ非反復的解法:
V(s) = ∑_{n=0}^∞ V_n s^n
s: 埋め込みパラメータ定理 (Newton-Kantorovich定理)
$\mathbf{f}: \mathbb{R}^n \to \mathbb{R}^n$ を $C^2$ 級関数とし、$\mathbf{x}^$ を $\mathbf{f}(\mathbf{x}^) = \mathbf{0}$ の解とする。以下の条件が満たされるとき、$\mathbf{x}^*$ の近傍で開始されたNewton反復は二次収束する:
| $\mathbf{J}(\mathbf{x})$ がLipschitz連続:$ | \mathbf{J}(\mathbf{x}) - \mathbf{J}(\mathbf{y}) | \leq L | \mathbf{x} - \mathbf{y} | $ |
収束率:
| $ | \mathbf{x}^{(k+1)} - \mathbf{x}^* | \leq C | \mathbf{x}^{(k)} - \mathbf{x}^* | ^2$ … (70) |
ここで、$C$ は定数である。
Newton法の収束領域を拡大する手法:
信頼領域法 (Trust Region Method)
各反復で以下の制約付き最適化問題を解く:
| $\min_{\mathbf{d}} | \mathbf{f}(\mathbf{x}^{(k)}) + \mathbf{J}(\mathbf{x}^{(k)})\mathbf{d} | ^2 \quad \text{s.t.} \quad | \mathbf{d} | \leq \Delta^{(k)}$ … (71) |
線探索法 (Line Search Method)
ステップサイズ $\alpha$ を調整:
$\mathbf{x}^{(k+1)} = \mathbf{x}^{(k)} + \alpha^{(k)} \mathbf{d}^{(k)}$ … (72)
ここで、$\mathbf{d}^{(k)} = -[\mathbf{J}(\mathbf{x}^{(k)})]^{-1}\mathbf{f}(\mathbf{x}^{(k)})$ はNewton方向。
Gauss-Seidel法は不動点反復 $\mathbf{x} = \mathbf{G}(\mathbf{x})$ として定式化:
$V_i^{(k+1)} = G_i(V_1^{(k+1)}, \ldots, V_{i-1}^{(k+1)}, V_{i+1}^{(k)}, \ldots, V_n^{(k)})$ … (73)
収束条件 (Banachの不動点定理)
$\mathbf{G}: D \subset \mathbb{C}^n \to D$ が縮小写像:
| $ | \mathbf{G}(\mathbf{x}) - \mathbf{G}(\mathbf{y}) | \leq \rho | \mathbf{x} - \mathbf{y} | , \quad 0 \leq \rho < 1$ … (74) |
このとき、唯一の不動点 $\mathbf{x}^*$ が存在し、任意の初期値から収束:
| $ | \mathbf{x}^{(k)} - \mathbf{x}^* | \leq \rho^k | \mathbf{x}^{(0)} - \mathbf{x}^* | $ … (75) |
線形収束率:
| $\rho = \lim_{k \to \infty} \frac{ | \mathbf{x}^{(k+1)} - \mathbf{x}^* | }{ | \mathbf{x}^{(k)} - \mathbf{x}^* | }$ … (76) |
Gauss-Seidel行列を $\mathbf{M} = (\mathbf{D} + \mathbf{L})^{-1}\mathbf{U}$ とすると:
収束条件: $\rho(\mathbf{M}) < 1$ (スペクトル半径 < 1)
最適SOR係数:
$\omega_{opt} = \frac{2}{1 + \sqrt{1 - \rho(\mathbf{B})^2}}$ … (77)
ここで、$\mathbf{B}$ はJacobi反復行列。
$\mathbf{B}’$, $\mathbf{B}’’$ 行列の条件数による収束性評価:
| $\kappa(\mathbf{B}’) = | \mathbf{B}’ | \cdot | (\mathbf{B}’)^{-1} | $, $\kappa(\mathbf{B}’’) = | \mathbf{B}’’ | \cdot | (\mathbf{B}’’)^{-1} | $ … (78) |
収束性能:
真のJacobian $\mathbf{J}$ と近似Jacobian $\mathbf{J}_{FD}$ の差:
| $\Delta\mathbf{J} = \mathbf{J} - \mathbf{J}{FD} = [\mathbf{0}, \mathbf{J}{P | V | }; \mathbf{J}_{Q\theta}, \mathbf{0}]$ … (79) |
摂動理論により、収束次数は以下に低下:
| $r_{FD} \approx 1 + \frac{\log( | \Delta\mathbf{J} | )}{\log( | \mathbf{f} | )}$ … (80) |
最大値ノルム:
| $ | \mathbf{f}(\mathbf{x}^{(k)}) | _\infty = \max_i | f_i(\mathbf{x}^{(k)}) | < \varepsilon$ … (81) |
ユークリッドノルム:
| $ | \mathbf{f}(\mathbf{x}^{(k)}) | _2 = \sqrt{\sum_i f_i(\mathbf{x}^{(k)})^2} < \varepsilon$ … (82) |
| $\frac{ | \mathbf{x}^{(k+1)} - \mathbf{x}^{(k)} | }{ | \mathbf{x}^{(k)} | } < \varepsilon_{rel}$ … (83) |
工学的実用性を考慮:
| $ | \Delta P_i | < \varepsilon_P$ かつ $ | \Delta Q_i | < \varepsilon_Q$ $\forall i$ … (84) |
典型値:$\varepsilon_P = \varepsilon_Q = 10^{-4}$ pu (工学精度)、$10^{-8}$ pu (研究精度)
機械精度 $\varepsilon_m$ の影響:
| $\delta\mathbf{x} \approx \kappa(\mathbf{J}) \varepsilon_m | \mathbf{x} | $ … (85) |
実用判定基準:
部分ピボット選択:
| $ | a_{kk}^{(k)} | = \max_{k \leq i \leq n} | a_{ik}^{(k)} | $ … (86) |
完全ピボット選択:
| $ | a_{pq}^{(k)} | = \max_{\substack{k \leq i \leq n \ k \leq j \leq n}} | a_{ij}^{(k)} | $ … (87) |
x^(k+1) = ∑_{i=0}^m θ_i G(x^(k-i))
x_accelerated = x_k - (x_{k+1} - x_k)² / (x_{k+2} - 2x_{k+1} + x_k)
PA = LU
最適化:
function [V, converged] = newtonpf(Ybus, Sbus, V0, tol)
% 高水準・簡潔な実装
% 豊富な数値ライブラリ
% プロトタイピングに最適
def newton_raphson_pf(Ybus, Sbus, V0, tol=1e-8):
# NumPy/SciPy活用
# 可読性と拡張性
# 機械学習との親和性
bool newton_raphson_pf(SparseMatrix& Ybus, Vector& Sbus,
Vector& V, double tol) {
// 高速・メモリ効率
// 商用ソフトウェア品質
// 大規模系統対応
}
function newtonRaphsonPF(Ybus, Sbus, V0, tolerance) {
// ブラウザ実行
// 教育・可視化
// インタラクティブ性
}
| 手法 | 単位反復 | 総反復数 | 総計算量 |
|---|---|---|---|
| Newton-Raphson | O(n³) | 3-6 | O(n³) |
| Fast Decoupled | O(n³/2) | 5-12 | O(n³) |
| Gauss-Seidel | O(n²) | 20-50 | O(n³) |
| DC Power Flow | O(n³) | 1 | O(n³) |
| 手法 | 主要行列 | スパース率 | メモリ量 |
|---|---|---|---|
| Newton-Raphson | Jacobian | 90-95% | O(n²) |
| Fast Decoupled | B’, B’’ | 95-98% | O(n²/2) |
| Gauss-Seidel | なし | - | O(n²) |
| DC Power Flow | B | 95-98% | O(n²/4) |
手法 反復数 時間[ms] メモリ[MB] 精度
Newton-Raphson 4 12.3 2.1 1e-10
Fast Decoupled 7 8.9 1.6 1e-8
Gauss-Seidel 23 15.7 1.2 1e-8
DC Power Flow 1 1.2 0.8 1e-2
| 母線数 | 反復数 | 時間[s] | メモリ[MB] | 備考 |
|---|---|---|---|---|
| 14 | 4 | 0.01 | 2.1 | 小規模 |
| 57 | 5 | 0.15 | 8.9 | 中規模 |
| 118 | 6 | 0.45 | 28.4 | 大規模 |
| 300 | 7 | 2.1 | 156 | 実用規模 |
| 2383 | 8 | 45 | 1024 | 実系統規模 |
適用場面:
避けるべき場面:
適用場面:
避けるべき場面:
適用場面:
避けるべき場面:
適用場面:
避けるべき場面:
| 記号 | 意味 |
|---|---|
| V_i | 母線iの複素電圧 |
| θ_i | 母線iの電圧位相角 |
| P_i, Q_i | 母線iの有効・無効電力注入 |
| Y_{ij} | アドミタンス行列要素 |
| G_{ij}, B_{ij} | コンダクタンス・サセプタンス |
| J | Jacobian行列 |
| ε | 収束判定値 |
詳細な系統データは各HTMLファイルに実装済み。
基本的なアルゴリズム実装のテンプレートコードを提供。
文書: power_flow_methods.md
作成者: Power Flow Visualization Project
作成日: 2024年12月30日
更新日: 2025年1月1日
バージョン: 2.0 HTML対応版
ライセンス: Educational Use