Solving system of four second order ODEs

1 Ansicht (letzte 30 Tage)
Jake Barlow
Jake Barlow am 22 Jan. 2022
Kommentiert: Jake Barlow am 23 Jan. 2022
Hi,
I am trying to solve a system of four second order ordinary differential equations with the initial conditions . However after running the following code
syms U1(t) V1(t) U2(t) V2(t)
%Constants
alpha = 0.5;
dU1 = diff(U1, t);
dU2 = diff(U2, t);
dV1 = diff(V1, t);
dV2 = diff(V2, t);
t0 = 0; %initial time
tf = 5; %final time
%initial conditions
c1 = 1; %U(0)
c2 = 0; %U'(0)
c3 = 1; %V(0)
c4 = 0; %V'(0)
y0 = [c1 c2 c3 c4];
eq1 = diff(U2, 2) == alpha*dV2;
eq2 = diff(V2, 2) == -alpha*dU2;
eq3 = diff(U1, 2) == alpha*dV1;
eq4 = diff(V1, 2) == -alpha*dU1;
vars = [U1(t); V1(t); U2(t); V2(t)];
V = odeToVectorField([eq1,eq2,eq3,eq4]);
M = matlabFunction(V, 'vars', {'t','Y'});
interval = [t0 tf]; %time interval
ySol = ode45(M,interval,y0);
tValues = linspace(interval(1),interval(2),1000);
yValues1 = deval(ySol,tValues,1); %U(t)
yValues2 = deval(ySol,tValues,3); %V(t)
plot(yValues1, yValues2)
I get this error
Index exceeds the number of array elements. Index must not exceed 4.
Error in symengine>@(t,Y)[Y(2);Y(4).*(-1.0./2.0);Y(4);Y(2)./2.0;Y(6);Y(8)./2.0;Y(8);Y(6).*(-1.0./2.0)]
Error in odearguments (line 90)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 106)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
Could you please help? If there is another way to solve this sytem then it will be much appreciated if you can share it. Thank you.

Akzeptierte Antwort

Torsten
Torsten am 22 Jan. 2022
4 second-order ODEs need 8 initial conditions, not 4.
  3 Kommentare
Torsten
Torsten am 22 Jan. 2022
Use
[V,S] = odeToVectorField([eq1,eq2,eq3,eq4])
instead of
V = odeToVectorField([eq1,eq2,eq3,eq4])
to see in which order the solution variables appear.
Specify your initial conditions in this order, too.
Jake Barlow
Jake Barlow am 23 Jan. 2022
Hi Torsten, thank you very much. I specified the initial conditions in the correct order.

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (0)

Kategorien

Mehr zu Programming finden Sie in Help Center und File Exchange

Produkte


Version

R2021b

Community Treasure Hunt

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

Start Hunting!

Translated by