0

我一直在尝试使用这个1解决二阶微分方程，但是我无法正确解决，并且在网上没有发现任何有用的信息，但我相信我已经取得了进展。

我用过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)

我已经设置好了，如果代码运行顺利，它会打印“第一站，第二站，第三站”和“第四站”，它只打印出第一站，这就是我得到错误的地方。

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)for position 和y(3*(N+k-1)+1:3*(N+k-1)+3)for the velocity for the k-th ofparticles来分离位置和速度N。

于 2015-05-18T20:22:58.400 回答