ODE function will not call the OutputFcn

I want to access specific variables which have been solved for in the ODE function but are not differential equations. I am setting an 'OutputFcn', however it is not being accessed while the program is being run. Here are the important parts of my code:

 options = odeset('OutputFcn',@(status)RotorOut(t,y,flag));
 [t,y] = ode45(@QR_solver,delt,y0,options,[],D);
--- 
 function dydt = QR_solver(t,y,D)
   ...
   rs = XX;
   ...
 end
---
 function status = RotorOut(t,y,flag)
   global rs k

   if flag ~= strcmp(flag,'init')|strcmp(flag,'done')
     rotspeed = evalin('base','rotspeed');
     k = k + 1;
     rotspeed(k) = rs;
     assignin('base','rotspeed',rotspeed)
     status = 0;
   end
   t    % just to verify the function is being called
 end

Am I calling the function wrong? Is there a better to access the variable rs?

Answer 1

[Brian B]

1. The OutputFcn should have the format

 status = myfun(t,y,flag)

but your anonymous function takes a single input (status). Try

 options = odeset('OutputFcn',@RotorOut);

2. Perhaps it is hidden in the ellipses, but you don't seem to declare rs as global in QR_solver().

Note, however, that the time points at which the output from ode45 is returned may be different from the points at which the functions themselves are evaluated. So if you want to save intermediate results for later use, you will probably have to interpolate afterward, or just loop through the ode45 output and recompute the intermediate results.

I definitely wish there were a way to have ode45 compute auxiliary outputs!

A somewhat cleaner way to share the data between QR_solver() and RotorOut() would be to use nested functions:

 [ode, out] = getHandles();
 options = odeset('OutputFcn',@(status)RotorOut(t,y,flag));
 [t,y] = ode45(@QR_solver,delt,y0,options,[],D);

 function [ode, out] = getHandles()
   ode = @QR_solver;
   out = @RotorOut;

   rs = []; % make this variable shared

   function dydt = QR_solver(t,y,D)
     ...
     rs = XX;
     ...
   end
   ---
   function status = RotorOut(t,y,flag)

     if flag ~= strcmp(flag,'init')|strcmp(flag,'done')
       rotspeed = evalin('base','rotspeed');
       k = k + 1;
       rotspeed(k) = rs;
       assignin('base','rotspeed',rotspeed)
       status = 0;
     end
     t    % just to verify the function is being called
   end
 end

[Norman Johnson]

Brian, Thanks for the quick response. I got my program to work. I like your suggestion 2) however I am not very familiar with nested functions and so I have not tried it. How should I make rs shared? It looks like I call the ode45 function before I set the handles, is that correct?

[Brian B]

Nested functions take some getting used to! The rule is that a variable exists in the workspace of the outer-most function that accesses it. Thus the declaration

 rs = [];

in getHandles() makes rs shared between getHandles(), QR_solver(), and RotorOut(), since the two nested functions also access that same variable.

I forgot to change the function handles in the solver call. It should be

 [ode, out] = getHandles();
 options = odeset('OutputFcn',out);
 [t,y] = ode45(ode,delt,y0,options,[],D);

Thus the handles ode and out are created before the call to ode45. For more details on variable and function scope, see http://www.mathworks.com/help/matlab/matlab_prog/nested-functions.html.