0

私のコードでは、ソリューション ベクトルのどこでもゼロが返されますが、その理由がわかりません。結合された 2 次 ODE を 4 つの 1 次 ODE に分解しました。

関数を xp.m として定義しています

function zprime = f(t,z)
a = 1;
b = 1;
c = 1.5;

zprime = zeros(4,1);

zprime(1) = z(2);
zprime(2) = -a*z(1) + b*(z(3) - z(1));
zprime(3) = z(4);
zprime(4) = -c*(z(3) - z(1));
end

次のコマンドを使用して、matlab で実行します。

[t,z] = ode45('xp',[1,100],[0 0 0 0]);

私の初期条件はすべて0なので、私の初期条件が0の解または何か他のものを与えるのですか? ic を変更すると、予想どおり、ソリューションが変更されます。

ありがとう

4

2 に答える 2

4

IC を使用した特定のケースではz_0 = [0,0,0,0]、解は の値で安定している必要がありますz_out = [0,0,0,0]。プラグインしz = z_0てODEソルバーで実行するときは、関数を見てください。

zprime(1) = z(2);                       % --->  0
zprime(2) = -a*z(1) + b*(z(3) - z(1));  % ---> -a*(0) + b*(0)
zprime(3) = z(4);                       % ---> 0
zprime(4) = -c*(z(3) - z(1));           % ---> -c*(0-0) = 0

数値解の一般的な前提を念頭に置いてください。初期条件を取り、導関数の式に入力します。これは、解いている関数の勾配を示します。その勾配を使用して、将来のある時点 (または近くの場所など) での関数の値を決定し、それを微分式にフィードバックして、最初からやり直します。

異なる方法 (前方/後方オイラー、RK4 など) の唯一の大きな違いは、現在の位置で勾配を決定するために使用する方法です。

于 2013-10-17T20:19:13.030 に答える
2

あなたの実際の質問は、Matlab やプログラミングよりも数学に関係しています。function にゼロをプラグインすると、ゼロf以外の答えが得られないことがすぐにわかります。平衡点または不動点を調べる必要があります。平衡が不安定であっても (丘の頂上を想像してみてください)、外乱 (外部入力、数値エラー) のない正確にその時点の状態は、永久に平衡を保ちます。微分方程式を扱う場合は、平衡点を見つける方法と、これらの点で評価されるヤコビ行列を計算して系の特性を決定する方法を知っておくとよいでしょう。この分野でさらに質問がある場合は、math.stackexchange.comにアクセスすることをお勧めします。

また、統合関数を呼び出すために古いスキームを使用しています。パラメータを渡すこともできます。統合関数をメイン関数のサブ関数として定義するか、別の関数ファイルとして定義します (関数名と同じファイル名 – 以外の名前が必要ですf) 。

function zprime = f(t,z,a,b,c)
zprime(1,1) = z(2);
zprime(2,1) = -a*z(1) + b*(z(3) - z(1));
zprime(3,1) = z(4);
zprime(4,1) = -c*(z(3) - z(1));

次に、メイン関数で呼び出します

a = 1;
b = 1;
c = 1.5;
[t,z] = ode45(@(t,z)f(t,z,a,b,c),[1 100],[0 0 0 0]);
于 2013-10-17T20:22:55.310 に答える