Problem using ode23, must return a column vector.

Hi, im trying to simulate 2 mass MDS system step response using ode23.

this is my state space:

K1, K2, M1, M2, C1, C2 are all known constants

global A B u
A = [0, 1, 0, 0; (-K1-K2)/M1, (-C1-C2)/M1, K2/M1, C2/M1;
    0, 0, 0, 1; K2/M2, C2/M2, -K2/M2, -C2/M2];

Unrecognized function or variable 'K1'.

B = [0; 1/M1; 0; 0];
C = [1, 0, 0, 0; 0, 0, 1, 0];
D = [0; 0];

and this is my ode23 call:

T_end = 10;
u = [1,1,1,1];  %(1, because step response)
x0 = [0, 0, 0, 0];  %(initial conditions)
[t,y] = ode23('racsimsemODE',[0,T_end],x0)

and my ode23 function:

function [dy] = racsimsemODE(t,x)
global A B u
dy = A*x + B*u;
end

however, this doesent work and I dont know why...

[Matt J]

It doesn't sound like you've actually checked whether racsimsemODE is returning a column vector...

Answer 1

[Star Strider]

Multiply ‘dy’ by an appropriate vector of ones and it works.

You need to determine if it produces the correct result.

% K1, K2, M1, M2, C1, C2 are all known constants
vv = num2cell(rand(6,1));                           % Define Constants
[K1, K2, M1, M2, C1, C2] = vv{:};                   % Assign Constants

A = [0, 1, 0, 0; (-K1-K2)/M1, (-C1-C2)/M1, K2/M1, C2/M1;
    0, 0, 0, 1; K2/M2, C2/M2, -K2/M2, -C2/M2];
B = [0; 1/M1; 0; 0];
C = [1, 0, 0, 0; 0, 0, 1, 0];
D = [0; 0];
% and this is my ode23 call:
T_end = 10;
u = [1,1,1,1];  % (1, because step response)
x0 = [0, 0, 0, 0];  % (initial conditions)
[t,y] = ode23(@(t,x)racsimsemODE(t,x,A,B,u),[0,T_end],x0);

figure
plot(t,y)
grid
legend(compose('$x_{%d}$',1:4), 'Location','best', 'Interpreter','latex')

% and my ode23 function:
function [dy] = racsimsemODE(t,x,A,B,u)
dy = A*x + B*u;
dy = dy*ones(4,1);
end

.

[Nino Dragicevic]

Thanks, but I've figured out my stupid mistake... u isnt a vector, it should be just a constant... So I've replaced u vector with u = 1, and it works.

[Star Strider]

As always, my pleasure!

(I went with the code as posted.)

[Nino Dragicevic]

Do you maybe have an idea how do I modify the code to show sine response? In a way that the input function u(t) = sin(omega * t)?

[Star Strider]

There are several possibilities.

This is one option:

% K1, K2, M1, M2, C1, C2 are all known constants
vv = num2cell(rand(6,1));                           % Define Constants
[K1, K2, M1, M2, C1, C2] = vv{:};                   % Assign Constants

A = [0, 1, 0, 0; (-K1-K2)/M1, (-C1-C2)/M1, K2/M1, C2/M1;
    0, 0, 0, 1; K2/M2, C2/M2, -K2/M2, -C2/M2];
B = [0; 1/M1; 0; 0];
C = [1, 0, 0, 0; 0, 0, 1, 0];
D = [0; 0];
% and this is my ode23 call:
T_end = 10;
% u = [1,1,1,1];  % (1, because step response)
omega = 2*pi*0.5;                                       % Define 'omega'
u = @(t) sin(omega*t);                                  % Define 'u(t)'
x0 = [0, 0, 0, 0];  % (initial conditions)
[t,y] = ode23(@(t,x)racsimsemODE(t,x,A,B,u),[0,T_end],x0);

figure
plot(t,y)
grid
legend(compose('$x_{%d}$',1:4), 'Location','best', 'Interpreter','latex')

% and my ode23 function:
function [dy] = racsimsemODE(t,x,A,B,u)
dy = A*x + B*u(t);                                      % Call 'u(t)'
end

If ‘u’ is defined as a vector instead of a function (as I did here), it would be necesary to interpolate it instead. See the documentation section on ODE with Time-Dependent Terms for details.

.

[Nino Dragicevic]

Thanks alot, you just saved my semester!

[Star Strider]

As always, my pleasure!