Pass extra variable parameters to ode15s function (MATLAB)

Question

I'm trying to solve a system of ordinary differential equations in MATLAB.

I have a simple equation: dy = -k/M *x - c/M *y+ F/M.

This is defined in my ode function test2.m, dependant on the values X and t. I want to trig 'F' with a signal, generated by my custom function squaresignal.m. The output hereof, is the variable u, spanding from 0 to 1, as it is a smooth heaviside function. - Think square wave. The inputs in squaresignal.m, is t and f.

u=squaresignal(t,f)

These values are to be used inside my function test2, in order to enable or disable variable 'F' with the value u==1 (enable). Disable for all other values.

My ode function test2.m reads:

function dX = test2(t ,X, u)
x = X (1) ;
y = X (2) ;

M = 10;
k = 50;
c = 10;
F = 300;

if u == 1
    F = F;
else
    F = 0,
end

dx = y ;
dy = -k/M *x - c/M *y+ F/M ;

dX = [ dx dy ]';
end

And my runscript reads:

clc
clear all

tstart = 0;
tend = 10;
tsteps = 0.01;

tspan = [0 10];
t = [tstart:tsteps:tend];
f = 2;

u = squaresignal(t,f)

for ii = 1:length(u)
    options=odeset('maxstep',tsteps,'outputfcn',@odeplot);
    [t,X]=ode15s(@(t,X)test2(t,X,u(ii)),[tstart tend],[0 0],u);
end


figure (1);
plot(t,X(:,1))

figure (2);
plot(t,X(:,2));

However, the for-loop does not seem to do it's magic. I still only get F=0, instead of F=F, at times when u==1. And i know, that u is equal to one at some times, because the output of squaresignal.m is visible to me.

So the real question is this. How do i properly pass my variable u, to my function test2.m, and use it there to trig F? Is it possible that the squaresignal.m should be inside the odefunction test2.m instead?

Answer 1

Here's an example where I pass a variable coeff to the differential equation:

function [T,Q] = main()
  t_span = [0 10];
  q0 = [0.1; 0.2];        % initial state
  ode_options = odeset(); % currently no options... You could add some here
  coeff = 0.3;            % The parameter we wish to pass to the differential eq.

  [T,Q] = ode15s(@(t,q)diffeq(t,q,coeff),t_span,q0, ode_options);
end

function dq = diffeq(t,q,coeff)
  % Preallocate vector dq
  dq = zeros(length(q),1);

  % Update dq:
  dq(1) = q(2);
  dq(2) = -coeff*sin(q(1));
end

EDIT: Could this be the problem?

tstart = 0;
tend = 10;
tsteps = 0.01;

tspan = [0 10];
t = [tstart:tsteps:tend];
f = 2;

u = squaresignal(t,f)

Here you create a time vector t which has nothing to do with the time vector returned by the ODE solver! This means that at first we have t[2]=0.01 but once you ran your ODE solver, t[2] can be anything. So yes, if you want to load an external signal source depending on time, then you need to call your squaresignal.m from within the differential equation and pass the solver's current time t! Your code should look like this (note that I'm passing f now as an additional argument to the diffeq):

function dX = test2(t ,X, f)
x = X (1) ;
y = X (2) ;

M = 10;
k = 50;
c = 10;
F = 300;

u = squaresignal(t,f)
if u == 1
    F = F;
else
    F = 0,
end

dx = y ;
dy = -k/M *x - c/M *y+ F/M ;

dX = [ dx dy ]';
end

Note however that matlab's ODE solvers do not like at all what you're doing here. You are drastically (i.e. non-smoothly) changing the dynamics of your system. What you should do is to use one of the following:

a) events if you want to trigger some behaviour (like termination) depending on the integrated variable x or

b) If you want to trigger the behaviour based on the time t, you should segment your integration into different parts where the differential equation does not vary during one segment. You can then resume your integration by using the current state and time as x0 and t0 for the next run of ode15s. Of course this only works of you're external signal source u is something simple like a step funcion or square wave. In case of the square wave you would only integrate for a timespan during which the wave does not jump. And then exactly at the time of the jump you start another integration with altered differential equations.