このページの内容は最新ではありません。最新版の英語を参照するには、ここをクリックします。

トランジスタのスティッフな微分代数方程式の求解

この例では ode23t を使用して電気回路を記述するスティッフな微分代数方程式 (DAE) を解く方法を説明します [1]。例のファイル amp1dae.m にコード記述されている 1 石トランジスタ アンプの問題は、半陽的な形式に書き換えることもできますが、この例では元の形式 $$M{u}^{\prime }=\varphi (u)$$ で解きます。この問題には、定数の特異質量行列 $\mathit{M}$ が含まれます。

このトランジスタ アンプ回路には 6 つの抵抗、3 つのコンデンサ、1 つのトランジスタが含まれます。

• 初期の電圧信号は $$U_e(t)=0.4\sin \left(200\pi t\right)$$ です。

• 動作電圧は ${\mathit{U}}_{\mathit{b}}=6$ です。

• ノードの電圧は ${\mathit{U}}_{\mathit{i}}\left(\mathit{t}\right)\left(\mathit{i}=1,2,3,4,5\right)$ によって与えられます。

• 抵抗器の値 ${\mathit{R}}_{\mathit{i}}\left(\mathit{i}=1,2,3,4,5,6\right)$ は定数であり、各抵抗器を通る電流は $\mathit{I}=\mathit{U}/\mathit{R}$ を満たします。

• コンデンサの値 ${\mathit{C}}_{\mathit{i}}\left(\mathit{i}=1,2,3\right)$ は定数であり、各コンデンサを通る電流は $\mathit{I}=\mathit{C}\cdot \mathrm{dU}/\mathrm{dt}$ を満たします。

目的はノード 5 にかかる出力電圧 ${\mathit{U}}_5\left(\mathit{t}\right)$ を求めることです。

MATLAB® でこの方程式を解くには、方程式をコード化し、質量行列をコード化し、初期条件と積分区間を設定してからソルバー ode23t を呼び出します。(ここで行ったように) 必要な関数をファイルの最後にローカル関数として含めることも、あるいは個別の名前付きファイルとして MATLAB パスのディレクトリに保存することもできます。

質量行列のコード化

キルヒホッフの法則を使用して各ノード (1 ～ 5) を通る電流を等しくすると、回路を記述する 5 つの方程式からなる系が得られます。

$$\begin{array}{l}\mathrm{node}1:\frac{{\mathit{U}}_{\mathit{e}}\left(\mathit{t}\right)}{{\mathit{R}}_0}-\frac{{\mathit{U}}_1}{{\mathit{R}}_0}+{\mathit{C}}_1\left({{\mathit{U}}_2}^{\prime }-{{\mathit{U}}_1}^{\prime }\right)=0,\\ \mathrm{node}2:\frac{{\mathit{U}}_{\mathit{b}}}{{\mathit{R}}_2}-{\mathit{U}}_2\left(\frac{1}{{\mathit{R}}_1}+\frac{1}{{\mathit{R}}_2}\right)+{\mathit{C}}_1\left({{\mathit{U}}_1}^{\prime }-{{\mathit{U}}_2}^{\prime }\right)-0.01\mathit{f}\left({\mathit{U}}_2-{\mathit{U}}_3\right)=0,\\ \mathrm{node}3:\mathit{f}\left({\mathit{U}}_2-{\mathit{U}}_3\right)-\frac{{\mathit{U}}_3}{{\mathit{R}}_3}-{\mathit{C}}_2{{\mathit{U}}_3}^{\prime }=0,\\ \mathrm{node}4:\frac{{\mathit{U}}_{\mathit{b}}}{{\mathit{R}}_4}-\frac{{\mathit{U}}_4}{{\mathit{R}}_4}+{\mathit{C}}_3\left({{\mathit{U}}_5}^{\prime }-{{\mathit{U}}_4}^{\prime }\right)-0.99\mathit{f}\left({\mathit{U}}_2-{\mathit{U}}_3\right)=0,\\ \mathrm{node}5:-\frac{{\mathit{U}}_5}{{\mathit{R}}_5}+{\mathit{C}}_3\left({{\mathit{U}}_4}^{\prime }-{{\mathit{U}}_5}^{\prime }\right)=0.\end{array}$$

この方程式系の質量行列は、方程式の左辺の微分項を集めることによって求められ、次の形式になります。

