0

t=0、すべての初期条件 t=0 y_0=(0,0,0) で一連の ODE (dy/dt) を解いています。さまざまな時点で y の値に数値を追加して (たとえば、t=10 では y1 をその数値に追加し、t=20 では y2 をその数値に追加する必要があります)、方程式を解くことはできますか?

4

2 に答える 2

2

あなたが提案する方法（および@macduffによって示される方法）でODEに大きな不連続を挿入すると、精度が低下し、計算時間が長くなる可能性があります（特にode45-ode15sを使用すると、より良いオプションになるか、少なくとも絶対許容誤差と相対許容誤差が確実になるようにしてください適切）。非常に剛性の高いシステムを効果的に作成しました。特定の時間から始まる ODE に数値を追加する場合は、ソルバーが独自に選択した特定の時点でのみこれらの方程式を評価することに注意してください。(3 つ以上の要素を指定することで、固定ステップ サイズの出力を取得できるという事実に惑わされないでくださいtspan。Matlab のソルバーはすべて可変ステップ サイズ ソルバーであり、エラー基準に基づいて実際のステップを選択します。)

より良いオプションは、システムを区分的に統合し、各実行からの結果の出力を一緒に追加することです。

% t = 0 to t = 10, pass parameter a = 0 to add to ODEs
a = 0;
tspan = [0 10];
[T,Y] = ode45(@(t,y)myfun(t,y,a),tspan,y_0);

% t = 10 to t = 20, pass parameter a = 10 to add to ODEs
a = 10;
[t,y] = ode45(@(t,y)myfun(t,y,a),tspan+T(end),Y(end,:));
T = [T;t(2:end)];
Y = [Y;y(2:end,:)];

% t = 20 to t = 30, pass parameter a = 20 to add to ODEs
a = 20;
[t,y] = ode45(@(t,y)myfun(t,y,a),tspan+T(end),Y(end,:));
T = [T;t(2:end)];
Y = [Y;y(2:end,:)];

TMatlab エディターは、配列が事前に割り当てられてYいないか、または成長していないことについて不平を言うことがありますが、この場合、大きなチャンクで数回しか成長していないため、問題ありません。または、固定の出力ステップ サイズが必要な場合は、次のようにすることができます。

dt = 0.01;
T = 0:dt:30;
Y = zeros(length(T),length(y_0));

% t = 0 to t = 10, pass parameter a = 0 to add to ODEs
a = 0;
[~,Y(1:10/dt+1,:)] = ode45(@(t,y)myfun(t,y,a),T(1:10/dt+1),y_0);

% t = 10 to t = 20, pass parameter a = 10 to add to ODEs
a = 10;
[~,Y(10/dt+1:20/dt+1,:)] = ode45(@(t,y)myfun(t,y,a),T(10/dt+1:20/dt+1),Y(10/dt+1,:));

% t = 20 to t = 30, pass parameter a = 20 to add to ODEs
a = 20;
[~,Y(20/dt+1:end,:)] = ode45(@(t,y)myfun(t,y,a),T(20/dt+1:end),Y(20/dt+1,:));

必要に応じて、上記の両方のコード ブロックをよりコンパクトなforループに簡単に変換できます。

どちらの場合も、ODE 関数は次のようmyfunにパラメーターを組み込みますa。

function ydot = myfun(t,y,a)
    y(1) = ... % add a however you like
    ...
于 2013-06-28T17:48:08.823 に答える