0

私はこの1を使用して 2 次微分方程式を解こうとしましたが、それを正しく行うことができず、オンラインで役立つものは何も見つかりませんでしたが、進歩したと思います。

私はdsolveを使用しました。

syms x(t) v(t) fi(t)

[x(t), v(t)] = dsolve(diff(x) == v, diff(v) == fi/m, x(0) == [-L, -L], v(0) == [5, 10] )

それは私に与えます;

   x(t) =



int(fi(x)/5, x, 0, t, 'IgnoreSpecialCases', true, 'IgnoreAnalyticConstraints', true) - 5





v(t) =

C2 + t*(int(fi(x)/5, x, 0, t, 'IgnoreSpecialCases', true, 'IgnoreAnalyticConstraints', true) - 5) + int(-(x*fi(x))/5, x, 0, t, 'IgnoreSpecialCases', true, 'IgnoreAnalyticConstraints', true)

今、結果を解釈するのに助けが必要で、この結果を使用して ode45 から何かを取得できるかどうか疑問に思っています。また、力場を移動する 500 個の粒子のシミュレーションの参照軌道として、解をプロットしたいと考えています。

ODE45 を使用:

function dxdt = solution(t,y0)
frprintf('Second stop')
.....
dxdt = [x, v]
end

メインファイルから呼び出します:

t:dt:t_f
y0 = [x0,v0]
fprintf('first stop')
[x, v] = ode45(@solution, y0, t)

コードがスムーズに実行される場合、「最初のストップ、2番目のストップ、3番目のストップ」と「4番目のストップ」が出力されるように設定しました。最初のストップのみが出力され、エラーが発生します。

4

1 に答える 1

1

あなたがどこかで定義されていると仮定します

function F = fi(t)
    F = ...
end

ODE関数を（mグローバル変数として使用して）として定義します

function doty = odefunc(t,y)
    doty = [ y(2); fi(t)/m ]
end

そして電話する

t = t0:dt:tf
y0 = [ x0, v0 ]
t,y = ode45(odefunc, t, y0)

plot( t, y(:,1) )

一般に、ベクトルyには、システム内のすべての粒子またはオブジェクトの位相空間 (位置および速度/インパルス) 内のポイントが含まれます。次に、そのodefunc時点でのこの特定の位相空間ポイントの力を計算し、それからt微分ベクトルを位相空間ポイントに構成します。

たとえば、3D シミュレーションでy(6*(k-1)+1:6*(k-1)+3)は、粒子の位置kとy(6*(k-1)+4:6*(k-1)+6)速度ベクトルを調整できます。または、粒子の番目のy(3*(k-1)+1:3*(k-1)+3)位置と速度で位置と速度を分離することもできます。y(3*(N+k-1)+1:3*(N+k-1)+3)kN

于 2015-05-18T20:22:58.400 に答える