Re: passing parameters to lsode

[Tatsuro MATSUOKA]

--- "address@hidden" <address@hidden> wrote:
> Hi -
> this is my first message here. Hope I do this right  ;)  Question: How do I 
> pass parameters to a
> differential equation? lsode wants me to declare the diff. eq. as function 
> xdot =
> something(x,t). i would rather like to declare this as 
> xdot=something(x,t,parameter1,...). I
> help myself by declaring global variables for every parameter. This doesn't 
> look very elegant
> and transparent to me.
> 
> Am I missing something?

For fsolve expression, 'xdot=something(x,t,parameter1,...). ' can be used now. 
Please see the following and its succeeding threads:

http://www.cae.wisc.edu/pipermail/help-octave/2007-August/005444.html

I do not not know, the same representation is allowed for lsode of current 
octave.
One of the solution, David showed in the succeeding thread.

>I help myself by declaring global variables for every parameter. This doesn't 
>look very elegant
> and transparent to me.

At least octave 2.1.xx,
your solution is the only way to pass the parameter to the function.
An example code is written at the end of mail.

This is the all I can answer now.

Hi, David
Do you know parameter treating extension to be applied for the lsode in the 
current octave?
  

*****************
%Octave script m file ; incompatible with MATLAB
% Analysis of non-constant temperature batch reactor
subplot(111);clg;
format;
%****function
function dT=NTCBRf(T,xA)
global p ;
% p : A parameter vector
dT(1)=p(1)+p(2)*exp(p(3)/T(1))/(1-xA)*(p(4)-T(1));
dT(2)=p(5)*exp(p(3)/T(1))/(1-xA);
end
%
%**** Main script
global p
% p : A parameter vector
% constants
p = [ -64.883 3.8828e-17 2.2396e+04 620 3.098e-14];
Ts=620;
xA=0:0.01:0.99; % A vector of conversions to be calclated.
T0(1)=Ts; T0(2)=0; % T0(1) : initial temperature, T0(2) : start time
% ODE Solver 'lsode'
T=lsode('NTCBRf', T0, xA);
% T is a Matrix
% first column T(:,1) is a temperature vector corresponding to the xA vector
% second column T(:,2) is a time vector corresponding to the xA vector
%
fprintf("Solutions of the problem are :\n");
fprintf(" Time required when xA reaches 0.97 : t_0.97 = %7.2f s \n", T(98,2));
fprintf(" Temperature at the time : T_0.97 = %7.2f K \n", T(98,1));
disp('');
%plot graph
msg = 'Time dependence of conversion (upper side) and temperature (lower side) 
are
plotted in graphs.';
msg =[msg ' (t_0.97, xA=0.97) and (t_0.97, T_0.97) are marked as points in the
graphs.'];
disp(msg);
axis('auto');
subplot(211);
ylabel('xA / -');
pmsg=sprintf("(%7.2f, 0.97)", T(98,2));
plot(T(:,2),xA,";xA;",T(98,2),xA(98),["@32;" pmsg ";"]);
subplot(212);
xlabel('time / s');
ylabel('T / K');
pmsg=sprintf("(%7.2f, %7.2f)", T(98,2), T(98,1));
plot(T(:,2),T(:,1),";T;",T(98,2),T(98,1),["@32;" pmsg ";"]);


--------------------------------------
Easy + Joy + Powerful = Yahoo! Bookmarks x Toolbar
http://pr.mail.yahoo.co.jp/toolbar/