Prev: Positioning of colorbar
Next: Product of eigenvalues - not equal to its determinant
From: Matthieu Guérinel on 25 Feb 2010 07:46
Hello,
I'm having troubles trying to use the 'outputfcn' command of the options of the ode45 solver. Indeed, I have a model made of several .m files all linked together, and I need to perform differential equations calculus at the end in order to get my result.
The thing is that it worked very well so far (without the outputfcn command), but I would like to improve its accuracy by re-calculating a function (which in this case correspond to a surface) after each integration step, and insert it in the solver at the next time step.
I have an m-file in which I define the dydt function, of the the following form:
function dydt = odefun(t, y, k1, k2, k3, k4, k5, k6), where k1,k2...k6 are matrices or structures.

Then I create my solver with options as follow:

options =odeset('RelTol',1e-3,'AbsTol',1e-6,'InitialStep',1e-5,...
'OutputFcn',@myfun);

[t,y] = ode45(@odefun,tspan,y0,options,k1,k2,k3,k4,k5,k6);

myfun is the function of the form: function S = myfun(t,y,k4), that I want to call after each integration step, which depends on t,y as well as on k4, and which will give me the surface that I want to insert after in the solver.
My problem is that each time I launch the model, an error occurs tellig me :
"??? Error using ==> myfun
Too many input arguments. "
I already read some threads about this but still couldn't find a way to solve it...

I also called the myfun function directly in the odefun function, which worked, but I'd rather do it that way since it seems more correct to me.

Thanks in for any help.

Matthieu
From: Steven Lord on 25 Feb 2010 09:30

"Matthieu Guérinel" <mamat.guerinel(a)gmail.com> wrote in message
news:hm5rec$9cd$1(a)fred.mathworks.com...
> Hello,
> I'm having troubles trying to use the 'outputfcn' command of the options
> of the ode45 solver. Indeed, I have a model made of several .m files all
> linked together, and I need to perform differential equations calculus at
> the end in order to get my result.
> The thing is that it worked very well so far (without the outputfcn
> command), but I would like to improve its accuracy by re-calculating a
> function (which in this case correspond to a surface) after each
> integration step, and insert it in the solver at the next time step.
> I have an m-file in which I define the dydt function, of the the following
> form:
> function dydt = odefun(t, y, k1, k2, k3, k4, k5, k6), where k1,k2...k6 are
> matrices or structures.
>
> Then I create my solver with options as follow:
>
> options =odeset('RelTol',1e-3,'AbsTol',1e-6,'InitialStep',1e-5,...
> 'OutputFcn',@myfun);
>
> [t,y] = ode45(@odefun,tspan,y0,options,k1,k2,k3,k4,k5,k6);

If you try to pass additional parameters into odefun like this, then ALL of
the functions whose handles are stored in the options structure MUST accept
all those additional parameters as well, even if they don't need them.

Use anonymous functions instead to control exactly what each function
accepts.

k1= 1;
k2 = 2;
k3 = 3;
k4 = 4;
k5 = 5;
k6 = 6;
options =odeset('RelTol',1e-3,'AbsTol',1e-6,'InitialStep',1e-5,...
'OutputFcn',@(t, y) myfun(t, y, k4));
[t, y] = ode45(@(t, y) odefun(t, y, k1, k2, k3, k4, k5, k6), tspan, y0,
options);

Alternately, you could combine all six of your additional parameters into a
vector/matrix/cell array and have odefun and myfun extract only those pieces
that they need, and pass that "collection" into the solver as an additional
input.

--
Steve Lord
slord(a)mathworks.com
comp.soft-sys.matlab (CSSM) FAQ: http://matlabwiki.mathworks.com/MATLAB_FAQ


From: Matthieu Guérinel on 25 Feb 2010 11:40
Thank you very much for your answer!
However (and unfortunately), I still have some troubles and the same mistake still appears in the Matlab prompt after running.
Indeed, I already tried something similar regarding the anonymous functions that I saw in a previous thread in which you participated, but I think I omitted something in my first post.

I forgot to mention that in the output function that I am calling, the last parameter (here k4) is also depending on the integration time step t. This parameter k4 is a 3 columns matrix, in which the first column correspond to integer time values (e.g. 1,2,3...n), and the second and 3 columns to time specific values. As I want the corresponding values of these 2 last columns to the integration step time, I have to perform an interpolation in the "dydt" function (for the 2nd column), as well as in the "outputfcn" function (for the 3rd column).

That is the reason why I was wondering if the problem might come from here, because I don't really know the way the differential equation solver, but I have seen in some threads that attention should be paid when asking "extra" outputs in the way that it might false the final result.

As I mentioned before, I also wrote a similar m-file without using the "outputfcn" function, but using the myfun function inside the definition of the function "dydt", and it was running properly, but as I don't have any way to verify the veracity of the results I got, even though they look quite good, I prefer to try to use this "outputfcn" command, since Matlab built-in functions are usually more efficient...

You'll find hereinafter the real m-file I'm using, as well as the function I want to call which will probably make things clearer... The whole file is called as a function and is intended to solve the state-space model of my whole problem.

=> Ode solver file:


"function [T Y PTO MOR] = SystemSol(time,gp,ss,RMass,TFexc,PTO,MOR);


ii=length(gp.dof);
MASS=blkdiag(eye(gp.ord(1)-ii),RMass);
y0=[zeros(1,length(ss.Bg))];
MASS=inv(MASS);
options =odeset('RelTol',1e-3,'AbsTol',1e-6,'InitialStep',1e-5,...
'OutputFcn',@Swp_func);

[t,y] = ode45(@testdiff,time,y0,options,gp,ss,MASS,TFexc,PTO,MOR);


function dydt = testdiff(t,y,gp,ss,MASS,TFexc,PTO,MOR)

IFe=interp1(TFexc(:,1),TFexc(:,2),t,'linear','extrap');

SSInp=vertcat(zeros(1,length(ss.Bg)-length(IFe)).',IFe.');
dydt=MASS*(ss.Ag*y+eval(MOR.poly)+ss.Bg.*SSInp); "


Function to call:

"function Sinst = Swp_func(t,y,TFexc)

Rc=5; Hc=10; g=9.81; wd=1000;
tfsei=interp1(TFexc(:,1),TFexc(:,3),t,'linear','extrap');
zt=y(3)-tfsei;
rc = Rc*(1-abs(zt)./Hc);
Swp = (pi().*rc.^2)*g*wd;
if zt < 0
Sw = ((2*pi()*Rc*sqrt((Rc^2+Hc^2)))-pi()*rc.*sqrt((rc.^2+(Hc-abs(zt)).^2)))*g*wd;
elseif zt >= 0
Sw = (pi()*rc.*sqrt((rc.^2+(Hc-zt).^2)))*g*wd;
end

Sinst=[Swp Sw]; "

My final objective is to insert the new value of Swp in the next integration step of the solver, in the matrix "ss.Ag"

Thank you very much for any help you could give me.

Best regards,

Matthieu
 | 
Pages: 1
Prev: Positioning of colorbar
Next: Product of eigenvalues - not equal to its determinant