0

2 つの方程式の ODE 系がありますが、1 つの方程式だけを使用して他の方程式の結果を最小化したいと考えています。

1)

t=linspace(0,2,3);
syms x(t) y(t);
inits='x(0)=2,y(0)=0';
[x,y]=dsolve('Dx=y','Dy=(y*2)-x', inits)

x = 2*exp(t) - 2*t*exp(t);
y = -2*t*exp(t)

xx=eval(vectorize(x));

xx = 2.0000; 0; -14.7781

yy=eval(vectorize(y));

yy = 0; -5.4366; -29.5562

結果が得られた後、1 つの方程式だけで解いて、Dy 方程式で xx 配列を使用しようとしました。

2)

inits='y(0)=0';
[y]=dsolve('Dy=(y*2)-xx', inits);

y = xx/2 - (xx*exp(2*t))/2

yy=eval(vectorize(y));

yy = 0; 0; 396.0397

値は、最初の例と同じではありません。配列を使用して同じ結果を得るには?

4

1 に答える 1

1

問題の 1 つは、変数 xx がシンボリックではないため、シンボリック ソルバーがそれを定数と見なしているように見えることです。

より大きな問題は、xx 値が単に 3 点のベクトルである場合に、matlab が xx 値を連続関数として処理する方法を正確に特定していないことです! 2 番目のケースでも出力が同じであると期待しているという事実は、ある種の誤解を示しています。

しかし、これを明確にするために、xx を ZOH (ゼロ次ホールド) 連続信号として扱いたいと仮定しましょう。これを象徴的に処理するには、Heaviside 関数を使用して明示的に ZOH 信号を作成する必要があると思います。

別の方法として、たとえば ode45 を使用して数値的に解決することもできます

t  = [0,1,2];
xx = [2, 0, -14.7781];
dydy = @(t,y) 2*y - xx(1+trunc(t));
y = ode45(dydt, [0,2], 0);

これは、それぞれ[0, 1, 2]の t 値で[0, -6.39, -47.21]の yy 値を返します。

これは、ZOH システムの[0, 1-e^2, e^2-e^4]の理論値 (手動で計算) とよく一致します。

ご覧のとおり、上記の回答は yy = [0, -5.4366, -29.5562]の元のソリューションとより一致しています。ただし、当然のことながら、2 つのシステムは異なります。最初のシステムには連続時間指数信号が供給されたのに対し、2 番目のシステムには非常に粗くサンプリングされた近似が供給されたからです。

より高速なレート (より細かい時間間隔) でサンプリングし、ZOH よりも優れたものでサンプル間ポイントを補間することによって、2 つをより類似させます。

更新: ありがとうございます。ZOH連続信号の作成を手伝ってもらえませんか？どうやってするか？

上記の例では、xx ベクトルで指定された 3 つのポイントを使用し、"xx(1+trunc(t))" を使用してこれらにアクセスすることにより、微分関数 (dydx) で ZOH を作成しました。これは、trunc (切り捨て) を使用して、サンプル間の (非整数) 時間中に入力定数を明示的に保持します。

ODE は線形であるため、matlab 関数 "lsim()" を使用して、時間ベクトルと入力ベクトルを直接指定したり、入力補間のタイプを直接指定したりすることもできます (実際には ZOH を含む)。デフォルト）。

例えば：

t=[0,1,2]
x=[2,0,-2*e^2]
num=-1
den=[1,-2]
mytf = tf(num,den)
y = lsim(mytf,x,t,0,'zoh');

以前の ode45 数値解と同様に、これは次の同じ解を与えます。

y = [0.00000; -6.38906; -47.20909]

更新（再） ありがとうございます。ZOH連続信号の作成を手伝ってもらえませんか？どうやってするか？

シンボリック ソルバーについて説明します。私は Matlab シンボリック ライブラリにアクセスできませんが、本当にシンボリック ソルバーを使用したい場合は、前に説明したように、ヘヴィサイド (単位ステップ) 関数を使用して連続時間 ZOH 信号を作成できます。次のようなものがそれを行うはずです：

syms xzoh(t)
xzoh = xx(1)*heaviside(t) + (xx(2)-xx(1))*heaviside(t-1) + (xx(3)-xx(2))*heaviside(t-2)
于 2013-03-27T14:17:22.293 に答える