$$M=\left(\begin{array}{ccccc}-c_1& c_1& 0& 0& 0\\ c_1& -c_1& 0& 0& 0\\ 0& 0& -c_2& 0& 0\\ 0& 0& 0& -c_3& c_3\\ 0& 0& 0& c_3& -c_3\end{array}\right),$$

ここで、$$k=1,2,3$$ に対して $$c_k=k\times 10^{-6}$$ です。

適切な定数 ${\mathit{c}}_{\mathit{k}}$ を使って質量行列を作成し、関数 odeset を使用して質量行列を指定します。質量行列が特異行列であることは明らかですが、ソルバーによる DAE 問題の自動検出をテストするために、'MassSingular' オプションを既定値 'maybe' のままにします。

c = 1e-6 * (1:3);
M = zeros(5,5);
M(1,1) = -c(1);
M(1,2) =  c(1);
M(2,1) =  c(1);
M(2,2) = -c(1);
M(3,3) = -c(2);
M(4,4) = -c(3);
M(4,5) =  c(3);
M(5,4) =  c(3);
M(5,5) = -c(3);
opts = odeset('Mass',M);

方程式のコード作成

関数 transampdae にはこの例の方程式系が含まれています。この関数は、すべての電圧および定数パラメーターの値を定義します。方程式の左辺に集められた導関数は質量行列にコード化され、transampdae は方程式の右辺をコード化します。

function dudt = transampdae(t,u)
% Define voltages and parameters
Ue = @(t) 0.4*sin(200*pi*t);
Ub = 6;
R0 = 1000;
R15 = 9000;
alpha = 0.99;
beta = 1e-6;
Uf = 0.026;

% Define system of equations
f23 = beta*(exp((u(2) - u(3))/Uf) - 1);
dudt = [ -(Ue(t) - u(1))/R0
    -(Ub/R15 - u(2)*2/R15 - (1-alpha)*f23)
    -(f23 - u(3)/R15)
    -((Ub - u(4))/R15 - alpha*f23)
    (u(5)/R15) ];
end

メモ: この関数は、ローカル関数として例の最後に記載されています。

初期条件のコード化

初期条件を設定します。この例では、[1] で計算された各ノードを通る電流に対し、整合性のある初期条件を使用します。

Ub = 6;
u0(1) = 0;
u0(2) = Ub/2;
u0(3) = Ub/2;
u0(4) = Ub;
u0(5) = 0;

方程式を解く

ode23t を使用して、時間区間 [0 0.05] での DAE 系の解を求めます。

tspan = [0 0.05];
[t,u] = ode23t(@transampdae,tspan,u0,opts);

結果のプロット

初期電圧 $$U_e(t)$$ と出力電圧 $$U_5(t)$$ をプロットします。

Ue = @(t) 0.4*sin(200*pi*t);
plot(t,Ue(t),'o',t,u(:,5),'.')
axis([0 0.05 -3 2]);
legend('Input Voltage U_e(t)','Output Voltage U_5(t)','Location','NorthWest');
title('One Transistor Amplifier DAE Problem Solved by ODE23T');
xlabel('t');

参照

[1] Hairer, E., and Gerhard Wanner. Solving Ordinary Differential Equations II: Stiff and Differential-Algebraic Problems. Springer Berlin Heidelberg, 1991, p. 377.

ローカル関数

ここでは、ODE ソルバー ode23t が解を計算するために呼び出すローカル補助関数を紹介しています。あるいは、この関数を独自のファイルとして MATLAB パスのディレクトリに保存することもできます。

function dudt = transampdae(t,u)
% Define voltages and parameters
Ue = @(t) 0.4*sin(200*pi*t);
Ub = 6;
R0 = 1000;
R15 = 9000;
alpha = 0.99;
beta = 1e-6;
Uf = 0.026;

% Define system of equations
f23 = beta*(exp((u(2) - u(3))/Uf) - 1);
dudt = [ -(Ue(t) - u(1))/R0
    -(Ub/R15 - u(2)*2/R15 - (1-alpha)*f23)
    -(f23 - u(3)/R15)
    -((Ub - u(4))/R15 - alpha*f23)
    (u(5)/R15) ];
end

参考

ode23t | ode15s

関連するトピック

• ODE ソルバーの選択
• ODE オプションの概要
• 微分代数方程式 (DAE) の求解
• 状態依存の強い質量行列による ODE の求解