9 views (last 30 days)

Show older comments

Hey guys, I need to solve a coupled, simultaneous and implicit system of ODE's. Here's my code and then I'll briefly explain it:

close all

clear all

time_length = 25;

g = 9.8; % Acceleration due to gravity (m/s^2)

rho0 = 1.20; % Density of air (kg/m^3)

d = 0.22; % Ball diameter (m)

m = 0.43; % Ball mass (kg)

A = pi*(d/2)^2; % Ball cross-sectional area (m^2)

Cd = 0.3; % Drag coefficient

Cl = 0.3; % Lift coefficient

Sx = 1; % x-component of S

Sy = 1; % y-component of S

Sz = 0; % z-component of S

S = [Sx, Sy, Sz]; % Spin vector

% Coupled ODE's. x = x(1), vx = x(2). Same for y and z

odex = @(t,x,y,z) -Cd*(A/(2*m))*rho0*sqrt(x(2).^2 + y(2).^2 + z(2).^2).*x(2) + Cl*(A/(2*m))*rho0*sqrt(x(2).^2 + y(2).^2 + z(2).^2).*(Sy.*z(2)-Sz.*y(2));

odey = @(t,x,y,z) -Cd*(A/(2*m))*rho0*sqrt(x(2).^2 + y(2).^2 + z(2).^2).*y(2) + Cl*(A/(2*m))*rho0*sqrt(x(2).^2 + y(2).^2 + z(2).^2).*(Sz.*x(2)-Sx.*z(2));

odez = @(t,x,y,z) -g-Cd*(A/(2*m))*rho0*sqrt(x(2).^2 + y(2).^2 + z(2).^2).*z(2) + Cl*(A/(2*m))*rho0*sqrt(x(2).^2 + y(2).^2 + z(2).^2).*(Sx.*y(2)-Sy.*x(2));

W = @(x,y,z) [x(1);x(2);y(1);y(2);z(1);z(2)];

dxyz = @(t,x,y,z) [x(2);odex;y(2);odey;z(2);odez];

[T,W] = ode45(dxyz,linspace(0,time_length,1e4),[0,30/sqrt(2),0,0,30/sqrt(2)]);

So basically:

The notation is as follows:

x(1) = x,

x(2) = derivative of x = dx/dt = vx,

odex = second derivative of x = derivative of vx = d2x/dt2 = d(vx)/dt = dx(2)/dt = ax. Same for y and z.

Horrible equations in odex/y/z but basically each depend on the other variables, i.e. odex doesn't only depend on x, but y and z too (and their derivatives)

This last point means I need to solve them via matrix algebra. W is the column matrix that I'm planning to solve for, and dxyz is the column matrix containing the derivatives, i.e. dW/dt = dxyz.

in the ode45 bracket I have 6 initial conditions, corresponding to x(1) -> z(2) in order.

And now when I run it, I get these outputs:

Error using @(t,x,y,z)[x(2);odex;y(2);odey;z(2);odez]

Not enough input arguments.

Error in odearguments (line 88)

f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.

Error in ode45 (line 114)

[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...

Error in code (line 37)

[T,W] = ode45(dxyz,linspace(0,time_length,1e4),[0,30/sqrt(2),0,0,30/sqrt(2)]);

It's telling me not enough input arguments, I do not understand. I have everything in there that I need, don't I?

Also I'm getting the idea that all these errors correspond to just the one error, i.e. if I fix the one thing that's wrong, I'll fix it all. Unfortunately, I don't know what to do about it because I don't understand why I'm getting an error.

Any and all help is much appreciated!

Geoff Hayes
on 29 Sep 2014

Hi Joshua - while I haven't used ode45, I will try to give some hints at what may be causing the errors.

The message Not enough input arguments means that "someone" is not supplying enough input arguments to "some" function. In this case, that would be the ode45 not providing enough inputs to the dxyz function. Look at the documentation for ode45 and in particular the details on the function handle input parameter, odefun, which states

A function handle that evaluates the right side of the differential equations. All solvers solve systems of equations in the form y′ = f(t,y) or problems that involve a mass matrix, M(t,y)y′ = f(t,y). The ode23s solver can solve only equations with constant mass matrices. ode15s and ode23t can solve problems with a mass matrix that is singular, i.e., differential-algebraic equations (DAEs).

So the function must have a signature that expects two input parameters, t and y. Note how your function signature for dxyz requires four input parameters. Hence the error. So either your function signature will need to change or you will have to adapt the code to something similar to the second example (from the ode45 link).

Now note how dxyz is defined

dxyz = @(t,x,y,z) [x(2);odex;y(2);odey;z(2);odez];

It produces a 6x1 vector, calling odex, odey, and odez. But the code is not supplying any input parameters to these functions so would generate an error. I observed the following errors

Error using vertcat

The following error occurred converting from double to function_handle:

Error using function_handle

Too many output arguments.

So you will need to change the above so that (whatever) inputs are passed into to these functions

dxyz = @(t,x,y,z) [x(2);odex(t,x,y,z);y(2);odey(t,x,y,z);z(2);odez(t,x,y,z)];

As well, the initial condition vector, as supplied to ode45, has only five elements rather than six

[0,30/sqrt(2),0,0,30/sqrt(2)]

-------------

The above addresses some concerns with the code, but not with whether it is being used as it should be. For example, why does the initial conditions vector have six elements when the three position elements (x,y,z) are unused? Why can't we just supply the velocities, and the dxyz return a 3x1 vector instead of the 6x1?

Kate
on 12 Feb 2018

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

Start Hunting!