ODE45 to solve a system of two coupled 2nd order ODEs

2 次查看（过去 30 天）
Ricardo Machado
Ricardo Machado 2019-8-31
评论： David Goodmanson 2019-9-2
So this is my code for a system of coupled oscillators
syms x1(t) x2(t) k1 k2 m
Dx1 = diff(x1);
D2x1 = diff(x1,2);
Dx2 = diff(x2);
D2x2 = diff(x2,2);
Eq1 = D2x1 == (-(k1+k2)*x1+(k2)*x2)/m;
Eq2 = D2x2 == ((k2*x1)+((k1+k2)*x2))/m;
[V,Subs] = odeToVectorField(Eq1, Eq2);
ftotal = matlabFunction(V, 'Vars',{'t','Y','k1','k2','m'});
interval = [0 5];
y0 = [1 0; 0 0]; %initial conditions
ySol = ode45( @(t,Y)ftotal(t,Y,1,1,1),interval,y0);
% ftotal(t,Y,k1,k2,m)
tValues = linspace(interval(1),interval(2),5);
yValues = deval(ySol,tValues);
plot(tValues,yValues)
% plot(x,y)
I'm trying to numerically integrate and use the ode45 function to find the solution to this system of equations (Eq1 and Eq2). But somehow, my graphical solution is wrong and i don't get oscillations as expected.

请先登录，再回答此问题。

采纳的回答

David Goodmanson
David Goodmanson 2019-8-31
编辑：David Goodmanson 2019-8-31
Hi Ricardo,
I am inferring that you have the following system with fixed points S:
S---k1---M---k2---M---k1---S
Then eqn 1 is correct, but eqn 2 should be
Eq2 = D2x2 == ((k2*x1) - ((k1+k2)*x2))/m
i.e. with a minus sign. Then you get oscillations, which look a lot better when the linspace statement for tValues goes to, say, 100 points instead of 5.
The minus sign also makes eqn 2 symmetric with eqn 1 under interchange of x1 and x2, which reflects what is going on in the diagram.
If the system is actually S---k1---M---k2---M
with the right end free, then
Eq2 = D2x2 == ((k2*x1) - (k2*x2))/m
  1 个评论
David Goodmanson
David Goodmanson 2019-9-2
Hi Ricardo, it's a matter of selecting the right subarray from yValues. What have you tried so far?

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by