Solving system of four second order ODEs

1 view (last 30 days)
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.

Accepted Answer

Torsten
Torsten on 22 Jan 2022
4 second-order ODEs need 8 initial conditions, not 4.
  3 Comments
Torsten
Torsten on 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 on 23 Jan 2022
Hi Torsten, thank you very much. I specified the initial conditions in the correct order.

Sign in to comment.

More Answers (0)

Categories

Find more on Programming in Help Center and File Exchange

Products


Release

R2021b

Community Treasure Hunt

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

Start Hunting!