ODE System with 4 equations

Hi all,

I have a system with 4 ODEs which I want to solve simultanously.Each equations are feeded with some variables. All derivatives are with respect to time (t) only. The variables are x,v,p and u.

dx/dt = v(t)

dv/dt = - 2*v(t) - 1000*x(t) - p(t)

dp/dt = v(t) - u(t)

du/dt = p(t) - abs(u(t) * u(t)

Initial conditions are all zero at t = 0, i.e. x(0) = 0; v(0) = 0; p(0) = 0; u(0) = 0.

Looking forward to get your help.

I don't have any preference over the integration scheme but an application of ode45 should help. I also have access to the symbolic toolbox.

Best regards,

Baris

[Malak Abuhusien]

how solution with for loop?

Answer 1

[Josh Meyer]

When you have a system of equations, each equation gets its own spot in the solution vector y.

With the conventions

y(1) = x, dydt(1) = dx/dt

y(2) = v, dydt(2) = dv/dt

y(3) = p, dydt(3) = dp/dt

y(4) = u, dydt(4) = du/dt

You can write the system of equations in an ODE function as

function dydt = ODEsystem(t,y)
dydt = zeros(4,1);
dydt(1) = y(2);
dydt(2) = - 2*y(2) - 1000*y(1) - y(3);
dydt(3) =  y(2) - y(4);
dydt(4) =  y(3) - abs(y(4) * y(4));
end

After you save the function in a file in your current directory, you can set the initial conditions and integrate with:

y0 = zeros(4,1);
tspan = [0 10]; 
[t,y] = ode45(@ODEsystem,tspan,y0);
plot(t,y,'-o')

For your problem, with the initial conditions all zero, this integration doesn't do much because all of the terms in the equations depend on x, v, y, or p, so the terms all remain zero.

[Baris Gungordu]

Thank you very much. The original equations are much longer and I want to embed the variables inside the ODEsystem. When I included that in the function, I get the error of "Not enough input arguments"

Say, dydt(2) = -2*epsilon*y(2) ... and I want to embed epsilon which is pre-defined in the code.

function dydt = ODEsystem(t,y,eps)
dydt = zeros(4,1);
dydt(1) = y(2);
dydt(2) = -2*eps*y(2);
end

[Josh Meyer]

You can add additional inputs to the ODE function to pass in parameters, but ode45 still requires the function input to have only two inputs. So adding inputs requires that you create an anonymous function to pass to ode45:

fcn = @(t,y) ODEsystem(t,y,epsilon);
[t,y] = ode45(fcn,tspan,y0);

When fcn is created, it uses the value for "epsilon" in your workspace.

[Baris Gungordu]

Thank you. Strangely, when I have one single input parameter it's working. However, when I pass the second parameter, it does not recognise it any more. I get the error:

Unrecognized function or variable 'F_ext'.

Error in analytical_model_v2_edit>ODEsystem (line 82)
dydt(2) = -2*eps*omega*y(2) - omega^2 * y(1) + (F_ext*sin(2*pi*freq*t) - Ad*y(3))/mt;

Error in odearguments (line 90)
f0 = feval0,y0,args{:});   % ODE15I sets args{1} to yp0.

Error in ode45 (line 115)
  odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);

Error in analytical_model_v2_edit (line 34)
[t,y] = ode45(@ODEsystem,tspan,y0)ode is:

my code is:

fcn = @(t,y) ODEsystem(t,y,eps,F_ext); %eps,omega,V_c,Ao,Ad,K,F_ext,l_eff,gamma,p_amb,rho_amb,mt,freq);

[t,y] = ode45(@ODEsystem,tspan,y0);

[Josh Meyer]

This line

fcn = @(t,y) ODEsystem(t,y,eps,F_ext); 

defines a function handle that calls ODEsystem with two parameters and leaves the t,y inputs available for ode45 to use. So, when you call ode45, you don't want to specify "@ODEsystem" as the function, since that function takes 4 inputs. Instead, specify the function as "fcn" which you defined with the previous line:

[t,y] = ode45(fcn,tspan,y0);

If you want to use more parameters, you can change the "ODEsystem" function to accept more inputs, and then you just need to edit the "fcn = " line to call ODEsystem with the appropriate number of extra inputs to pass in the parameter values from your workspace.

[Baris Gungordu]

Great, many thanks.

[RITIKA Jaiswal]

what do do if we have odes of dimension 100.Since it was of order 4 we can easily write that but what if have order of 100 how can we implement that in our code?

please help.

[Torsten]

If there are regularities in the dydt terms, you can usually use a loop to set them up.

If not, you will have to write them down one by one.