<itemizedlist>
<listitem>
<para>In Fortran the calling sequence must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine fsub(x,zu,f)
double precision zu(*), f(*),x
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>In C the function prototype must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void fsub(double *x, double *zu, double *f)
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>And in Scilab:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-function f=fsub(x,zu,parameters)
- ]]></programlisting>
+ <screen><![CDATA[
+function f = fsub(x,zu,parameters)
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<itemizedlist>
<listitem>
<para>In Fortran the calling sequence must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine dfsub(x,zu,df)
double precision zu(*), df(*),x
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>In C the function prototype must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void dfsub(double *x, double *zu, double *df)
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>And in Scilab:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-function df=dfsub(x,zu,parameters)
- ]]></programlisting>
+ <screen><![CDATA[
+function df = dfsub(x,zu,parameters)
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<itemizedlist>
<listitem>
<para>In Fortran the calling sequence must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine gsub(i,zu,g)
double precision zu(*), g(*)
integer i
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>In C the function prototype must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void gsub(int *i, double *zu, double *g)
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>And in Scilab:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-function g=gsub(i,zu,parameters)
- ]]></programlisting>
+ <screen><![CDATA[
+function g = gsub(i,zu,parameters)
+]]></screen>
<para>
Note that in contrast to <literal>f</literal> in <literal>fsub</literal>, here only one value per call is returned in <literal>g</literal>.
</para>
<itemizedlist>
<listitem>
<para>In Fortran the calling sequence must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine dgsub(i,zu,dg)
double precision zu(*), dg(*)
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>In C the function prototype must be</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void dgsub(int *i, double *zu, double *dg)
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>And in Scilab </para>
- <programlisting role="no-scilab-exec"><![CDATA[
-function dg=dgsub(i,zu,parameters)
- ]]></programlisting>
+ <screen><![CDATA[
+function dg = dgsub(i,zu,parameters)
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<itemizedlist>
<listitem>
<para>In Fortran the calling sequence must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine guess(x,zu,dmval)
double precision x,z(*), dmval(*)
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>In C the function prototype must be</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void fsub(double *x, double *zu, double *dmval)
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>And in Scilab </para>
- <programlisting role="no-scilab-exec"><![CDATA[
-function [dmval,zu]=fsub(x,parameters)
- ]]></programlisting>
+ <screen><![CDATA[
+function [dmval,zu] = fsub(x,parameters)
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
</imageobject>
</mediaobject>
<programlisting role="example"><![CDATA[
-N=1;// just one differential equation
-m=4;//a fourth order differential equation
-M=sum(m);
+N = 1; // just one differential equation
+m = 4; // a fourth order differential equation
+M = sum(m);
-x_low=1;
-x_up=2; // the x limits
-zeta=[x_low,x_low,x_up,x_up]; //two constraints (on the value of u and its second derivative) on each bound.
+x_low = 1;
+x_up = 2; // the x limits
+zeta = [x_low,x_low,x_up,x_up]; // two constraints (on the value of u and its second derivative) on each bound.
-//The external functions
-//These functions are called by the solver with zu=[u(x);u'(x);u''(x);u'''(x)]
+// The external functions
+// These functions are called by the solver with zu=[u(x);u'(x);u''(x);u'''(x)]
// - The function which computes the right hand side of the differential equation
-function f=fsub(x,zu)
- f=(1-6*x^2*zu(4)-6*x*zu(3))/x^3
+function f = fsub(x,zu)
+ f = (1-6*x^2*zu(4)-6*x*zu(3))/x^3
endfunction
// - The function which computes the derivative of fsub with respect to zu
-function df=dfsub(x,zu)
- df=[0,0,-6/x^2,-6/x]
+function df = dfsub(x,zu)
+ df = [0,0,-6/x^2,-6/x]
endfunction
// - The function which computes the ith constraint for a given i
-function g=gsub(i,zu),
+function g = gsub(i,zu),
select i
- case 1 then //x=zeta(1)=1
- g=zu(1) //u(1)=0
- case 2 then //x=zeta(2)=1
- g=zu(3) //u''(1)=0
- case 3 then //x=zeta(3)=2
- g=zu(1) //u(2)=0
- case 4 then //x=zeta(4)=2
- g=zu(3) //u''(2)=0
+ case 1 then // x = zeta(1) = 1
+ g = zu(1) // u(1) = 0
+ case 2 then // x = zeta(2) = 1
+ g = zu(3) // u''(1) = 0
+ case 3 then // x = zeta(3) = 2
+ g = zu(1) // u(2) = 0
+ case 4 then // x = zeta(4) = 2
+ g = zu(3) // u''(2) = 0
end
endfunction
// - The function which computes the derivative of gsub with respect to z
-function dg=dgsub(i,z)
+function dg = dgsub(i,z)
select i
- case 1 then //x=zeta(1)=1
- dg=[1,0,0,0]
- case 2 then //x=zeta(2)=1
- dg=[0,0,1,0]
- case 3 then //x=zeta(3)=2
- dg=[1,0,0,0]
- case 4 then //x=zeta(4)=2
- dg=[0,0,1,0]
+ case 1 then // x = zeta(1) = 1
+ dg = [1,0,0,0]
+ case 2 then // x = zeta(2) = 1
+ dg = [0,0,1,0]
+ case 3 then // x = zeta(3) = 2
+ dg = [1,0,0,0]
+ case 4 then // x = zeta(4) = 2
+ dg = [0,0,1,0]
end
endfunction
// - The function which computes the initial guess, unused here
-function [zu,mpar]=guess(x)
- zu=0;
- mpar=0;
+function [zu,mpar] = guess(x)
+ zu = 0;
+ mpar = 0;
endfunction
-//define the function which computes the exact value of u for a given x ( for testing purposes)
-function zu=trusol(x)
- zu=0*ones(4,1)
+// define the function which computes the exact value of u for a given x ( for testing purposes)
+function zu = trusol(x)
+ zu = 0*ones(4,1)
zu(1) = 0.25*(10*log(2)-3)*(1-x) + 0.5 *( 1/x + (3+x)*log(x) - x)
zu(2) = -0.25*(10*log(2)-3) + 0.5 *(-1/x^2 + (3+x)/x + log(x) - 1)
zu(3) = 0.5*( 2/x^3 + 1/x - 3/x^2)
zu(4) = 0.5*(-6/x^4 - 1/x/x + 6/x^3)
endfunction
-fixpnt=[ ];//All boundary conditions are located at x_low and x_up
+fixpnt = []; // All boundary conditions are located at x_low and x_up
-// nonlin collpnt n ntol ndimf ndimi iprint iread iguess rstart nfxpnt
-ipar=[0 0 1 2 2000 200 1 0 0 0 0 ]
+// nonlin collpnt n ntol ndimf ndimi iprint iread iguess rstart nfxpnt
+ipar = [0 0 1 2 2000 200 1 0 0 0 0 ]
-ltol=[1,3];//set tolerance control on zu(1) and zu(3)
-tol=[1.e-11,1.e-11];//set tolerance values for these two controls
-xpoints=x_low:0.01:x_up;
+ltol = [1,3]; // set tolerance control on zu(1) and zu(3)
+tol = [1.e-11,1.e-11]; // set tolerance values for these two controls
+xpoints = x_low:0.01:x_up;
-zu=bvode(xpoints,N,m,x_low,x_up,zeta,ipar,ltol,tol,fixpnt,...
- fsub,dfsub,gsub,dgsub,guess)
-//check the constraints
-zu([1,3],[1 $]) //should be zero
+zu = bvode(xpoints,N,m,x_low,x_up,zeta,ipar,ltol,tol,fixpnt,...
+ fsub,dfsub,gsub,dgsub,guess)
+// check the constraints
+zu([1,3],[1 $]) // should be zero
plot(xpoints,zu(1,:)) // the evolution of the solution u
-zu1=[];
-for x=xpoints
- zu1=[zu1,trusol(x)];
-end;
+zu1 = [];
+for x = xpoints
+ zu1 = [zu1,trusol(x)];
+end
norm(zu-zu1)
- ]]></programlisting>
+ ]]></programlisting>
</listitem>
<listitem>
<para>
Same problem using <literal>bvodeS</literal> and an initial guess.
</para>
<programlisting role="no-scilab-exec"><![CDATA[
-function [z,lhS]=zstart(x)
- z=zeros(5,1);z(5)=1;
- lhS=[0;1];
+function [z,lhS] = zstart(x)
+ z = zeros(5,1); z(5) = 1;
+ lhS = [0;1];
endfunction
-zu=bvode(xpoints,N,m,x_low,x_up,zeta,ltol=ltol,tol=tol,guess=zstart)
+zu = bvode(xpoints, N, m, x_low,x_up, zeta, ltol=ltol,tol=tol, guess=zstart)
]]></programlisting>
</listitem>
<listitem>
thickness subject to a quadratically varying axisymmetric external pressure
distribution. Here <emphasis>φ</emphasis> is the meridian angle change of the
deformed shell and <emphasis>ψ</emphasis> is a stress function.
- For <emphasis>ε = μ = 10<superscript>-3</superscript></emphasis>, two
+ For <emphasis>ε = μ = 10<superscript>-3</superscript></emphasis>, two
different solutions may found depending on the starting point
</para>
<mediaobject>
for <literal>x=0</literal> and <literal>x=1</literal>
</para>
<programlisting role="example"><![CDATA[
-N=2;// two differential equations
-m=[2 2];//each differential equation is of second order
-M=sum(m);
+N = 2; // two differential equations
+m = [2 2]; // each differential equation is of second order
+M = sum(m);
-x_low=0;x_up=1; // the x limits
-zeta=[x_low,x_low, x_up x_up]; //two constraints on each bound.
+x_low = 0;x_up=1; // the x limits
+zeta = [x_low,x_low, x_up x_up]; // two constraints on each bound.
-//The external functions
-//These functions are called by the solver with zu=[u1(x);u1'(x);u2(x);u2'(x)]
+// The external functions
+// These functions are called by the solver with zu=[u1(x);u1'(x);u2(x);u2'(x)]
// - The function which computes the right hand side of the differential equation
-function f=fsub2(x,zu,eps,dmu,eps4mu,gam,xt),
- f=[zu(1)/x^2-zu(2)/x+(zu(1)-zu(3)*(1-zu(1)/x)-gam*x*(1-x^2/2))/eps4mu //phi''
- zu(3)/x^2-zu(4)/x+zu(1)*(1-zu(1)/(2*x))/dmu];//psi''
+function f = fsub2(x,zu,eps,dmu,eps4mu,gam,xt),
+ f = [zu(1)/x^2-zu(2)/x+(zu(1)-zu(3)*(1-zu(1)/x)-gam*x*(1-x^2/2))/eps4mu //phi''
+ zu(3)/x^2-zu(4)/x+zu(1)*(1-zu(1)/(2*x))/dmu];//psi''
endfunction
// - The function which computes the derivative of fsub with respect to zu
-function df=dfsub2(x,zu,eps,dmu,eps4mu,gam,xt),
- df=[1/x^2+(1+zu(3)/x)/eps4mu, -1/x, -(1-zu(1)/x)/eps4mu, 0
- (1-zu(1)/x)/dmu 0 1/x^2 -1/x];
+function df = dfsub2(x,zu,eps,dmu,eps4mu,gam,xt),
+ df = [1/x^2+(1+zu(3)/x)/eps4mu, -1/x, -(1-zu(1)/x)/eps4mu, 0
+ (1-zu(1)/x)/dmu 0 1/x^2 -1/x];
endfunction
// - The function which computes the ith constraint for a given i
-function g=gsub2(i,zu),
+function g = gsub2(i,zu),
select i
- case 1 then //x=zeta(1)=0
- g=zu(1) //u(0)=0
- case 2 then //x=zeta(2)=0
- g=-0.3*zu(3) //x*psi'-0.3*psi+0.7x=0
- case 3 then //x=zeta(3)=1
- g=zu(1) //u(1)=0
- case 4 then //x=zeta(4)=1
- g=1*zu(4)-0.3*zu(3)+0.7*1 //x*psi'-0.3*psi+0.7x=0
+ case 1 then // x = zeta(1) = 0
+ g = zu(1) // u(0) = 0
+ case 2 then // x = zeta(2) = 0
+ g = -0.3*zu(3) // x*psi'-0.3*psi+0.7x = 0
+ case 3 then // x = zeta(3) = 1
+ g = zu(1) // u(1) = 0
+ case 4 then // x = zeta(4) = 1
+ g = 1*zu(4)-0.3*zu(3)+0.7*1 // x*psi'-0.3*psi+0.7x = 0
end
endfunction
// - The function which computes the derivative of gsub with respect to z
-function dg=dgsub2(i,z)
+function dg = dgsub2(i,z)
select i
- case 1 then //x=zeta(1)=1
- dg=[1,0,0,0]
- case 2 then //x=zeta(2)=1
- dg=[0,0,-0.3,0]
- case 3 then //x=zeta(3)=2
- dg=[1,0,0,0]
- case 4 then //x=zeta(4)=2
- dg=[0,0,-0.3,1]
+ case 1 then // x = zeta(1) = 1
+ dg = [1,0,0,0]
+ case 2 then // x = zeta(2) = 1
+ dg = [0,0,-0.3,0]
+ case 3 then // x = zeta(3) = 2
+ dg = [1,0,0,0]
+ case 4 then // x = zeta(4) = 2
+ dg = [0,0,-0.3,1]
end
endfunction
-gam=1.1
-eps=1d-3
-dmu=eps
-eps4mu=eps^4/dmu
-xt=sqrt(2*(gam-1)/gam)
+[gam, eps, dmu] = (1.1, 1d-3, 1d-3);
+eps4mu = eps^4/dmu
+xt = sqrt(2*(gam-1)/gam)
-fixpnt=[ ];//All boundary conditions are located at x_low and x_up
-collpnt=4;
-nsizef=4+3*M+(5+collpnt*N)*(collpnt*N+M)+(2*M-2)*2*M ;
-nsizei=3 + collpnt*N+M;
-nmax=200;
-// nonlin collpnt n ntol ndimf ndimi iprint iread iguess rstart nfxpnt
-ipar=[1 collpnt 10 4 nmax*nsizef nmax*nsizei -1 0 0 0 0 ]
+fixpnt = []; // All boundary conditions are located at x_low and x_up
+collpnt = 4;
+nsizef = 4 + 3*M + (5+collpnt*N)*(collpnt*N+M)+(2*M-2)*2*M ;
+nsizei = 3 + collpnt*N+M;
+nmax = 200;
+// nonlin collpnt n ntol ndimf ndimi iprint iread iguess rstart nfxpnt
+ipar = [1 collpnt 10 4 nmax*nsizef nmax*nsizei -1 0 0 0 0 ];
-ltol=1:4;//set tolerance control on zu(1), zu(2), zu(3) and zu(4)
-tol=[1.e-5,1.e-5,1.e-5,1.e-5];//set tolreance values for these four controls
-xpoints=x_low:0.01:x_up;
+ltol = 1:4; // set tolerance control on zu(1), zu(2), zu(3) and zu(4)
+tol = [1.e-5,1.e-5,1.e-5,1.e-5]; // set tolreance values for these four controls
+xpoints = x_low:0.01:x_up;
// - The function which computes the initial guess, unused here
-function [zu,dmval]=guess2(x,gam),
- cons=gam*x*(1-x^2/2)
- dcons=gam*(1-3*x^2/2)
- d2cons=-3*gam*x
- dmval=zeros(2,1)
+function [zu,dmval] = guess2(x,gam),
+ cons = gam*x*(1-x^2/2)
+ dcons = gam*(1-3*x^2/2)
+ d2cons = -3*gam*x
+ dmval = zeros(2,1)
if x>xt then
- zu=[0 0 -cons -dcons]
- dmval(2)=-d2cons
+ zu = [0 0 -cons -dcons]
+ dmval(2) = -d2cons
else
- zu=[2*x;2;-2*x+cons;-2*dcons]
- dmval(2)=d2cons
+ zu = [2*x;2; -2*x+cons; -2*dcons]
+ dmval(2) = d2cons
end
endfunction
-zu=bvode(xpoints,N,m,x_low,x_up,zeta,ipar,ltol,tol,fixpnt,...
- fsub2,dfsub2,gsub2,dgsub2,guess2);
-scf(1);clf();plot(xpoints,zu([1 3],:)) // the evolution of the solution phi and psi
+zu = bvode(xpoints,N,m,x_low,x_up,zeta,ipar,ltol,tol,fixpnt,...
+ fsub2,dfsub2,gsub2,dgsub2,guess2);
+scf(1); clf(); plot(xpoints,zu([1 3],:)) // the evolution of the solution phi and psi
-//using an initial guess
-ipar(9)=1;//iguess
+// using an initial guess
+ipar(9) = 1; // iguess
-zu2=bvode(xpoints,N,m,x_low,x_up,zeta,ipar,ltol,tol,fixpnt,...
- fsub2,dfsub2,gsub2,dgsub2,guess2);
-scf(2);clf();plot(xpoints,zu2([1 3],:)) // the evolution of the solution phi and psi
- ]]></programlisting>
+zu2 = bvode(xpoints,N,m,x_low,x_up,zeta,ipar,ltol,tol,fixpnt,...
+ fsub2,dfsub2,gsub2,dgsub2,guess2);
+scf(2); clf(); plot(xpoints,zu2([1 3],:)) // the evolution of the solution phi and psi
+ ]]></programlisting>
</listitem>
<listitem>
<para>
// BV: y(0)=y'(0), y(0)=1; y(1)=0
// z=[y(x) ; y'(x) ; la]
-function rhs=fsub(x,z)
- rhs=[-z(3)*z(1);0]
+function rhs = fsub(x,z)
+ rhs = [-z(3)*z(1);0]
endfunction
-function g=gsub(i,z)
- g=[z(1)-z(2) z(1)-1 z(1)]
- g=g(i)
+function g = gsub(i,z)
+ g = [z(1)-z(2) z(1)-1 z(1)]
+ g = g(i)
endfunction
// The following start function is good for the first 8 eigenfunctions.
-function [z,lhs]=ystart(x,z,la0)
- z=[1;0;la0]
- lhs=[0;0]
+function [z,lhs] = ystart(x,z,la0)
+ z = [1;0;la0]
+ lhs = [0;0]
endfunction
-a=0;b=1;
-m=[2;1];
-n=2;
-zeta=[a a b];
-N=101;
-x=linspace(a,b,N)';
+[a, b, m, n] =(0, 1, [2;1], 2);
+zeta = [a a b];
+N = 101;
+x = linspace(a,b,N)';
// We have s(n)-(n+1/2)*pi -> 0 for n to infinity.
-la0=evstr(x_dialog('n-th eigenvalue: n= ?','10'));
-la0=(%pi/2+la0*%pi)^2;
+la0 = evstr(x_dialog('n-th eigenvalue: n= ?','10'));
+la0 = (%pi/2+la0*%pi)^2;
-z=bvodeS(x,m,n,a,b,fsub,gsub,zeta,ystart=list(ystart,la0));
+z = bvodeS(x,m,n,a,b,fsub,gsub,zeta,ystart=list(ystart,la0));
// The same call without any display
-z=bvodeS(x,m,n,a,b,fsub,gsub,zeta,ystart=list(ystart,la0),iprint=1);
+z = bvodeS(x,m,n,a,b,fsub,gsub,zeta,ystart=list(ystart,la0),iprint=1);
// The same with a lot of display
-z=bvodeS(x,m,n,a,b,fsub,gsub,zeta,ystart=list(ystart,la0),iprint=-1);
+z = bvodeS(x,m,n,a,b,fsub,gsub,zeta,ystart=list(ystart,la0),iprint=-1);
clf()
plot(x,[z(1,:)' z(2,:)'])
xtitle(['Startvalue = '+string(la0);'Eigenvalue = '+string(z(3,1))],'x',' ')
legend(['y(x)';'y''(x)']);
- ]]></programlisting>
+ ]]></programlisting>
<scilab:image localized="true">
function rhs=fsub(x,z)
rhs=[-z(3)*z(1);0]
<emphasis role="bold">A boundary value problem with more than one solution.</emphasis>
</para>
<programlisting role="example"><![CDATA[
-// DE: y''(x)=-exp(y(x))
+// DE: y''(x) = -exp(y(x))
// BV: y(0)=0; y(1)=0
// This boundary value problem has more than one solution.
// It is demonstrated how to find two of them with the help of
// some preinformation of the solutions y(x) to build the function ystart.
-// z=[y(x);y'(x)]
-
-a=0;
-b=1;
-m=2;
-n=1;
-zeta=[a b];
-N=101;
-tol=1e-8*[1 1];
-x=linspace(a,b,N);
-
-function rhs=fsub(x,z)
- rhs=-exp(z(1));
+// z = [y(x);y'(x)]
+
+[a, b, m, n] = (0, 1, 2, 1);
+zeta = [a b];
+N = 101;
+tol = 1e-8*[1 1];
+x = linspace(a,b,N);
+
+function rhs = fsub(x,z)
+ rhs = -exp(z(1))
endfunction
-function g=gsub(i,z)
- g=[z(1) z(1)]
- g=g(i)
+function g = gsub(i,z)
+ g = [z(1) z(1)]
+ g = g(i)
endfunction
-function [z,lhs]=ystart(x,z,M)
- //z=[4*x*(1-x)*M ; 4*(1-2*x)*M]
- z=[M;0]
- //lhs=[-exp(4*x*(1-x)*M)]
- lhs=0
+function [z,lhs] = ystart(x,z,M)
+ // z = [4*x*(1-x)*M ; 4*(1-2*x)*M]
+ z = [M;0]
+ //lhs = [-exp(4*x*(1-x)*M)]
+ lhs = 0
endfunction
-for M=[1 4]
+for M = [1 4]
if M==1
- z=bvodeS(x,m,n,a,b,fsub,gsub,zeta,ystart=list(ystart,M),tol=tol);
+ z = bvodeS(x,m,n,a,b,fsub,gsub,zeta,ystart=list(ystart,M),tol=tol);
else
- z1=bvodeS(x,m,n,a,b,fsub,gsub,zeta,ystart=list(ystart,M),tol=tol);
+ z1 = bvodeS(x,m,n,a,b,fsub,gsub,zeta,ystart=list(ystart,M),tol=tol);
end
end
// Integrating the ode yield e.g. the two solutions yex and yex1.
-function y=f(c)
- y=c.*(1-tanh(sqrt(c)/4).^2)-2;
+function y = f(c)
+ y = c.*(1-tanh(sqrt(c)/4).^2)-2;
endfunction
-c=fsolve(2,f);
+c = fsolve(2,f);
-function y=yex(x,c)
- y=log(c/2*(1-tanh(sqrt(c)*(1/4-x/2)).^2))
+function y = yex(x,c)
+ y = log(c/2*(1-tanh(sqrt(c)*(1/4-x/2)).^2))
endfunction
-function y=f1(c1), y=2*c1^2+tanh(1/4/c1)^2-1;endfunction
-c1=fsolve(0.1,f1);
+function y = f1(c1), y=2*c1^2+tanh(1/4/c1)^2-1; endfunction
+c1 = fsolve(0.1,f1);
-function y=yex1(x,c1)
- y=log((1-tanh((2*x-1)/4/c1).^2)/2/c1/c1)
+function y = yex1(x,c1)
+ y = log((1-tanh((2*x-1)/4/c1).^2)/2/c1/c1)
endfunction
disp('norm(yex(x)-z(1,:))= ', norm(z(1,:)-yex(x)))
subplot(2,1,2)
plot2d(x,z1(1,:),style=[5])
xtitle(' ','x',' ')
- ]]></programlisting>
+ ]]></programlisting>
<scilab:image localized="true">
a=0;
b=1;
</para>
<programlisting role="example"><![CDATA[
// DE y'''(x)=1
-// z=[y(x);y'(x);y''(x)]
+// z = [y(x);y'(x);y''(x)]
// BV: y(-1)=2 y(1)=2
// Side condition: y(0)=1
-a=-1;b=1;c=0;
+a = -1; b = 1; c = 0;
// The side condition point c must be included in the array fixpnt.
-n=1;
-m=[3];
+n = 1;
+m = 3;
-function rhs=fsub(x,z)
- rhs=1
+function rhs = fsub(x,z)
+ rhs = 1
endfunction
-function g=gsub(i,z)
- g=[z(1)-2 z(1)-1 z(1)-2]
- g=g(i)
+function g = gsub(i,z)
+ g = [z(1)-2 z(1)-1 z(1)-2]
+ g = g(i)
endfunction
-N=10;
-zeta=[a c b];
-x=linspace(a,b,N);
+N = 10;
+zeta = [a c b];
+x = linspace(a,b,N);
-z=bvodeS(x,m,n,a,b,fsub,gsub,zeta,fixpnt=c);
+z = bvodeS(x,m,n,a,b,fsub,gsub,zeta,fixpnt=c);
-function y=yex(x)
-y=x.^3/6+x.^2-x./6+1
+function y = yex(x)
+ y = x.^3/6 + x.^2 - x/6 + 1
endfunction
disp(norm(yex(x)-z(1,:)),'norm(yex(x)-z(1,:))= ')
- ]]></programlisting>
+ ]]></programlisting>
</listitem>
<listitem>
<para>
// the differential equation at 1.0 and y, z_3=1, z_4=1/(2*y-1)
// Ref: http://arxiv.org/pdf/hep-th/0407005
-y= 1.9d0;
-eigens=zeros(3,40); // To store the results
+y = 1.9d0;
+eigens = zeros(3,40); // To store the results
// General setup for bvode
// f(x)=2*x-1 and constant c_2, c_3, so dmval=0. Notice that the z-vector
// has mstar = 4 components, while dmval has ncomp = 3 components.
-function [z,dmval]=guess(x)
- z=[2*x-1, 2., 1., 1/(2*y-1)]
- dmval=[0,0,0]
+function [z,dmval] = guess(x)
+ z = [2*x-1, 2., 1., 1/(2*y-1)]
+ dmval = [0,0,0]
endfunction
// First execution has ipar(9)=1 and uses the guess
// Subsequent executions have ipar(9)=3 and use continuation. This is
// run in tight closed loop to not disturb the stack
-for i=1:40
-v=valv(i);
-sol=bvode(res,ncomp,m,aleft,aright,zeta,ipar,ltol,tol,fixpnt,...
- fsub,dfsub,gsub,dgsub,guess);
-eigens(:,i)=[v;sol(3,101);sol(4,101)]; // c_2 and c_3 are constant!
-ipar(9)=3;
+for i = 1:40
+ v = valv(i);
+ sol = bvode(res,ncomp,m,aleft,aright,zeta,ipar,ltol,tol,fixpnt,...
+ fsub,dfsub,gsub,dgsub,guess);
+ eigens(:,i) = [v;sol(3,101);sol(4,101)]; // c_2 and c_3 are constant!
+ ipar(9) = 3;
end
// To see the evolution of the eigenvalues with v, disp(eigens)
// This is markedly different from the case at small v.
// The continuation procedure allows to explore these exponential behaviours
// without skipping to other eigenstates.
- ]]></programlisting>
+ ]]></programlisting>
</listitem>
</itemizedlist>
</refsection>
- <refsection role="see also">
- <title>See also</title>
- <simplelist type="inline">
- <member>
- <link linkend="link">link</link>
- </member>
- <member>
- <link linkend="external">external</link>
- </member>
- <member>
- <link linkend="ode">ode</link>
- </member>
- <member>
- <link linkend="dassl">dassl</link>
- </member>
- </simplelist>
- </refsection>
<refsection>
<title>Used Functions</title>
<para>This function is based on the Fortran routine
</listitem>
</orderedlist>
</refsection>
+ <refsection role="see also">
+ <title>See also</title>
+ <simplelist type="inline">
+ <member>
+ <link linkend="link">link</link>
+ </member>
+ <member>
+ <link linkend="external">external</link>
+ </member>
+ <member>
+ <link linkend="ode">ode</link>
+ </member>
+ <member>
+ <link linkend="dassl">dassl</link>
+ </member>
+ </simplelist>
+ </refsection>
</refentry>
<para>This form of external is used to pass parameters to the
function. It must be as follows:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(res,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>res</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r=res(t,y,ydot,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>
<literal>res</literal> still returns the residual value
as a function of <literal>(t,x,xdot,x1,x2,...)</literal>, and
<para>This form of external is used to pass parameters to the
function. It must be as follows:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(jac,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>jac</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=jac(t,x,xdot,p1,p2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = jac(t,x,xdot,p1,p2,...)
+]]></screen>
<para>
<literal>jac</literal> still returns
<literal>dg/dx+cj*dg/dxdot</literal> as a function of
This form of <link linkend="external">external</link> is used to pass parameters to the
function. It must be as follows:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(surface,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>surface</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=surface(t,x,p1,p2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = surface(t,x,p1,p2,...)
+]]></screen>
</listitem>
</varlistentry>
<varlistentry>
Option <literal>"root"</literal> calls the <link linkend="dasrt">dasrt</link> routine,
and <literal>"root2"</literal> calls <link linkend="dasrt">daskr</link>.
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
g(t, x, xdot) = 0
x(t0) = x0 and xdot(t0) = xdot0
- ]]></programlisting>
+]]></screen>
<para>
If <literal>xdot0</literal> is not given in the <emphasis>initial</emphasis>
argument, the <literal>dae</literal> function tries to compute it solving
%DAEOPTIONS = list([], 1, [], [], [], 0, 0); // Ask dae mesh points to be returned
y = dae([x0, xd0], 0, 4d10, chemres); // Without jacobian
y = dae([x0, xd0], 0, 4d10, chemres, chemjac); // With jacobian
- ]]></programlisting>
+ ]]></programlisting>
<para>
Example #2: dasrt ("root")
</para>
plot(yy(1, 1), yy(2, 1), 'r+')
xstring(yy(1, 1)+0.1, yy(2, 1), string(nn(1)));
cd(previous_dir);
- ]]></programlisting>
+ ]]></programlisting>
<scilab:image localized="false"><![CDATA[
code = ['#include <math.h>'
'void res22(double *t, double *y, double *yd, double *res, int *ires, double *rpar, int *ipar)'
plot(yy(1, 1), yy(2, 1), 'r+')
xstring(yy(1, 1)+0.1, yy(2, 1), string(nn(1)));
cd(previous_dir);
- ]]></programlisting>
+ ]]></programlisting>
<scilab:image><![CDATA[
code = ['#include <math.h>'
'void res22(double *t, double *y, double *yd, double *res, int *ires, double *rpar, int *ipar)'
<para>
The variable <literal>%DAEOPTIONS</literal> is a <link linkend="list">list</link> with the following elements:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(tstop,imode,band,maxstep,stepin,nonneg,isest)
- ]]></programlisting>
+]]></screen>
<para>The default value is:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list([],0,[],[],[],0,0)
- ]]></programlisting>
+]]></screen>
<para>The meaning of the elements is described below.</para>
<variablelist>
<varlistentry>
<para>This form allows to pass parameters other than t, y, ydot to the function.
It must be as follows:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(res, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>res</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r = res(t, y, ydot, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>
<literal>res</literal> still returns
<literal>r = g(t, y, ydot)</literal> as a function of
linked with Scilab.
</para>
<para>In C, the syntax must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void res(double *t, double y[], double yd[], double r[],
int *ires, double rpar[], int ipar[])
- ]]></programlisting>
+]]></screen>
<para>In Fortran, it must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine res(t, y, yd, r, ires, rpar, ipar)
double precision t, y(*), yd(*),r(*),rpar(*)
integer ires, ipar(*)
- ]]></programlisting>
+]]></screen>
<para>
The <literal>rpar</literal> and <literal>ipar</literal>
arrays must be present but cannot be used.
<listitem>
<para>A list.</para>
<para>It must be as follows</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(jac, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>jac</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=jac(t, y, ydot, cj, x1, x2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = jac(t, y, ydot, cj, x1, x2,...)
+]]></screen>
<para>
<literal>jac</literal> still returns
<literal>dg/dy + cj*dg/dydot</literal> as a function of
<para>It must refer to the name of a C function or a Fortran subroutine linked with Scilab
</para>
<para>In C, the syntax must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void jac(double *t, double y[], double yd[], double pd[],
double *cj, double rpar[], int ipar[])
- ]]></programlisting>
+]]></screen>
<para>In Fortran, it must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine jac(t, y, yd, pd, cj, rpar, ipar)
double precision t, y(*), yd(*), pd(*), cj, rpar(*)
integer ipar(*)
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<listitem>
<para>A list.</para>
<para>It must be as follows</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(surf, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>surf</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r = surf(t, y, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>A character string.</para>
<para>It must refer to the name of a C function or a Fortran subroutine linked with Scilab.
</para>
<para>In C, the syntax must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void surf(int *ny, double *t, double y[], int *ng, double gout[])
- ]]></programlisting>
+]]></screen>
<para>In Fortran, it must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine surf(ny, t, y, ng, gout)
double precision t, y(*), gout(*)
integer ny, ng
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<listitem>
<para>A list.</para>
<para>It must be as follows:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(psol, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>
where the syntax of <literal>psol</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
psol(wp, iwp, b, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>
<literal>psol</literal> still returns the solution in <literal>r</literal>.
</para>
<para>It must refer to the name of a C function or a Fortran subroutine linked with Scilab
</para>
<para>In C, the syntax must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void psol (int*neq, double*t, double*y, double*ydot, double*savr,
double*wk, double*cj, double*wght, double*wp, int*iwp, double*b, double*eplin, int*ier, double*rpar, int*ipar)
- ]]></programlisting>
+]]></screen>
where the arrays <literal>wp</literal> and <literal>iwp</literal> contain matrix elements of LU-factored preconditioner
<literal>P</literal>, <literal>wp</literal> being the values and
<literal>iwp</literal> the pivots used in the factorization.
<para>In Fortran, it must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine psol (neq, t, y, ydot, savr, wk, cj, wght,
wp, iwp, b, eplin, ier, rpar, ipar)
double precision t,y(*), ydot(*), savr(*), wk(*), cj, wght(*), wp(*),
b(*), eplin, rpar(*)
integer neq, iwp(*), ier, ipar(*)
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<listitem>
<para>A list.</para>
<para>It must be as follows</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(pjac, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>
where the syntax of <literal>pjac</literal> is
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
pjac(neq, t, y, ydot, h, cj, rewt, savr, x1, x2,...)
- ]]></programlisting>
+]]></screen>
<para>
<literal>pjac</literal> still returns factorized
<literal>dg/dy + cj*dg/dydot</literal> as a function of
<para>It must refer to the name of a C function or a Fortran subroutine linked with Scilab
</para>
<para>In C, the syntax must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void pjac (double*res, int*ires, int*neq, double*t, double*y, double*ydot, double*rewt, double*savr,
double*wk, double*h, double*cj, double*wp, int*iwp, int*ier, double*rpar, int*ipar)
- ]]></programlisting>
+]]></screen>
<para>In Fortran, it must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine pjac (res, ires, neq, t, y, ydot, rewt, savr,
wk, h, cj, wp, iwp, ier, rpar, ipar)
double precision res(*), t, y(*), ydot(*), rewt(*), savr(*),
wk(*), h, cj, wp(*), rpar(*)
integer ires, neq, iwp(*), ier, ipar(*)
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<refsection>
<title>Description</title>
<para>Solution of the implicit differential equation:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
g(t, y, ydot) = 0
y(t0) = y0 and ydot(t0) = ydot0
- ]]></programlisting>
+]]></screen>
<para>Returns the surface crossing instants and the number of the surface
reached in <literal>nn</literal>.
</para>
nn
// Should return nn = [2.4698972 2]
- ]]></programlisting>
+ ]]></programlisting>
</refsection>
<refsection role="see also">
<title>See also</title>
</refnamediv>
<refsynopsisdiv>
<title>Syntax</title>
- <synopsis>[r,nn,[,hd]]=dasrt(x0,t0,t [,atol,[rtol]],res [,jac],ng, surf [,info] [,hd])</synopsis>
+ <synopsis>[r,nn,[,hd]] = dasrt(x0,t0,t [,atol,[rtol]],res [,jac],ng, surf [,info] [,hd])</synopsis>
</refsynopsisdiv>
<refsection>
<title>Arguments</title>
<para>This form allows to pass parameters other than t,y,ydot to
the function. It must be as follows:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(res,x1,x2,...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>res</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=res(t,y,ydot,x1,x2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = res(t,y,ydot,x1,x2,...)
+]]></screen>
<para>
<literal>res</literal> still returns
<literal>r=g(t,y,ydot)</literal> as a function of
linked with Scilab.
</para>
<para>In C The syntax must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void res(double *t, double y[], double yd[], double r[],
int *ires, double rpar[], int ipar[])
- ]]></programlisting>
+]]></screen>
<para>In Fortran it must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine res(t,y,yd,r,ires,rpar,ipar)
double precision t, y(*),yd(*),r(*),rpar(*)
integer ires,ipar(*)
- ]]></programlisting>
+]]></screen>
<para>
The <literal>rpar</literal> and <literal>ipar</literal> arrays must be present but cannot be
used.
<listitem>
<para>A list.</para>
<para>It must be as follows</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(jac,x1,x2,...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>jac</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=jac(t,y,ydot,cj,x1,x2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = jac(t,y,ydot,cj,x1,x2,...)
+]]></screen>
<para>
<literal>jac</literal> still returns
<literal>dg/dy + cj*dg/dydot</literal> as a function of
with scilab
</para>
<para>In C The syntax must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void jac(double *t, double y[], double yd[], double pd[],
double *cj, double rpar[], int ipar[])
- ]]></programlisting>
+]]></screen>
<para>In Fortran it must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine jac(t,y,yd,pd,cj,rpar,ipar)
double precision t, y(*),yd(*),pd(*),cj,rpar(*)
integer ipar(*)
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<listitem>
<para>A list.</para>
<para>It must be as follows</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(surf,x1,x2,...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>surf</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=surf(t,y,x1,x2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = surf(t,y,x1,x2,...)
+]]></screen>
</listitem>
<listitem>
<para>A character string.</para>
with scilab.
</para>
<para>In C the syntax must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void surf(int *ny, double *t, double y[], int *ng, double gout[])
- ]]></programlisting>
+]]></screen>
<para>In Fortran it must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine surf(ny,t,y,ng,gout)
double precision t, y(*),gout(*)
integer ny,ng
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<refsection>
<title>Description</title>
<para>Solution of the implicit differential equation.</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-g(t,y,ydot)=0
-y(t0)=y0 and ydot(t0)=ydot0
- ]]></programlisting>
+ <screen><![CDATA[
+g(t,y,ydot) = 0
+y(t0) = y0 and ydot(t0) = ydot0
+]]></screen>
<para>Returns the surface crossing instants and the number of the surface
reached in <literal>nn</literal>.
</para>
<refsection>
<title>Examples</title>
<programlisting role="example"><![CDATA[
-//dy/dt = ((2*log(y)+8)/t -5)*y, y(1) = 1, 1<=t<=6
-//g1 = ((2*log(y)+8)/t - 5)*y
-//g2 = log(y) - 2.2491
-y0=1;t=2:6;t0=1;y0d=3;
-atol=1.d-6;rtol=0;ng=2;
+// dy/dt = ((2*log(y)+8)/t -5)*y, y(1) = 1, 1<=t<=6
+// g1 = ((2*log(y)+8)/t - 5)*y
+// g2 = log(y) - 2.2491
+y0 = 1; t = 2:6; t0 = 1; y0d = 3;
+atol = 1.d-6; rtol = 0; ng = 2;
-deff('[delta,ires]=res1(t,y,ydot)','ires=0;delta=ydot-((2*log(y)+8)/t-5)*y')
-deff('[rts]=gr1(t,y)','rts=[((2*log(y)+8)/t-5)*y;log(y)-2.2491]')
+deff('[delta,ires] = res1(t,y,ydot)', 'ires=0; delta=ydot-((2*log(y)+8)/t-5)*y')
+deff('rts = gr1(t,y)', 'rts=[((2*log(y)+8)/t-5)*y;log(y)-2.2491]')
-[yy,nn]=dasrt([y0,y0d],t0,t,atol,rtol,res1,ng,gr1);
+[yy,nn] = dasrt([y0,y0d],t0,t,atol,rtol,res1,ng,gr1);
//(Should return nn=[2.4698972 2])
- ]]></programlisting>
+ ]]></programlisting>
</refsection>
<refsection role="see also">
<title>See also</title>
</refnamediv>
<refsynopsisdiv>
<title>Syntax</title>
- <synopsis>[r [,hd]]=dassl(x0,t0,t [,atol,[rtol]],res [,jac] [,info] [,hd])</synopsis>
+ <synopsis>[r [,hd]] = dassl(x0,t0,t [,atol,[rtol]],res [,jac] [,info] [,hd])</synopsis>
</refsynopsisdiv>
<refsection>
<title>Arguments</title>
<literal>y</literal>, <literal>ydot</literal> to
the function. It must be as follows:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(res,x1,x2,...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>res</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=res(t,y,ydot,x1,x2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = res(t,y,ydot,x1,x2,...)
+]]></screen>
<para>
<literal>res</literal> still returns
<literal>r=g(t,y,ydot)</literal> as a function of
linked with Scilab.
</para>
<para>In C the syntax must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void res(double *t, double y[], double yd[], double r[],
int *ires, double rpar[], int ipar[])
- ]]></programlisting>
+]]></screen>
<para>In Fortran it must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine res(t,y,yd,r,ires,rpar,ipar)
double precision t, y(*),yd(*),r(*),rpar(*)
integer ires,ipar(*)
- ]]></programlisting>
+]]></screen>
<para>
The <literal>rpar</literal> and <literal>ipar</literal> arrays must be present but cannot be
used.
<listitem>
<para>A list.</para>
<para>It must be as follows</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(jac,x1,x2,...)
- ]]></programlisting>
+]]></screen>
<para>where the syntax of the function
<literal>jac</literal> is now
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=jac(t,y,ydot,cj,x1,x2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = jac(t,y,ydot,cj,x1,x2,...)
+]]></screen>
<para>
<literal>jac</literal> still returns
<literal>dg/dy+cj*dg/dydot</literal> as a function of
with Scilab.
</para>
<para>In C the syntax must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void jac(double *t, double y[], double yd[], double pd[],
double *cj, double rpar[], int ipar[])
- ]]></programlisting>
+]]></screen>
<para>In Fortran it must be:</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine jac(t,y,yd,pd,cj,rpar,ipar)
double precision t, y(*),yd(*),pd(*),cj,rpar(*)
integer ipar(*)
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
The <literal>dassl</literal> function integrate the differential algebraic equation and
returns the evolution of <literal>y</literal> a given time points
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-g(t,y,ydot)=0
-y(t0)=y0 and ydot(t0)=ydot0
- ]]></programlisting>
+ <screen><![CDATA[
+g(t,y,ydot) = 0
+y(t0) = y0 and ydot(t0) = ydot0
+ ]]></screen>
</refsection>
<refsection>
<title>Examples</title>
<programlisting role="example"><![CDATA[
-function [r,ires]=chemres(t,y,yd)
+function [r,ires] = chemres(t,y,yd)
r=[-0.04*y(1)+1d4*y(2)*y(3)-yd(1)
0.04*y(1)-1d4*y(2)*y(3)-3d7*y(2)*y(2)-yd(2)
y(1)+y(2)+y(3)-1];
ires=0
endfunction
-function pd=chemjac(x,y,yd,cj)
+function pd = chemjac(x,y,yd,cj)
pd=[-0.04-cj , 1d4*y(3) , 1d4*y(2);
0.04 ,-1d4*y(3)-2*3d7*y(2)-cj ,-1d4*y(2);
1 , 1 , 1 ]
endfunction
-y0=[1;0;0];
-yd0=[-0.04;0.04;0];
-t=[1.d-5:0.02:.4,0.41:.1:4,40,400,4000,40000,4d5,4d6,4d7,4d8,4d9,4d10];
+y0 = [1;0;0];
+yd0 = [-0.04;0.04;0];
+t = [1.d-5:0.02:.4,0.41:.1:4,40,400,4000,40000,4d5,4d6,4d7,4d8,4d9,4d10];
-y=dassl([y0,yd0],0,t,chemres);
+y = dassl([y0,yd0],0,t,chemres);
-info=list([],0,[],[],[],0,0);
-info(2)=1;
-y1=dassl([y0,yd0],0,4d10,chemres,info);
-y2=dassl([y0,yd0],0,4d10,chemres,chemjac,info);
+info = list([],0,[],[],[],0,0);
+info(2) = 1;
+y1 = dassl([y0,yd0],0,4d10,chemres,info);
+y2 = dassl([y0,yd0],0,4d10,chemres,chemjac,info);
//Using extra argument for parameters
//-----------------------------------
-function [r,ires]=chemres(t,y,yd ,a,b,c)
+function [r,ires] = chemres(t,y,yd ,a,b,c)
r=[-a*y(1)+b*y(2)*y(3)-yd(1)
a*y(1)-b*y(2)*y(3)-c*y(2)*y(2)-yd(2)
y(1)+y(2)+y(3)-1];
- ires=0
+ ires = 0
endfunction
-function pd=chemjac(x,y,yd,cj, a,b,c)
+function pd = chemjac(x,y,yd,cj, a,b,c)
pd=[-a-cj , b*y(3) , b*y(2);
a ,-b*y(3)-2*c*y(2)-cj ,-b*y(2);
1 , 1 , 1 ]
endfunction
-y3=dassl([y0,yd0],0,t,list(chemres,0.04,1d4,3d7),list(chemjac,0.04,1d4,3d7));
+y3 = dassl([y0,yd0], 0, t, list(chemres,0.04,1d4,3d7),list(chemjac,0.04,1d4,3d7));
-//using C code
-//------------
+// using C code
+// ------------
// - create the C code
cd TMPDIR
rescode=['void chemres(double *t, double y[], double yd[], double r[], int *ires, double rpar[], int ipar[])'
' pd[7] = -1.0e4*y[1];'
' pd[8] = 1.0;'
' }'];
-mputl([rescode;jaccode],'mycode.c') //create the C file
+mputl([rescode;jaccode], 'mycode.c') //create the C file
// - compile it
ilib_for_link(['chemres','chemjac'],'mycode.c',[],'c','','loader.sce');//compile
exec('loader.sce') //incremental linking
// - call dassl
-y4=dassl([y0,yd0],0,t,'chemres','chemjac');
- ]]></programlisting>
+y4 = dassl([y0,yd0], 0, t, 'chemres', 'chemjac');
+ ]]></programlisting>
</refsection>
<refsection role="see also">
<title>See also</title>
<refsynopsisdiv>
<title>Syntax</title>
<synopsis>
- y=diff(x)
- y=diff(x [,n [,dim]])
+ y = diff(x)
+ y = diff(x, n)
+ y = diff(x, n, dim)
</synopsis>
</refsynopsisdiv>
<refsection>
<title>Examples</title>
<programlisting role="example">
<![CDATA[
-v=(1:8)^3;
+v = (1:8)^3;
diff(v)
diff(v,3)
diff(A,3,2)
-//approximate differentiation
-step=0.001
-t=0:step:10;
-y=sin(t);
-dy=diff(sin(t))/step; //approximate differentiation of sine function
-norm(dy-cos(t(1:$-1)),%inf)
+// approximate differentiation
+step = 0.001
+t = 0:step:10;
+y = sin(t);
+dy = diff(sin(t)) / step; //approximate differentiation of sine function
+norm(dy-cos(t(1:$-1)), %inf)
]]>
</programlisting>
</refsection>
</refnamediv>
<refsynopsisdiv>
<title>Syntax</title>
- <synopsis>y=impl([type],y0,ydot0,t0,t [,atol, [rtol]],res,adda [,jac])</synopsis>
+ <synopsis>y = impl([type,] y0, ydot0, t0, t [,atol [,rtol]], res, adda [,jac])</synopsis>
</refsynopsisdiv>
<refsection>
<title>Arguments</title>
If <literal>res</literal> is a function, its syntax must be as
follows:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r = res(t,y,ydot)
- ]]></programlisting>
+]]></screen>
<para>
where <literal>t</literal> is a real scalar (time) and
<literal>y</literal> and <literal>ydot</literal> are real vector (state
If <literal>adda</literal> is a function, its syntax must be as
follows:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r = adda(t,y,p)
- ]]></programlisting>
+]]></screen>
<para>
and it must return <literal>r=A(t,y)+p</literal> where
<literal>p</literal> is a matrix to be added to
If <literal>jac</literal> is a function, its syntax must be as
follows:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
j = jac(t,y,ydot)
- ]]></programlisting>
+]]></screen>
<para>and it must return the Jacobian of
<literal>r=g(t,y)-A(t,y)*ydot</literal> with respect to
<literal>y</literal>.
<refsection>
<title>Examples</title>
<programlisting role="example"><![CDATA[
-y=impl([1;0;0],[-0.04;0.04;0],0,0.4,'resid','aplusp');
+y = impl([1;0;0], [-0.04;0.04;0], 0, 0.4, 'resid', 'aplusp');
// Using hot restart
-//[x1,w,iw]=impl([1;0;0],[-0.04;0.04;0],0,0.2,'resid','aplusp');
+//[x1,w,iw] = impl([1;0;0], [-0.04;0.04;0], 0, 0.2, 'resid', 'aplusp');
// hot start from previous call
-//[x1]=impl([1;0;0],[-0.04;0.04;0],0.2,0.4,'resid','aplusp',w,iw);
+//x1 = impl([1;0;0],[-0.04;0.04;0], 0.2, 0.4, 'resid', 'aplusp', w, iw);
//max(abs(x1-x))
]]></programlisting>
</refsection>
</refnamediv>
<refsynopsisdiv>
<title>Syntax</title>
- <synopsis>[result,err]=int3d(X,Y,Z,f [,nf[,params]])</synopsis>
+ <synopsis>
+ [result, err] = int3d(X, Y, Z, f)
+ [result, err] = int3d(X, Y, Z, f, nf)
+ [result, err] = int3d(X, Y, Z, f, nf, params)
+ </synopsis>
</refsynopsisdiv>
<refsection>
<title>Arguments</title>
<para>The function calculates an approximation to a given vector of
definite integrals
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
I I I (F ,F ,...,F ) dx(3)dx(2)dx(1),
1 2 numfun
- ]]></programlisting>
+]]></screen>
<para>
where the region of integration is a collection of <literal>NUMTET</literal>
tetrahedrons and where
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
F = F (X(1),X(2),X(3)), J = 1,2,...,NUMFUN.
J J
- ]]></programlisting>
+]]></screen>
<para>A globally adaptive strategy is applied in order to compute
approximations <literal>result(k)</literal> hopefully satisfying, for each
component of <literal>I</literal>, the following claim for accuracy:
<refsection>
<title>Examples</title>
<programlisting role="example"><![CDATA[
-X=[0;1;0;0];
-Y=[0;0;1;0];
-Z=[0;0;0;1];
-[RESULT,ERROR]=int3d(X,Y,Z,'int3dex')
+X = [0;1;0;0];
+Y = [0;0;1;0];
+Z = [0;0;0;1];
+[RESULT, ERROR] = int3d(X, Y, Z, 'int3dex')
// computes the integrand exp(x*x+y*y+z*z) over the
//tetrahedron (0.,0.,0.),(1.,0.,0.),(0.,1.,0.),(0.,0.,1.)
//integration over a cube -1<=x<=1;-1<=y<=1;-1<=z<=1
// bottom -top- right -left- front -rear-
-X=[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0;
+X =[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0;
-1,-1, -1,-1, 1, 1, -1,-1, -1,-1, -1,-1;
1,-1, 1,-1, 1, 1, -1,-1, 1,-1, 1,-1;
1, 1, 1, 1, 1, 1, -1,-1, 1, 1, 1, 1];
-Y=[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0;
+Y =[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0;
-1,-1, -1,-1, -1, 1, -1, 1, -1,-1, 1, 1;
-1, 1, -1, 1, 1, 1, 1, 1, -1,-1, 1, 1;
1, 1, 1, 1, -1,-1, -1,-1, -1,-1, 1, 1];
-Z=[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0;
+Z =[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0;
-1,-1, 1, 1, -1, 1, -1, 1, -1,-1, -1,-1;
-1,-1, 1, 1, -1,-1, -1,-1, -1, 1, -1, 1;
-1,-1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1];
-function v=f(xyz,numfun),v=exp(xyz'*xyz),endfunction
-[result,err]=int3d(X,Y,Z,f,1,[0,100000,1.d-5,1.d-7])
+function v = f(xyz,numfun),v=exp(xyz'*xyz), endfunction
+[result, err] = int3d(X, Y, Z, f, 1, [0,100000,1.d-5,1.d-7])
-function v=f(xyz,numfun),v=1,endfunction
-[result,err]=int3d(X,Y,Z,f,1,[0,100000,1.d-5,1.d-7])
- ]]></programlisting>
+function v = f(xyz,numfun), v=1, endfunction
+[result, err] = int3d(X, Y, Z, f, 1, [0,100000,1.d-5,1.d-7])
+ ]]></programlisting>
</refsection>
<refsection role="see also">
<title>See also</title>
x1 = -10:0.1:10;
Y = integrate(['if x==0 then 1,';
'else sin(x)/x,end'], 'x', 0, x1)
- ]]></programlisting>
+ ]]></programlisting>
</refsection>
<refsection role="see also">
<title>See also</title>
<title>Examples</title>
<programlisting role="example"><![CDATA[
// Function written in the Scilab language
-function y=f(x),y=x*sin(30*x)/sqrt(1-((x/(2*%pi))^2)),endfunction
-exact=-2.5432596188;
-I=intg(0,2*%pi,f)
+function y = f(x), y = x*sin(30*x)/sqrt(1-((x/(2*%pi))^2)), endfunction
+exact = -2.5432596188;
+I = intg(0,2*%pi,f)
abs(exact-I)
// Function with an argument written in the Scilab language
</refnamediv>
<refsynopsisdiv>
<title>Syntax</title>
- <synopsis>[y]=intl(a,b,z0,r,f)</synopsis>
+ <synopsis>y = intl(a, b, z0, r, f)</synopsis>
</refsynopsisdiv>
<refsection>
<title>Arguments</title>
<refsection>
<title>Examples</title>
<programlisting role="example"><![CDATA[
-function y=f(z)
+function y = f(z)
y = z^(3 + %pi * %i)
endfunction
-intl(1,2,1+%i, 3, f)
+intl(1, 2, 1+%i, 3, f)
]]></programlisting>
</refsection>
</refnamediv>
<refsynopsisdiv>
<title>Syntax</title>
- <synopsis>v = inttrap([x,] s)</synopsis>
+ <synopsis>
+ v = inttrap(s)
+ v = inttrap(x, s)
+ </synopsis>
</refsynopsisdiv>
<refsection>
<title>Arguments</title>
<refsection>
<title>Examples</title>
<programlisting role="example"><![CDATA[
-t=0:0.1:%pi
+t = 0:0.1:%pi
inttrap(t,sin(t))
]]></programlisting>
</refsection>
<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:id="ode_discrete" xml:lang="en">
<refnamediv>
<refname>ode_discrete</refname>
- <refpurpose>ordinary differential equation solver, discrete time
- simulation
+ <refpurpose>ordinary differential equation solver, discrete time simulation
</refpurpose>
</refnamediv>
<refsynopsisdiv>
<title>Syntax</title>
- <synopsis>y=ode("discrete",y0,k0,kvect,f)</synopsis>
+ <synopsis>y = ode("discrete", y0, k0,kvect,f)</synopsis>
</refsynopsisdiv>
<refsection>
<title>Arguments</title>
<refsection>
<title>Examples</title>
<programlisting role="example"><![CDATA[
-y1=[1;2;3]; deff("yp=a_function(k,y)","yp=A*y+B*u(k)")
-A=diag([0.2,0.5,0.9]); B=[1;1;1];u=1:10;n=5;
-y=ode("discrete",y1,1,1:n,a_function);
-y(:,2)-(A*y1+B*u(1))
+y1 = [1;2;3]; deff("yp = a_function(k,y)", "yp = A*y+B*u(k)")
+A = diag([0.2,0.5,0.9]); B = [1;1;1]; u = 1:10; n = 5;
+y = ode("discrete", y1, 1, 1:n, a_function);
+y(:,2) - (A*y1+B*u(1))
// Now y evaluates at [y3,y5,y7,y9]
-y=ode("discrete",y1,1,3:2:9,a_function)
+y = ode("discrete", y1, 1, 3:2:9, a_function)
]]></programlisting>
</refsection>
<refsection role="see also">
If <literal>g</literal> is a function the syntax should be as
follows:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
z = g(t,y)
- ]]></programlisting>
+]]></screen>
<para>
where <literal>t</literal> is a real scalar (time) and
<literal>y</literal> a real vector (state). It returns a vector of size
<programlisting role="example"><![CDATA[
// Integration of the differential equation
// dy/dt=y , y(0)=1, and finds the minimum time t such that y(t)=2
-deff("[ydot]=f(t,y)","ydot=y")
-deff("[z]=g(t,y)","z=y-2")
-y0=1;ng=1;
-[y,rd]=ode("root",y0,0,2,f,ng,g)
+deff("ydot = f(t,y)", "ydot=y")
+deff("z = g(t,y)", "z=y-2")
+y0 = 1; ng = 1;
+[y,rd] = ode("root", y0, 0, 2, f, ng, g)
-deff("[z]=g(t,y)","z=y-[2;2;33]")
-[y,rd]=ode("root",1,0,2,f,3,g)
+deff("z = g(t,y)", "z = y-[2;2;33]")
+[y,rd] = ode("root", 1, 0, 2, f, 3, g)
]]></programlisting>
</refsection>
<refsection role="see also">
<title>Examples</title>
<programlisting role="example"><![CDATA[
//Linear system with switching input
-deff('xdu=phis(t,x,u,flag)','if flag==0 then xdu=A*x+B*u; else xdu=1-u;end');
-x0=[1;1];
-A=[-1,2;-2,-1];
-B=[1;2];
-u=0;
-nu=1;
-stdel=[1,0];
-u0=0;
-t=0:0.05:10;
-xu=odedc([x0;u0],nu,stdel,0,t,phis);
-x=xu(1:2,:);
-u=xu(3,:);
-nx=2;
+deff('xdu = phis(t,x,u,flag)', 'if flag==0 then xdu=A*x+B*u; else xdu=1-u; end');
+x0 = [1;1];
+A = [-1,2;-2,-1];
+B = [1;2];
+u = 0;
+nu = 1;
+stdel = [1,0];
+u0 = 0;
+t = 0:0.05:10;
+xu = odedc([x0;u0], nu, stdel, 0, t, phis);
+x = xu(1:2,:);
+u = xu(3,:);
+nx = 2;
plot2d(t, x', [1:nx], '161')
-plot2d2('onn',t',u',[nx+1:nx+nu],'000');
+plot2d2('onn', t', u', [nx+1:nx+nu], '000');
//Fortran external (see fydot2.f):
-norm(xu-odedc([x0;u0],nu,stdel,0,t,'phis'),1)
+norm(xu - odedc([x0;u0], nu, stdel, 0, t, 'phis'), 1)
]]></programlisting>
<scilab:image>
deff('xdu=phis(t,x,u,flag)','if flag==0 then xdu=A*x+B*u; else xdu=1-u;end');
// (feedback) |
// | u=hd(t,xd)
//
-deff('xcd=f(t,xc,xd,iflag)',...
+deff('xcd = f(t,xc,xd,iflag)',...
['if iflag==0 then '
' xcd=fc(t,xc,e(t)-hd(t,xd));'
'else '
' xcd=fd(xd,hc(t,xc));'
'end']);
-A=[-10,2,3;4,-10,6;7,8,-10];
-B=[1;1;1];
-C=[1,1,1];
-Ad=[1/2,1;0,1/20];
-Bd=[1;1];
-Cd=[1,1];
-deff('st=e(t)','st=sin(3*t)')
-deff('xdot=fc(t,x,u)','xdot=A*x+B*u')
-deff('y=hc(t,x)','y=C*x')
-deff('xp=fd(x,y)','xp=Ad*x + Bd*y')
-deff('u=hd(t,x)','u=Cd*x')
-h=0.1;t0=0;t=0:0.1:2;
-x0c=[0;0;0];
-x0d=[0;0];
-nd=2;
-xcd=odedc([x0c;x0d],nd,h,t0,t,f);
-norm(xcd-odedc([x0c;x0d],nd,h,t0,t,'fcd1')) // Fast calculation (see fydot2.f)
-plot2d([t',t',t'],xcd(1:3,:)');
+A = [-10,2,3;4,-10,6;7,8,-10];
+B = [1;1;1];
+C = [1,1,1];
+Ad = [1/2,1;0,1/20];
+Bd = [1;1];
+Cd = [1,1];
+deff('st = e(t)', 'st = sin(3*t)')
+deff('xdot = fc(t,x,u)', 'xdot = A*x+B*u')
+deff('y = hc(t,x)', 'y = C*x')
+deff('xp = fd(x,y)', 'xp = Ad*x + Bd*y')
+deff('u = hd(t,x)', 'u = Cd*x')
+h = 0.1; t0=0; t=0:0.1:2;
+x0c = [0;0;0];
+x0d = [0;0];
+nd = 2;
+xcd = odedc([x0c;x0d], nd, h, t0, t, f);
+norm(xcd - odedc([x0c;x0d], nd, h, t0, t, 'fcd1')) // Fast calculation (see fydot2.f)
+plot2d([t',t',t'], xcd(1:3,:)');
scf(2);
-plot2d2("gnn",[t',t'],xcd(4:5,:)');
+plot2d2("gnn", [t',t'], xcd(4:5,:)');
scf(0);
]]></programlisting>
<scilab:image>
condition <literal>y(0)=0</literal> claiming the solution be stored at each mesh value.
</para>
<programlisting role="example"><![CDATA[
-function ydot=f(t,y)
- ydot=y^2-y*sin(t)+cos(t)
+function ydot = f(t,y)
+ ydot = y^2 - y*sin(t) + cos(t)
endfunction
-%ODEOPTIONS=[2,0,0,%inf,0,2,500,12,5,0,-1,-1];
-y=ode(0,0,%pi,f);
+%ODEOPTIONS = [2,0,0,%inf,0,2,500,12,5,0,-1,-1];
+y = ode(0,0,%pi,f);
plot(y(1,:),y(2,:))
clear %ODEOPTIONS
]]></programlisting>
Cette forme d'<link linkend="external">external</link> sert à passer des paramètres à la fonction.
Elle doit se présenter comme suit:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(res,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>où la séquence d'appel de la fonction
<literal>res</literal> est
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r=res(t,y,ydot,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>
<literal>res</literal> retourne toujours le résidu comme fonction de
<literal>(t,x,xdot,x1,x2,...)</literal>, et
Cette forme d'<link linkend="external">external</link> est utilisée pour passer des paramètres à la fonction.
Elle doit se présenter comme suit :
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(jac,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>où la séquence d'appel de la fonction
<literal>jac</literal> est
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=jac(t,x,xdot,p1,p2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = jac(t,x,xdot,p1,p2,...)
+]]></screen>
<para>
<literal>jac</literal> retourne toujours
<literal>dg/dx+cj*dg/dxdot</literal> comme fonction de
passer des paramètres à la fonction.
Elle doit se présenter comme suit :
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(surface,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>où la qéquence d'appel de la fonction
<literal>surface</literal> est
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=surface(t,x,p1,p2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = surface(t,x,p1,p2,...)
+]]></screen>
</listitem>
</varlistentry>
<varlistentry>
L'option <literal>"root"</literal> appelle la routine <link linkend="dasrt">dasrt</link>,
et <literal>"root2"</literal> appelle <link linkend="dasrt">daskr</link>.
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
g(t, x, xdot) = 0
x(t0) = x0 and xdot(t0) = xdot0
- ]]></programlisting>
+]]></screen>
<para>
Si <literal>xdot0</literal> n'est pas donné en paramètre initial,
<literal>dae</literal> tente de la calculer en résolvant
%DAEOPTIONS = list([], 1, [], [], [], 0, 0); // Demande à dae les points à retourner
y = dae([x0, xd0], 0, 4d10, chemres); // Sans jacobian
y = dae([x0, xd0], 0, 4d10, chemres, chemjac); // Avec jacobien
- ]]></programlisting>
+ ]]></programlisting>
<para>
Exemple #2: dasrt ("root")
</para>
plot(yy(1, 1), yy(2, 1), 'r+')
xstring(yy(1, 1)+0.1, yy(2, 1), string(nn(1)));
cd(previous_dir);
- ]]></programlisting>
+ ]]></programlisting>
<scilab:image localized="false"><![CDATA[
code = ['#include <math.h>'
'void res22(double *t, double *y, double *yd, double *res, int *ires, double *rpar, int *ipar)'
plot(yy(1, 1), yy(2, 1), 'r+')
xstring(yy(1, 1)+0.1, yy(2, 1), string(nn(1)));
cd(previous_dir);
- ]]></programlisting>
+ ]]></programlisting>
<scilab:image><![CDATA[
code = ['#include <math.h>'
'void res22(double *t, double *y, double *yd, double *res, int *ires, double *rpar, int *ipar)'
<para>Cette forme permet de passer des paramètres autres que t, y, ydot à la fonction.
Elle doit se présenter comme suit :
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(res, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>où la séquence d'appel de la fonction
<literal>res</literal> est maintenant
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r = res(t, y, ydot, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>
<literal>res</literal> retourne toujours
<literal>r = g(t, y, ydot)</literal> comme fonction de
reliée à Scilab.
</para>
<para>En C, la séquence d'appel doit être :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void res(double *t, double y[], double yd[], double r[],
int *ires, double rpar[], int ipar[])
- ]]></programlisting>
+]]></screen>
<para>En Fortran, elle doit être :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine res(t, y, yd, r, ires, rpar, ipar)
double precision t, y(*), yd(*),r(*),rpar(*)
integer ires, ipar(*)
- ]]></programlisting>
+]]></screen>
<para>
Les tableaux <literal>rpar</literal> et <literal>ipar</literal>
doivent être présents mais ne peuvent pas être utilisés.
<listitem>
<para>Une liste.</para>
<para>Elle doit se présenter comme suit :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(jac, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>où la séquence d'appel de la fonction
<literal>jac</literal> est désormais
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=jac(t, y, ydot, cj, x1, x2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = jac(t, y, ydot, cj, x1, x2,...)
+]]></screen>
<para>
<literal>jac</literal> retourne toujours
<literal>dg/dy + cj*dg/dydot</literal> comme fonction de
<para>Elle doit se référer au nom d'une fonction C ou une routine Fortran reliée à Scilab.
</para>
<para>En C, la séquence d'appel doit être :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void jac(double *t, double y[], double yd[], double pd[],
double *cj, double rpar[], int ipar[])
- ]]></programlisting>
+]]></screen>
<para>En Fortran, elle doit être :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine jac(t, y, yd, pd, cj, rpar, ipar)
double precision t, y(*), yd(*), pd(*), cj, rpar(*)
integer ipar(*)
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<listitem>
<para>Une liste.</para>
<para>Elle doit se présenter comme suit :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(surf, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>où la séquence d'appel de la fonction
<literal>surf</literal> est maintenant
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r = surf(t, y, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
</listitem>
<listitem>
<para>Une chaîne de caractères.</para>
<para>Elle doit se référer au nom d'une fonction C ou une routine Fortran reliée à Scilab.
</para>
<para>En C, la séquence d'appel doit être :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void surf(int *ny, double *t, double y[], int *ng, double gout[])
- ]]></programlisting>
+]]></screen>
<para>En Fortran, elle doit être :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine surf(ny, t, y, ng, gout)
double precision t, y(*), gout(*)
integer ny, ng
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<listitem>
<para>Une liste.</para>
<para>Elle doit se présenter comme suit :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(psol, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>
où la séquence d'appel de <literal>psol</literal> est
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
psol(wp, iwp, b, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>
<literal>psol</literal> retourne toujours la solution dans <literal>r</literal>.
</para>
<para>Elle doit se référer au nom d'une fonction C ou une routine Fortran reliée à Scilab
</para>
<para>En C, la séquence d'appel doit être :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void psol (int*neq, double*t, double*y, double*ydot, double*savr,
double*wk, double*cj, double*wght, double*wp, int*iwp, double*b, double*eplin, int*ier, double*rpar, int*ipar)
- ]]></programlisting>
+]]></screen>
où les vecteurs <literal>wp</literal> et <literal>iwp</literal> contiennent le préconditionneur LU-factorisé
<literal>P</literal>, <literal>wp</literal> représentant les valeurs et
<literal>iwp</literal> les pivots utilisés dans la factorisation.
<para>En Fortran, elle doit être :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine psol (neq, t, y, ydot, savr, wk, cj, wght,
wp, iwp, b, eplin, ier, rpar, ipar)
double precision t,y(*), ydot(*), savr(*), wk(*), cj, wght(*), wp(*),
b(*), eplin, rpar(*)
integer neq, iwp(*), ier, ipar(*)
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<listitem>
<para>Une liste.</para>
<para>Elle doit se présenter comme suit :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(pjac, x1, x2, ...)
- ]]></programlisting>
+]]></screen>
<para>
où la séquence d'appel de <literal>pjac</literal> est
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
pjac(neq, t, y, ydot, h, cj, rewt, savr, x1, x2,...)
- ]]></programlisting>
+]]></screen>
<para>
<literal>pjac</literal> retourne toujours
<literal>dg/dy + cj*dg/dydot</literal> comme fonction de
<para>Elle doit se référer au nom d'une fonction C ou une routine Fortran reliée à Scilab
</para>
<para>En C, la séquence d'appel doit être :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void pjac (double*res, int*ires, int*neq, double*t, double*y, double*ydot, double*rewt, double*savr,
double*wk, double*h, double*cj, double*wp, int*iwp, int*ier, double*rpar, int*ipar)
- ]]></programlisting>
+]]></screen>
<para>En Fortran, elle doit être :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine pjac (res, ires, neq, t, y, ydot, rewt, savr,
wk, h, cj, wp, iwp, ier, rpar, ipar)
double precision res(*), t, y(*), ydot(*), rewt(*), savr(*),
wk(*), h, cj, wp(*), rpar(*)
integer ires, neq, iwp(*), ier, ipar(*)
- ]]></programlisting>
+]]></screen>
</listitem>
</itemizedlist>
</listitem>
<refsection>
<title>Description</title>
<para>Solution de l'équation différentielle implicite :</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
g(t, y, ydot) = 0
y(t0) = y0 et ydot(t0) = ydot0
- ]]></programlisting>
+]]></screen>
<para>Retourne les temps de traversée de surface et le nombre de surfaces traversées
dans <literal>nn</literal>.
</para>
nn
// Retourne nn = [2.4698972 2]
- ]]></programlisting>
+ ]]></programlisting>
</refsection>
<refsection role="see also">
<title>Voir aussi</title>
// integration de la fonction
I = intg(0, 2*%pi, 'cfun')
abs(exact - I)
- ]]></programlisting>
+ ]]></programlisting>
</refsection>
<refsection role="see also">
<title>Voir aussi</title>
</refnamediv>
<refsynopsisdiv>
<title>Séquence d'appel</title>
- <synopsis>v = inttrap([x,] s)</synopsis>
+ <synopsis>
+ v = inttrap(s)
+ v = inttrap(x, s)
+ </synopsis>
</refsynopsisdiv>
<refsection>
<title>Paramètres</title>
<refsection>
<title>Exemples</title>
<programlisting role="example"><![CDATA[
-t=0:0.1:%pi
+t = 0:0.1:%pi
inttrap(t,sin(t))
]]></programlisting>
</refsection>
</refnamediv>
<refsynopsisdiv>
<title>Syntaxe</title>
- <synopsis>y=ode(y0,t0,t,f)
- [y,w,iw]=ode([type],y0,t0,t [,rtol [,atol]],f [,jac] [,w,iw])
- [y,rd,w,iw]=ode("root",y0,t0,t [,rtol [,atol]],f [,jac],ng,g [,w,iw])
- y=ode("discrete",y0,k0,kvect,f)
+ <synopsis>y = ode(y0,t0,t,f)
+ [y,w,iw] = ode([type],y0,t0,t [,rtol [,atol]],f [,jac] [,w,iw])
+ [y,rd,w,iw] = ode("root",y0,t0,t [,rtol [,atol]],f [,jac],ng,g [,w,iw])
+ y = ode("discrete",y0,k0,kvect,f)
</synopsis>
</refsynopsisdiv>
<refsection role="parameters">
<para>Si c'est une subroutine Fortran, sa liste d'appel doit
être
</para>
- <programlisting role=""><![CDATA[
+ <screen><![CDATA[
subroutine fex(n,t,y,ydot)
integer n
double precision t,y(*),ydot(*)
- ]]></programlisting>
+]]></screen>
<para>Si c'est une fonction C son prototype doit être:</para>
- <programlisting role=""><![CDATA[
+ <screen><![CDATA[
void fex(int *n,double *t,double *y,double *ydot)
- ]]></programlisting>
+]]></screen>
<para>
Cet external peut être compilé par l'utilitaire <link linkend="ilib_for_link">ilib_for_link</link> et chargé
dynamiquement par la fonction <link linkend="link">link</link>.
<para>En Fortran, Cette routine doit avoir la liste d'appel suivante
:
</para>
- <programlisting role=""><![CDATA[
+ <screen><![CDATA[
subroutine fex(n,t,y,ml,mu,J,nrpd)
integer n,ml,mu,nrpd
double precision t,y(*),J(*)
- ]]></programlisting>
+]]></screen>
<para>Si c'est une fonction C son prototype doit être:</para>
- <programlisting role=""><![CDATA[
+ <screen><![CDATA[
void fex(int *n,double *t,double *y,int *ml,int *mu,double *J,int *nrpd,)
- ]]></programlisting>
+]]></screen>
<para>Dans la plupart des cas il n'est pas nécessaire d'utiliser
<literal>ml</literal>, <literal>mu</literal> et
<literal>nrpd</literal>, qui sont relatifs à la possibilité de
関数にパラメータを指定する際に使用されます.
以下のような形式とします:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(res,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>
ただし,ここで関数<literal>res</literal>の呼び出し手順は以下のようになります
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=res(t,y,ydot,p1,p2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = res(t,y,ydot,p1,p2,...)
+]]></screen>
<para>
この場合も<literal>res</literal> は
<literal>(t,x,xdot,x1,x2,...)</literal>の関数として残差の値を返し,
関数にパラメータを指定する際に使用されます.
以下のような形式とします:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(jac,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>
ただしこの場合の関数<literal>jac</literal>の呼び出し手順は
以下となります
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=jac(t,x,xdot,p1,p2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = jac(t,x,xdot,p1,p2,...)
+]]></screen>
<para>
この場合でも<literal>jac</literal> は,
<literal>(t,x,xdot,cj,p1,p2,...)</literal>の関数として
関数にパラメータを指定する際に使用されます.
以下のような形式とします:
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(surface,p1,p2,...)
- ]]></programlisting>
+]]></screen>
<para>
ただしこの場合の関数<literal>surface</literal>の呼び出し手順は
以下となります
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-r=surface(t,x,p1,p2,...)
- ]]></programlisting>
+ <screen><![CDATA[
+r = surface(t,x,p1,p2,...)
+]]></screen>
</listitem>
</varlistentry>
<varlistentry>
<link linkend="dasrt">dasrt</link> ルーチンをコールし,
<literal>"root2"</literal> は <link linkend="dasrt">daskr</link>をコールします.
</para>
- <programlisting role="no-scilab-exec"><![CDATA[
-g(t,x,xdot)=0
-x(t0)=x0 および xdot(t0)=xdot0
- ]]></programlisting>
+ <screen><![CDATA[
+g(t,x,xdot) = 0
+x(t0) = x0 および xdot(t0)=xdot0
+]]></screen>
<para>
<emphasis>initial</emphasis>引数に<literal>xdot0</literal> が指定されない場合,
<literal>dae</literal>関数は <literal>g(t,x0,xdot0)=0</literal>
1 , 1 , 1 ]
endfunction
-x0=[1; 0; 0];
-xd0=[-0.04; 0.04; 0];
-t=[1.d-5:0.02:.4, 0.41:.1:4, 40, 400, 4000, 40000, 4d5, 4d6, 4d7, 4d8, 4d9, 4d10];
+x0 = [1; 0; 0];
+xd0 = [-0.04; 0.04; 0];
+t = [1.d-5:0.02:.4, 0.41:.1:4, 40, 400, 4000, 40000, 4d5, 4d6, 4d7, 4d8, 4d9, 4d10];
-y=dae([x0,xd0],0,t,chemres);// 指定した観測時刻での点を返す
+y = dae([x0,xd0],0,t,chemres);// 指定した観測時刻での点を返す
-%DAEOPTIONS=list([],1,[],[],[],0,0); // dae メッシュ点を返すかどうかを確認
-y=dae([x0,xd0],0,4d10,chemres); // ヤコビアン指定なし
-y=dae([x0,xd0],0,4d10,chemres,chemjac); // ヤコビアン指定あり
- ]]></programlisting>
+%DAEOPTIONS = list([],1,[],[],[],0,0); // dae メッシュ点を返すかどうかを確認
+y = dae([x0,xd0],0,4d10,chemres); // ヤコビアン指定なし
+y = dae([x0,xd0],0,4d10,chemres,chemjac); // ヤコビアン指定あり
+ ]]></programlisting>
<para>
例 #2: dasrt ("root")
</para>
ilib_for_link(['res22' 'jac22' 'gr22'],'t22.c',[],'c',[],TMPDIR+'/t22loader.sce');
exec(TMPDIR+'/t22loader.sce')
//-3- 実行
-rtol=[1.d-6;1.d-6];atol=[1.d-6;1.d-4];
-t0=0;y0=[2;0];y0d=[0;-2];t=[20:20:200];ng=1;
+rtol = [1.d-6;1.d-6];atol=[1.d-6;1.d-4];
+t0 = 0; y0 = [2;0]; y0d = [0;-2]; t = [20:20:200]; ng = 1;
//簡単なシミュレーション
-t=0:0.003:300;
-yy=dae([y0,y0d],t0,t,atol,rtol,'res22','jac22');
-clf();plot(yy(1,:),yy(2,:))
+t = 0:0.003:300;
+yy = dae([y0,y0d],t0,t,atol,rtol,'res22','jac22');
+clf(); plot(yy(1,:),yy(2,:))
//find first point where yy(1)=0
-[yy,nn,hotd]=dae("root",[y0,y0d],t0,300,atol,rtol,'res22','jac22',ng,'gr22');
+[yy,nn,hotd] = dae("root",[y0,y0d],t0,300,atol,rtol,'res22','jac22',ng,'gr22');
plot(yy(1,1),yy(2,1),'r+')
xstring(yy(1,1)+0.1,yy(2,1),string(nn(1)))
//次の点までホットリスタート
-t01=nn(1);
-[pp,qq]=size(yy);
-y01=yy(2:3,qq);y0d1=yy(3:4,qq);
-[yy,nn,hotd]=dae("root",[y01,y0d1],t01,300,atol,rtol,'res22','jac22',ng,'gr22',hotd);
+t01 = nn(1);
+[pp,qq] = size(yy);
+y01 = yy(2:3,qq); y0d1 = yy(3:4,qq);
+[yy,nn,hotd] = dae("root",[y01,y0d1],t01,300,atol,rtol,'res22','jac22',ng,'gr22',hotd);
plot(yy(1,1),yy(2,1),'r+')
xstring(yy(1,1)+0.1,yy(2,1),string(nn(1)));
cd(previous_dir);
- ]]></programlisting>
+ ]]></programlisting>
<scilab:image localized="false"><![CDATA[
code = ['#include <math.h>'
'void res22(double *t, double *y, double *yd, double *res, int *ires, double *rpar, int *ipar)'
plot(yy(1, 1), yy(2, 1), 'r+')
xstring(yy(1, 1)+0.1, yy(2, 1), string(nn(1)));
cd(previous_dir);
- ]]></programlisting>
+ ]]></programlisting>
<scilab:image><![CDATA[
code = ['#include <math.h>'
'void res22(double *t, double *y, double *yd, double *res, int *ires, double *rpar, int *ipar)'
<?xml version="1.0" encoding="UTF-8"?>
-
<!--
* Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
* Copyright (C) 2008 - INRIA
* along with this program.
*
-->
-
<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns4="http://www.w3.org/1999/xhtml" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:id="daeoptions" xml:lang="ja">
-
<refnamediv>
-
<refname>daeoptions</refname>
-
<refpurpose>daeソルバのオプションを設定</refpurpose>
-
</refnamediv>
-
<refsynopsisdiv>
-
<title>呼出し手順</title>
-
<synopsis>daeoptions()</synopsis>
-
</refsynopsisdiv>
-
<refsection>
-
<title>説明</title>
-
<para>
-
dae関数の呼び出しコンテキストに変数
-
<literal>%DAEOPTIONS</literal>が存在した場合,
-
<literal>dae</literal>関数はこれをオプションとして設定します.
-
</para>
-
<para>
-
この daeoptions 関数は,<link linkend="dae">dae</link>
-
ソルバの様々なオプションを設定する際に
-
実行されるコマンドを対話的に表示します.
-
</para>
-
<para>
-
<caution>
-
注意: <literal>dae</literal>関数は
-
この変数が存在するかどうかを確認し,存在する場合にはこれを使用します.
-
デフォルト値を使用する場合にはこの変数を消去する必要があります.
-
<literal>daeoptions</literal> は
-
この変数を作成しないことに注意してください.この変数を作成するには,
-
<literal>daeoptions</literal>で表示されるコマンド行で実行する
-
必要があります.
-
</caution>
-
</para>
-
<para>
-
変数 <literal>%DAEOPTIONS</literal> は
-
以下の要素を有する <link linkend="list">リスト</link>です:
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
-list(tstop,imode,band,maxstep,stepin,nonneg,isest)
- ]]></programlisting>
-
+ <screen><![CDATA[
+list(tstop, imode, band, maxstep, stepin, nonneg, isest)
+]]></screen>
<para>デフォルト値は:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list([],0,[],[],[],0,0)
- ]]></programlisting>
-
+]]></screen>
<para>要素の意味を以下に示します.</para>
-
<variablelist>
-
<varlistentry>
-
<term>solver</term>
-
<listitem>
-
<para>
-
0の場合, <literal>dae</literal> は dassl/dasrtソルバを使用します,
-
1の場合, <literal>dae</literal> は daskrを使用します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>tstop</term>
-
<listitem>
-
<para>実数スカラーまたは空の行列で, 許容される
-
<literal>g </literal>の評価の最大回数を指定します.
-
空の行列は回数の"制限なし"を意味します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>imode</term>
-
<listitem>
-
<para>0 を指定した場合, dae は
-
ユーザが指定した時間での値のみを返します.
-
1を指定した場合, dae は計算した経過値を返します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>band</term>
-
<listitem>
-
<para>2個の要素を有するベクトルで,
-
<literal>jac</literal> で計算された帯行列の定義
-
<literal>[ml,mu]</literal>を指定します.
-
</para>
-
<para>
-
<literal>r(i - j + ml + mu + 1,j)</literal> =
-
<literal>dg(i)/dy(j)+cj*dg(i)/dydot(j)</literal> .
-
<literal>jac</literal> が完全な行列を返す場合,
-
<literal>band=[]</literal>を指定します
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>maxstep</term>
-
<listitem>
-
<para>スカラーまたは空の行列で, ステップの最大値.
-
空の行列は"指定なし"を意味します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>stepin</term>
-
<listitem>
-
<para>スカラーまたは空の行列で, ステップの最小値.
-
空の行列は"指定なし"を意味します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>nonneg</term>
-
<listitem>
-
<para>
-
スカラーで, <literal>解が
-
非負であることが既知の場合
-
</literal>
-
には 0 とする
-
必要があります.
-
その他の場合には 1 に設定する必要があります.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>isest</term>
-
<listitem>
-
<para>スカラーで,
-
指定した初期条件が以下と互換の場合には0 とする
-
必要があります
-
: <literal>g(t0,x0,xdot0)=0</literal>.
-
<literal>xdot0</literal>が単なる推定値である場合には
-
1 に設定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>method</term>
-
<listitem>
-
<para>
-
スカラー, 直接 / Krylov 法.
-
ソルバでKrylov反復を使用したい場合,
-
<literal>1</literal>に設定し,ルーチン
-
<literal>psol</literal>を<literal>dae()</literal>に指定します.
-
それ以外(daskrの直接法)は<literal>0</literal>を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>Kry_params</term>
-
<listitem>
-
<para>
-
ベクトル,
-
method=0を設定した場合はダミー変数として扱われます.
-
その他の場合,<literal>[]</literal>または
-
<literal>[maxl kmp nrmax epli]</literal>に設定します.
-
ただし:
-
- <literal>maxl </literal> = GMResアルゴリズムの最大反復回数 (デフォルト min(5, neq)),
-
- <literal>kmp </literal> = GMResアルゴリズムにおいて直交化を行うベクトルの数 (デフォルト maxl),
-
- <literal>nrmax</literal> = 非線形反復毎のGMResアルゴリズムの再スタートの最大回数(デフォルト 5),
-
- <literal>epli </literal> = GMResアルゴリズムの収束確認定数 (デフォルト 0.05).
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>init</term>
-
<listitem>
-
<para>
-
スカラー, isest=0と設定した場合はダミー引数として扱われます.
-
初期条件の計算後にソルバを停止する場合は1, それ以外は0を設定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>precond</term>
-
<listitem>
-
<para>
-
スカラー, <literal>psol</literal>用のプリコンディショナ計算および
-
LU分解ルーチン.
-
method=0の場合,ダミー引数として扱われます.
-
外部ルーチン<literal>psol</literal>が特定のルーチンを
-
使用する必要がある場合,
-
1に設定し, <literal>dae</literal>に
-
<literal>pjac</literal>ルーチンを提供します.
-
それ以外は 0を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>control</term>
-
<listitem>
-
<para>
-
ベクトル, isest=0を指定した場合はダミー引数として扱われます.
-
全ての変数のエラーをローカルに制御したい場合,
-
[] を指定します.
-
そうでない場合, [+-1, ..., +-1]を設定します.
-
ただし,
-
<literal>y(i)</literal> が微分変数の場合は 1,
-
<literal>y(i)</literal> が代数変数
-
(この微分が関数<literal>g(t, y, ydot)</literal>に明示的に
-
現れない)の場合は -1とします.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>heuristic</term>
-
<listitem>
-
<para>
-
ベクトル, 経験的パラメータ.
-
isest=0の場合はダミー引数として扱われます.
-
そうでない場合, [] または
-
<literal>[mxnit mxnj mxnh lsoff stptol epinit]</literal>に
-
設定します, ただし:
-
- <literal>mxnit</literal> = ヤコビアンまたはプリコンディショナ評価毎のニュートン反復の最大回数
-
(デフォルト: info(8)=0の場合は 5, それ以外は 15),
-
- <literal>mxnj</literal> = ヤコビアンまたはプリコンディショナ評価の最大回数
-
(デフォルト: info(8)=0の場合は 6, それ以外は 2),
-
- <literal>mxnh</literal> = info(7) ≠ [] の時に試行する人工的ステップサイズパラメータ h の最大値
-
(デフォルト 5),
-
- <literal>lsoff</literal> = ライン探索アルゴリズムを無効にするフラグ
-
(lsoff = 0 はライン探索をオン, lsoff = 1 はオフにすることを意味します) (デフォルト: 0),
-
- <literal>stptol</literal> = ライン探索アルゴリズムの最小ステップ (デフォルト (unit roundoff)^(2/3)),
-
- <literal>epinit</literal> = ニュートン反復収束確認のスイング係数(デフォルト 0.01).
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>verbosity</term>
-
<listitem>
-
<para>
-
スカラー, 冗長度. 外部出力を最小にする場合は 1, 完全に出力する場合は 2,
-
それ以外は 0に設定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</refsection>
-
<refsection role="see also">
-
<title>参照</title>
-
<simplelist type="inline">
-
<member>
-
<link linkend="dae">dae</link>
-
</member>
-
</simplelist>
-
</refsection>
-
-</refentry>
-
+</refentry>
\ No newline at end of file
<?xml version="1.0" encoding="UTF-8"?>
-
<!--
* Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
* Copyright (C) 2013 - Scilab Enterprises
* along with this program.
*
-->
-
<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:id="daskr" xml:lang="ja">
-
<refnamediv>
-
<refname>daskr</refname>
-
<refpurpose>ゼロ交差を有するDAEソルバー</refpurpose>
-
</refnamediv>
-
<refsynopsisdiv>
-
<title>呼び出し手順</title>
-
<synopsis>[r, nn [, hd]] = daskr(x0, t0, t [, atol [, rtol]], res [, jac], ng, surf [, info [, psol] [, pjac]] [, hd])</synopsis>
-
</refsynopsisdiv>
-
<refsection>
-
<title>引数</title>
-
<variablelist>
-
<varlistentry>
-
<term>x0</term>
-
<listitem>
-
<para>
-
は <literal>y0</literal> (<literal>ydot0</literal> は
-
ゼロを初期推定値として <literal>dassl</literal> により推定)
-
または行列 <literal>[y0 ydot0]</literal>のどちらかです.
-
<literal>g(t,y0,ydot0)</literal> はゼロに等しい必要があります.
-
<literal>ydot0</literal> の推定値のみ既知の場合,
-
<literal>info(7)=1</literal> を指定してください.
-
</para>
-
<variablelist>
-
<varlistentry>
-
<term>y0</term>
-
<listitem>
-
<para>初期条件の実数列ベクトル.</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>ydot0</term>
-
<listitem>
-
<para>
-
<literal>y</literal> at <literal>t0</literal> (推定値でも可)における<literal>y</literal>の時間微分の実数列ベクトル.
-
</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>t0</term>
-
<listitem>
-
<para>実数, 初期時間.</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>t</term>
-
<listitem>
-
<para>実数スカラーまたはベクトル.
-
解を得る時刻を指定します.
-
<literal>info(2)=1</literal>と指定することにより,
-
dasslの各ステップ点で解を得ることができます.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>nn</term>
-
<listitem>
-
<para>
-
2個のエントリ <literal>[times num]</literal>を有するベクトル.
-
<literal>times</literal> は面が交差する点における時間の値,
-
<literal>num</literal> は交差する面の数です.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>atol, rtol</term>
-
<listitem>
-
<para>
-
<literal>y</literal>と同じ大きさの実数スカラーまたは列ベクトル.
-
<literal>atol, rtol</literal> はそれぞれ解の絶対及び相対許容誤差を
-
表します. この許容誤差は
-
<literal>y</literal>の各要素について指定されます.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>res</term>
-
<listitem>
-
<para>
-
<link linkend="external">external</link> (関数またはリストまたは文字列).
-
<literal>g(t,y,ydot)</literal>の値を計算します. 以下のようになります :
-
</para>
-
<itemizedlist>
-
<listitem>
-
<para>Scilab関数.</para>
-
<para>呼び出し手順は以下とします
-
<literal>[r,ires]=res(t,y,ydot)</literal> かつ
-
<literal>res</literal> は残差
-
<literal>r=g(t,y,ydot)</literal> かつ エラーフラグ
-
<literal>ires</literal>を返すこと.
-
<literal>res</literal> は, <literal>r</literal>の計算に成功した場合は<literal>ires = 0</literal>,
-
残差が<literal>(t,y,ydot)</literal> において局所的に定義されないには
-
<literal>=-1</literal>, パラメータが許容範囲を逸脱している場合には <literal>=-2</literal>
-
となります.
-
</para>
-
</listitem>
-
<listitem>
-
<para>リスト.</para>
-
<para>この形式により, t, y, ydot 以外のパラメータを関数に指定出来ます.
-
以下のような形式とする必要があります:
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(res,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>ただし関数
-
<literal>res</literal> の呼び出し形式は以下とします
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r=res(t,y,ydot,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>
-
<literal>res</literal> は
-
<literal>(t,y,ydot,x1,x2,...)</literal>の関数として<literal>r=g(t,y,ydot)</literal>を返します.
-
</para>
-
<para>
-
警告: この形式は<literal>関数</literal>に渡す別の引数が存在しない場合には
-
使用するべきではありません.
-
</para>
-
</listitem>
-
<listitem>
-
<para>文字列.</para>
-
<para>ScilabにリンクされたCまたはFortranサブルーチンの名前を指定します.
-
</para>
-
<para>Cの場合, 呼び出し手順は次のようにします:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void res(double *t, double y[], double yd[], double r[],
int *ires, double rpar[], int ipar[])
- ]]></programlisting>
-
+]]></screen>
<para>Fortranの場合, 以下とします:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine res(t,y,yd,r,ires,rpar,ipar)
double precision t, y(*),yd(*),r(*),rpar(*)
integer ires,ipar(*)
- ]]></programlisting>
-
+]]></screen>
<para>
-
<literal>rpar</literal> および <literal>ipar</literal> 配列は存在する必要がありますが使用できません.
-
</para>
-
</listitem>
-
</itemizedlist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>jac</term>
-
<listitem>
-
<para>
-
<link linkend="external">external</link> (関数またはリストまたは文字列). 指定した<literal>cj</literal>
-
の値について, <literal>dg/dy + cj*dg/dydot</literal> の値を計算します.
-
</para>
-
<itemizedlist>
-
<listitem>
-
<para>Scilab関数.</para>
-
<para>
-
呼び出し手順は <literal>r=jac(t,y,ydot,cj)</literal> で,
-
<literal>jac</literal> 関数は
-
<literal>r=dg(t,y,ydot)/dy+cj*dg(t,y,ydot)/dydot</literal>を返す必要があります.
-
ただし,
-
<literal>cj</literal> は実数スカラーです.
-
</para>
-
</listitem>
-
<listitem>
-
<para>リスト.</para>
-
<para>以下のように指定します</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(jac,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>ただし関数
-
<literal>jac</literal> の呼び出し手順は以下のようになります
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r=jac(t,y,ydot,cj,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>
-
さらに<literal>jac</literal> は
-
<literal>(t,y,ydot,cj,x1,x2,...)</literal>の関数として<literal>dg/dy + cj*dg/dydot</literal> を返します.
-
</para>
-
</listitem>
-
<listitem>
-
<para>文字列.</para>
-
<para>scilabにリンクされたCまたはFortran サブルーチンの名前を指定します
-
</para>
-
<para>Cの場合, 呼び出し手順は以下となります:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void jac(double *t, double y[], double yd[], double pd[],
double *cj, double rpar[], int ipar[])
- ]]></programlisting>
-
+]]></screen>
<para>Fortranの場合は以下となります:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine jac(t,y,yd,pd,cj,rpar,ipar)
double precision t, y(*),yd(*),pd(*),cj,rpar(*)
integer ipar(*)
- ]]></programlisting>
-
+]]></screen>
</listitem>
-
</itemizedlist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>surf</term>
-
<listitem>
-
<para>
-
<link linkend="external">external</link> (関数またはリストまたは文字列).
-
k
-
<literal>ng</literal>個の要素を有する列ベクトル <literal>surf(t,y)</literal> の値を計算します.
-
各要素は局面を定義します.
-
これは以下のように定義できます:
-
</para>
-
<itemizedlist>
-
<listitem>
-
<para>Scilab関数.</para>
-
<para>呼び出し手順は
-
<literal>surf(t,y)</literal>とします
-
</para>
-
</listitem>
-
<listitem>
-
<para>リスト.</para>
-
<para>以下のように定義します</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(surf,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>ただし関数
-
<literal>surf</literal>の呼び出し手順は以下とします
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r=surf(t,y,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
</listitem>
-
<listitem>
-
<para>文字列.</para>
-
<para>ScilabにリンクされたC または Fortran サブルーチンの名前を指定します.
-
</para>
-
<para>Cの場合, 呼び出し手順は以下とします:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void surf(int *ny, double *t, double y[], int *ng, double gout[])
- ]]></programlisting>
-
+]]></screen>
<para>Fortranの場合, 以下とします:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine surf(ny,t,y,ng,gout)
double precision t, y(*),gout(*)
integer ny,ng
- ]]></programlisting>
-
+]]></screen>
</listitem>
-
</itemizedlist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info</term>
-
<listitem>
-
<para>
-
<literal>14</literal> 個の要素を有するリスト. デフォルト値は
-
<literal>list([],0,[],[],[],0,[],0,[],0,[],[],[],1)</literal>です.
-
</para>
-
<variablelist>
-
<varlistentry>
-
<term>info(1)</term>
-
<listitem>
-
<para>
-
<literal>g</literal> を評価できる時刻の最大値を実数スカラーで指定します.
-
時間を制限しない場合は空の行列<literal>[]</literal>を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(2)</term>
-
<listitem>
-
<para>
-
<literal>dassl</literal> が中間計算値を返すか(<literal>flag=1</literal>)
-
またはユーザが指定した時刻のみを返すか
-
(<literal>flag=0</literal>)を指定するフラグ.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(3)</term>
-
<listitem>
-
<para>
-
<literal>2</literal> 要素のベクトルで,
-
<literal>jac</literal>; <literal>r(i - j + ml + mu + 1,j) =
-
"dg(i)/dy(j)+cj*dg(i)/dydot(j)"
-
</literal>
-
で計算するバンド行列の定義域<literal>[ml,mu]</literal>を指定します.
-
<literal>info(3)=[]</literal>の場合, <literal>jac</literal> は通常の行列
-
を返します.
-
<literal>info(8)=1</literal>の場合, ダミーとして扱われます.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(4)</term>
-
<listitem>
-
<para>最大ステップサイズを指定する実数スカラー値.
-
制限しない場合は <literal>info(4)=[]</literal>を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(5)</term>
-
<listitem>
-
<para>初期ステップサイズを指定する実数スカラー値.
-
指定しない場合は <literal>info(5)=[]</literal>を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(6)</term>
-
<listitem>
-
<para>
-
解が非負であることが既知な場合は <literal>info(6)=1</literal>,
-
それ以外は <literal>info(6)=0</literal> を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(7)</term>
-
<listitem>
-
<para>
-
ydot0 が<literal>g(t0,y0,ydot0)=0</literal>のように指定された場合,
-
<literal>info(7)=[]</literal> と指定します.
-
それ以外の場合は <literal>info(7)=[+-1, ..., +-1]</literal>とします.
-
ここで,y(i)が微分変数の場合は<literal>info(7)(i)=1</literal>,
-
y(i)が代数変数(関数g(t,y,ydot)に微分値が陽に現れない)の場合は
-
<literal>info(7)(i)=-1</literal>とします.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(8)</term>
-
<listitem>
-
<para>
-
直接法 / Krylov.
-
ソルバーでKrylov反復法を使用したい場合,<literal>info(8)=1</literal>と指定し,
-
<literal>psol</literal>にルーチンを指定します.
-
それ以外の場合(dasslの直接法)は<literal>info(8)=0</literal>を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(9)</term>
-
<listitem>
-
<para>
-
Krylovパラメータ.
-
<literal>info(8)=0</literal>とした場合はダミー引数として処理されます.
-
それ以外の場合は
-
<literal>info(9)=[maxl kmp nrmax epli]</literal>を指定します, ただし:
-
</para>
-
<para>
-
- maxl = mSPIGMRアルゴリズムの最大反復回数 (デフォルト
-
<literal>min(5,neq)</literal>),
-
</para>
-
<para>
-
- kmp = SPIGMRアルゴリズムで直交化を行うベクトルの数
-
(デフォルト maxl),
-
</para>
-
<para>
-
- nrmax = 各非線形反復におけるSPIGMRアルゴリズムのリスタート最大数
-
(デフォルト <literal>5</literal>),
-
</para>
-
<para>
-
- epli = SPIGMRアルゴリズムの収束安定定数 (デフォルト <literal>0.05</literal>).
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(10)</term>
-
<listitem>
-
<para>
-
初期条件.
-
<literal>info(7)=[]</literal>の場合はダミーとして処理されます.
-
初期条件の計算直後にソルバーを停止させる場合は
-
<literal>info(10)=1</literal>と指定し,それ以外は
-
<literal>info(10)=0</literal>と指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(11)</term>
-
<listitem>
-
<para>
-
<literal>psol</literal>用のプリコンディショナ計算およびLU分解ルーチン.
-
<literal>info(8)=0</literal>の場合はダミーとして扱われます.
-
<link linkend="external">external</link> <literal>psol</literal>が特定の
-
ルーチンを使用する必要がある場合は,<literal>info(11)=1</literal>を指定し,
-
<literal>pjac</literal>ルーチンを指定します.
-
それ以外は<literal>info(11)=0</literal>とします.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(12)</term>
-
<listitem>
-
<para>
-
全変数についてエラー制御をローカルに行いたい場合,
-
<literal>info(12)=[]</literal>を指定します. それ以外は,
-
<literal>info(12)=[+-1, ..., +-1]</literal>を指定します.
-
ただし, y(i)が微分変数の場合は <literal>info(12)(i)=1</literal>,
-
y(i)が代数変数(微係数が関数g(t,y,ydot)に陽に現れない) の場合は
-
<literal>info(12)(i)=-1</literal>とします.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(13)</term>
-
<listitem>
-
<para>
-
経験的なパラメータ.
-
<literal>info(7)=[]</literal>の場合はダミーとして扱われます.
-
それ以外の場合,
-
<literal>info(13)=[mxnit mxnj mxnh lsoff stptol epinit]</literal>を指定します,
-
ただし:
-
</para>
-
<para>
-
- mxnit = ヤコビアンまたはプリコンディショナ評価毎のニュートン反復の最大回数
-
(デフォルト:<literal>info(8)=0</literal>の場合は <literal>5</literal>,
-
それ以外は <literal>15</literal>),
-
</para>
-
<para>
-
- mxnj = ヤコビアンまたはプリコンディショナ評価の最大回数
-
(デフォルト:<literal>info(8)=0</literal>の場合は <literal>6</literal>,
-
それ以外は <literal>2</literal>),
-
</para>
-
<para>
-
- mxnh = info(7) ≠ [] の場合に試行する人工的なステップサイズ
-
パラメータHの最大値
-
(デフォルト <literal>5</literal>),
-
</para>
-
<para>
-
- lsoff = ライン探索アルゴリズムをオフにするフラグ
-
(lsoff = 0 はライン探索有効を意味します, lsoff = 1 はオフを意味します)
-
(デフォルト: <literal>0</literal>),
-
</para>
-
<para>
-
- stptol = ライン探索アルゴリズムの最小スケール刻み (デフォルト: <literal>(丸め誤差単位)^(2/3)</literal>),
-
</para>
-
<para>
-
- epinit = ニュートン反復収束判定の振幅係数 (デフォルト <literal>0.01</literal>).
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(14)</term>
-
<listitem>
-
<para>
-
冗長性.
-
出力を最小化する場合は <literal>info(14)=1</literal>,
-
最大化する場合は <literal>info(14)=2</literal>, それ以外は
-
<literal>info(14)=0</literal>を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>psol</term>
-
<listitem>
-
<para>
-
<link linkend="external">external</link> (関数またはリストまたは文字列).
-
線形システム<literal>P*x = b</literal>を解きます.
-
ただし, P はプリコンディショナで, ルーチン <literal>pjac</literal>
-
が事前に計算して, wp 及び iwpに保存します.
-
</para>
-
<itemizedlist>
-
<listitem>
-
<para>Scilab関数.</para>
-
<para>呼び出し手順は
-
<literal>[r, ier] = psol(wp, iwp, b)</literal> とします.
-
<literal>psol</literal> 関数はシステムの解<literal>r</literal>と
-
エラーフラグ<literal>ier</literal>を返す必要があります.
-
</para>
-
</listitem>
-
<listitem>
-
<para>リスト.</para>
-
<para>以下のような形式とします</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(psol,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>ただし, 関数
-
<literal>psol</literal>の呼び出し手順は以下とします
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
psol(wp,iwp,b,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>
-
<literal>psol</literal> は <literal>r</literal>の解を返します.
-
</para>
-
</listitem>
-
<listitem>
-
<para>文字列.</para>
-
<para>ScilabにリンクされたCまたはFortranサブルーチンの名前を指定します
-
</para>
-
<para>Cの場合,呼び出し手順は以下とします:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void psol (int*neq, double*t, double*y, double*ydot, double*savr,
double*wk, double*cj, double*wght, double*wp, int*iwp, double*b, double*eplin, int*ier, double*rpar, int*ipar)
- ]]></programlisting>
-
+]]></screen>
ただし,配列 <literal>wp</literal> および <literal>iwp</literal> はプリコンディショナ
-
<literal>P</literal>の行列要素を行圧縮形式で有します.
-
<para>Fortranの場合は以下とします:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine psol (neq, t, y, ydot, savr, wk, cj, wght,
wp, iwp, b, eplin, ier, rpar, ipar)
double precision t,y(*),ydot(*),savr(*),wk(*),cj,wght(*),wp(*),
b(*),eplin,rpar(*)
integer neq,iwp(*),ier,ipar(*)
- ]]></programlisting>
-
+]]></screen>
</listitem>
-
</itemizedlist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>pjac</term>
-
<listitem>
-
<para>
-
<link linkend="external">external</link> (関数またはリストまたは文字列).
-
指定したパラメータ<literal>cj</literal>について
-
<literal>dg/dy + cj*dg/dydot</literal>を計算し,
-
doubleおよびint配列にLU分解を行います.
-
</para>
-
<itemizedlist>
-
<listitem>
-
<para>Scilab関数.</para>
-
<para>呼び出し手順は
-
<literal>[wp, iwp, ires] = pjac(neq, t, y, ydot, h, cj, rewt, savr)</literal>
-
となり, 結果として,
-
配列 wp および iwp には全てのプリコンディショナ情報が
-
圧縮スパース行形式で返されます.
-
</para>
-
</listitem>
-
<listitem>
-
<para>リスト.</para>
-
<para>以下の形式とします</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(pjac,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>ただし,関数
-
<literal>pjac</literal>の呼び出し手順は以下とします
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
pjac(neq,t,y,ydot,h,cj,rewt,savr,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>
-
<literal>pjac</literal> は,
-
<literal>(neq,t,y,ydot,h,cj,rewt,savr,x1,x2,...)</literal>の関数として
-
分解した <literal>dg/dy + cj*dg/dydot</literal> を返します.
-
</para>
-
</listitem>
-
<listitem>
-
<para>文字列.</para>
-
<para>ScilabにリンクされたCまたはFortran形式のサブルーチン名を指定します
-
</para>
-
<para>Cの場合,呼び出し手順は以下とします:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void pjac (double*res, int*ires, int*neq, double*t, double*y, double*ydot, double*rewt, double*savr,
double*wk, double*h, double*cj, double*wp, int*iwp, int*ier, double*rpar, int*ipar)
- ]]></programlisting>
-
+]]></screen>
<para>Fortranの場合, 以下のとします:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine pjac (res, ires, neq, t, y, ydot, rewt, savr,
wk, h, cj, wp, iwp, ier, rpar, ipar)
double precision res(*), t, y(*), ydot(*), rewt(*), savr(*),
wk(*), h, cj, wp(*), rpar(*)
integer ires, neq, iwp(*), ier, ipar(*)
- ]]></programlisting>
-
+]]></screen>
</listitem>
-
</itemizedlist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>hd</term>
-
<listitem>
-
<para>
-
<literal>dassl</literal>コンテキストを保存し,積分を再開できるようにする実数ベクトル.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>r</term>
-
<listitem>
-
<para>
-
実数行列 .
-
各列はベクトル<literal>[t;x(t);xdot(t)]</literal>で,
-
<literal>t</literal> は解が計算される添字の時刻です.
-
</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</refsection>
-
<refsection>
-
<title>説明</title>
-
<para>陰的微分方程式の解.</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
g(t,y,ydot) = 0
y(t0) = y0 and ydot(t0) = ydot0
- ]]></programlisting>
-
+]]></screen>
<para>
-
局面交差した時刻と面交差した数を<literal>nn</literal>に返します.
-
</para>
-
<para>詳細な例は SCI/modules/differential_equations/tests/unit_tests/daskr.tst にあります.</para>
-
</refsection>
-
<refsection>
-
<title>例</title>
-
<programlisting role="example"><![CDATA[
// dy/dt = ((2*log(y)+8)/t -5)*y, y(1) = 1, 1 <= t <=6
// g1 = ((2*log(y)+8)/t - 5)*y
[yy,nn] = daskr([y0,y0d],t0,t,atol,rtol,res1,ng,gr1);
nn
// nn = [2.4698972 2] を返します
- ]]></programlisting>
-
+ ]]></programlisting>
</refsection>
-
<refsection role="see also">
-
<title>参照</title>
-
<simplelist type="inline">
-
<member>
-
<link linkend="ode">ode</link>
-
</member>
-
<member>
-
<link linkend="dasrt">dasrt</link>
-
</member>
-
<member>
-
<link linkend="dassl">dassl</link>
-
</member>
-
<member>
-
<link linkend="impl">impl</link>
-
</member>
-
<member>
<link linkend="call">call</link>
</member>
-
<member>
-
<link linkend="link">link</link>
-
</member>
-
<member>
-
<link linkend="external">external</link>
-
</member>
-
</simplelist>
-
</refsection>
-
<refsection>
-
<title>履歴</title>
-
<revhistory>
-
<revision>
-
<revnumber>5.5.0</revnumber>
-
<revdescription>daskr ソルバーが追加されました</revdescription>
-
</revision>
-
</revhistory>
-
</refsection>
-
-</refentry>
-
+</refentry>
\ No newline at end of file
<?xml version="1.0" encoding="UTF-8"?>
-
<!--
* Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
* Copyright (C) 2008 - INRIA
* along with this program.
*
-->
-
<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:id="dasrt" xml:lang="ja">
-
<refnamediv>
-
<refname>dasrt</refname>
-
<refpurpose>ゼロ交差するDAE</refpurpose>
-
</refnamediv>
-
<refsynopsisdiv>
-
<title>呼び出し手順</title>
-
- <synopsis>[r,nn,[,hd]]=dasrt(x0,t0,t [,atol,[rtol]],res [,jac],ng, surf [,info] [,hd])</synopsis>
-
+ <synopsis>[r,nn,[,hd]] = dasrt(x0,t0,t [,atol,[rtol]],res [,jac],ng, surf [,info] [,hd])</synopsis>
</refsynopsisdiv>
-
<refsection>
-
<title>引数</title>
-
<variablelist>
-
<varlistentry>
-
<term>x0</term>
-
<listitem>
-
<para>
-
<literal>y0</literal> (<literal>ydot0</literal> は
-
0 を初期推定値として<literal>dassl</literal>で推定されたもの)
-
または行列<literal>[y0 ydot0]</literal>のどちらかです.
-
<literal>g(t,y0,ydot0)</literal>は 0 に等しい必要があります.
-
<literal>ydot0</literal>の推定値のみが既知の場合,
-
<literal>info(7)=1</literal> と設定してください.
-
</para>
-
<variablelist>
-
<varlistentry>
-
<term>y0</term>
-
<listitem>
-
<para>初期条件の実数列ベクトル.</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>ydot0</term>
-
<listitem>
-
<para>
-
<literal>t0</literal>における<literal>y</literal>の時間微分の
-
実数列ベクトル (推定値の場合もある).
-
</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>t0</term>
-
<listitem>
-
<para>初期時間を表す実数.</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>t</term>
-
<listitem>
-
<para>実数のスカラーおよびベクトル. 解を求める時間を指定します.
-
<literal>info(2)=1</literal>を設定することにより,
-
dasslの各ステップ点での解を得ることができることに注意してください.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>nn</term>
-
<listitem>
-
<para>
-
2個のエントリ <literal>[times num]</literal>を有するベクトル.
-
<literal>times</literal> は面が交差した時刻の値,
-
<literal>num</literal> は交差した面の数です
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>atol, rtol</term>
-
<listitem>
-
<para>
-
実数スカラーまたは<literal>y</literal>と同じ大きさの列ベクトル.
-
<literal>atol, rtol</literal> はそれぞれ絶対および相対許容誤差を指定します.
-
ベクトルの場合,<literal>y</literal>の各要素毎に許容誤差が指定されます.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>res</term>
-
<listitem>
-
<para>外部 (関数またはリストまたは文字列).
-
<literal>g(t,y,ydot)</literal>の値を計算します.次のようになります :
-
</para>
-
<itemizedlist>
-
<listitem>
-
<para>Scilab関数.</para>
-
<para>この場合, その呼出し手順を
-
<literal>[r,ires]=res(t,x,xdot)</literal> とする
-
必要があり,<literal>res</literal> は
-
残差<literal>r=g(t,x,xdot)</literal> とエラーフラグ
-
<literal>ires</literal>を返す必要があります.
-
<literal>res</literal>が<literal>r</literal>の計算に
-
成功した場合には<literal>ires = 0</literal>,
-
<literal>(t,x,xdot)</literal>のローカルな残差が定義されない
-
場合には <literal>=-1</literal>,
-
パラメータが許容範囲外の場合は <literal>=-2</literal> となります.
-
</para>
-
</listitem>
-
<listitem>
-
<para>リスト.</para>
-
<para>外部ルーチンのこの形式は
-
関数にt,y,ydot以外のパラメータを指定する際に使用されます.
-
以下のような形式とします:
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(res,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>
-
ただし,ここで関数<literal>res</literal>の呼び出し手順は
-
以下のようになります
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
-r=res(t,y,ydot,x1,x2,...)
- ]]></programlisting>
-
+ <screen><![CDATA[
+r = res(t,y,ydot,x1,x2,...)
+]]></screen>
<para>
-
<literal>res</literal> はこの場合も
-
<literal>(t,y,ydot,x1,x2,...)</literal>の関数として
-
<literal>r=g(t,y,ydot)</literal>を返します.
-
</para>
-
<para>警告: 関数にオプションの引数を指定しない場合,
-
この形式を使用する必要はありません
-
</para>
-
</listitem>
-
<listitem>
-
<para>文字列.</para>
-
<para>
-
ScilabにリンクされるCまたはFortranサブルーチンの名前を指します.
-
</para>
-
<para>Cの場合, 呼出し手順は以下のようにします:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void res(double *t, double y[], double yd[], double r[],
int *ires, double rpar[], int ipar[])
- ]]></programlisting>
-
+]]></screen>
<para>Fortranの場合, 呼出し手順は以下のようにします:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine res(t,y,yd,r,ires,rpar,ipar)
double precision t, y(*),yd(*),r(*),rpar(*)
integer ires,ipar(*)
- ]]></programlisting>
-
+]]></screen>
<para>
-
<literal>rpar</literal> および <literal>ipar</literal> 配列は存在しますが使用できません.
-
</para>
-
</listitem>
-
</itemizedlist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>jac</term>
-
<listitem>
-
<para>外部ルーチン (関数またはリストまたは文字列).
-
パラメータ<literal>cj</literal>の指定した値に関して
-
<literal>dg/dy+cj*dg/dydot</literal>の値を計算します.
-
</para>
-
<itemizedlist>
-
<listitem>
-
<para>Scilab関数.</para>
-
<para>その呼び出し手順は
-
<literal>r=jac(t,y,ydot,cj)</literal> とする必要があり,
-
<literal>jac</literal>関数は
-
<literal>r=dg(t,y,ydot)/dy+cj*dg(t,y,ydot)/dydot</literal>を
-
返す必要があります.
-
ただし, <literal>cj</literal>は実数スカラーです
-
</para>
-
</listitem>
-
<listitem>
-
<para>リスト.</para>
-
<para>以下のようにする必要があります</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(jac,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>
-
ただし,この場合の関数<literal>jac</literal>の呼び出し手順は
-
以下のようになります
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
-r=jac(t,y,ydot,cj,x1,x2,...)
- ]]></programlisting>
-
+ <screen><![CDATA[
+r = jac(t,y,ydot,cj,x1,x2,...)
+]]></screen>
<para>
-
この場合も<literal>jac</literal> は
-
<literal>(t,y,ydot,cj,x1,x2,...)</literal>の関数として
-
<literal>dg/dy+cj*dg/dydot</literal>を返します.
-
</para>
-
</listitem>
-
<listitem>
-
<para>文字列.</para>
-
<para>
-
ScilabにリンクされるFortranサブルーチンの名前を指す
-
必要があります
-
</para>
-
<para>Cの場合,呼び出し手順は次のようになります:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void jac(double *t, double y[], double yd[], double pd[],
double *cj, double rpar[], int ipar[])
- ]]></programlisting>
-
+]]></screen>
<para>Fortranの場合, 呼び出し手順は次のようになります:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine jac(t,y,yd,pd,cj,rpar,ipar)
double precision t, y(*),yd(*),pd(*),cj,rpar(*)
integer ipar(*)
- ]]></programlisting>
-
+]]></screen>
</listitem>
-
</itemizedlist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>surf</term>
-
<listitem>
-
<para>外部ルーチン (関数またはリストまたは文字列).
-
<literal>ng</literal>個の要素を有する列ベクトル
-
<literal>surf(t,y)</literal>の値を計算します.
-
各要素は一つの面を定義します.
-
以下のように定義されます:
-
</para>
-
<itemizedlist>
-
<listitem>
-
<para>Scilab関数.</para>
-
<para>呼び出し手順は
-
<literal>surf(t,y)</literal>となります
-
</para>
-
</listitem>
-
<listitem>
-
<para>リスト.</para>
-
<para>以下のようになります</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(surf,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>
-
ただし,この場合の関数<literal>surf</literal>の
-
呼び出し手順は以下のようになります
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
-r=surf(t,y,x1,x2,...)
- ]]></programlisting>
-
+ <screen><![CDATA[
+r = surf(t,y,x1,x2,...)
+]]></screen>
</listitem>
-
<listitem>
-
<para>文字列.</para>
-
<para>ScilabにリンクされるCまたはFortranサブルーチンの名前を指します</para>
-
<para>Cの場合,呼び出し手順は次のようになります:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void surf(int *ny, double *t, double y[], int *ng, double gout[])
- ]]></programlisting>
-
+]]></screen>
<para>Fortranの場合,呼び出し手順は次のようになります:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine surf(ny,t,y,ng,gout)
double precision t, y(*),gout(*)
integer ny,ng
- ]]></programlisting>
-
+]]></screen>
</listitem>
-
</itemizedlist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info</term>
-
<listitem>
-
<para>
-
<literal>7</literal> 個の要素を含むオプションのリスト.
-
デフォルト値は list([],0,[],[],[],0,0);
-
</para>
-
<variablelist>
-
<varlistentry>
-
<term>info(1)</term>
-
<listitem>
-
<para>
-
<literal>g</literal>を評価できる時間の最大値
-
を指定する実数スカラー.
-
または,時間に制限を課さない場合には,空の行列
-
<literal>[]</literal>
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(2)</term>
-
<listitem>
-
<para>
-
<literal>dasrt</literal> が
-
その中間的な計算値を返す(<literal>flag=1</literal>),
-
もしくはユーザが指定した時間においてのみ値を返す
-
(<literal>flag=0</literal>)
-
かどうかを指定するフラグ.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(3)</term>
-
<listitem>
-
<para>
-
<literal>2</literal> 個の要素のベクトルで,
-
<literal>jac</literal>; <literal>r(i - j + ml + mu + 1,j) =
-
"dg(i)/dy(j)+cj*dg(i)/dydot(j)"
-
</literal>
-
により計算された
-
帯行列の定義<literal>[ml,mu]</literal>を指定します.
-
<literal>jac</literal> が完全な行列を返す場合,
-
<literal>info(3)=[]</literal> を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(4)</term>
-
<listitem>
-
<para>実数スカラーで,
-
ステップの大きさの最大値を指定します.
-
制限をしない場合には,<literal>info(4)=[]</literal>
-
を指定してください.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(5)</term>
-
<listitem>
-
<para>実数スカラーで,ステップの大きさの初期値を指定します.
-
指定しない場合には,<literal>info(5)=[]</literal>
-
を指定してください.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(6)</term>
-
<listitem>
-
<para>
-
解が非負であることが既知の場合に<literal>info(6)=1</literal>に設定,
-
そうでない場合は <literal>info(6)=0</literal> に設定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(7)</term>
-
<listitem>
-
<para>
-
<literal>ydot0</literal>が単なる推定値の場合は
-
<literal>info(7)=1</literal>,
-
<literal>g(t0,y0,ydot0)=0</literal>の場合は
-
<literal>info(7)=0</literal> を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>hd</term>
-
<listitem>
-
<para>
-
実数ベクトルで,<literal>dasrt</literal>のコンテキストを保存し,
-
積分を再開することを可能にします.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>r</term>
-
<listitem>
-
<para>
-
実数行列. 各列はベクトル<literal>[t;x(t);xdot(t)]</literal>です.
-
ただし <literal>t</literal> は解を計算した時間です.
-
</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</refsection>
-
<refsection>
-
<title>説明</title>
-
<para>陰的な微分方程式を解きます.</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
-g(t,y,ydot)=0
-y(t0)=y0 and ydot(t0)=ydot0
- ]]></programlisting>
-
+ <screen><![CDATA[
+g(t,y,ydot) = 0
+y(t0) = y0 and ydot(t0) = ydot0
+]]></screen>
<para>
-
面を横切る時間および
-
到達した面の数を<literal>nn</literal>に返します.
-
</para>
-
<para>詳細な例は SCIDIR/tests/unit_tests/dassldasrt.tst で見つけることができます</para>
-
</refsection>
-
<refsection>
-
<title>例</title>
-
<programlisting role="example"><![CDATA[
//dy/dt = ((2*log(y)+8)/t -5)*y, y(1) = 1, 1<=t<=6
//g1 = ((2*log(y)+8)/t - 5)*y
//g2 = log(y) - 2.2491
-y0=1;t=2:6;t0=1;y0d=3;
-atol=1.d-6;rtol=0;ng=2;
-deff('[delta,ires]=res1(t,y,ydot)','ires=0;delta=ydot-((2*log(y)+8)/t-5)*y')
-deff('[rts]=gr1(t,y)','rts=[((2*log(y)+8)/t-5)*y;log(y)-2.2491]')
-[yy,nn]=dasrt([y0,y0d],t0,t,atol,rtol,res1,ng,gr1);
+y0 = 1; t = 2:6; t0 = 1; y0d = 3;
+atol = 1.d-6; rtol = 0; ng = 2;
+deff('[delta, ires] = res1(t,y,ydot)', 'ires=0; delta=ydot-((2*log(y)+8)/t-5)*y')
+deff('rts = gr1(t,y)', 'rts = [((2*log(y)+8)/t-5)*y;log(y)-2.2491]')
+[yy, nn] = dasrt([y0,y0d], t0, t, atol, rtol, res1, ng, gr1);
//(Should return nn=[2.4698972 2])
- ]]></programlisting>
-
+ ]]></programlisting>
</refsection>
-
<refsection role="see also">
-
<title>参照</title>
-
<simplelist type="inline">
-
<member>
-
<link linkend="ode">ode</link>
-
</member>
-
<member>
-
<link linkend="dassl">dassl</link>
-
</member>
-
<member>
-
<link linkend="impl">impl</link>
-
</member>
-
<member>
<link linkend="call">call</link>
</member>
-
<member>
-
<link linkend="link">link</link>
-
</member>
-
<member>
-
<link linkend="external">external</link>
-
</member>
-
</simplelist>
-
</refsection>
-
-</refentry>
-
+</refentry>
\ No newline at end of file
<?xml version="1.0" encoding="UTF-8"?>
-
<!--
* Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
* Copyright (C) 2008 - INRIA
* along with this program.
*
-->
-
<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:id="dassl" xml:lang="ja">
-
<refnamediv>
-
<refname>dassl</refname>
-
<refpurpose>微分代数方程式</refpurpose>
-
</refnamediv>
-
<refsynopsisdiv>
-
<title>呼出し手順</title>
-
- <synopsis>[r [,hd]]=dassl(x0,t0,t [,atol,[rtol]],res [,jac] [,info] [,hd])</synopsis>
-
+ <synopsis>[r [,hd]] = dassl(x0,t0,t [,atol,[rtol]],res [,jac] [,info] [,hd])</synopsis>
</refsynopsisdiv>
-
<refsection>
-
<title>引数</title>
-
<variablelist>
-
<varlistentry>
-
<term>x0</term>
-
<listitem>
-
<para>
-
<literal>y0</literal> (<literal>ydot0</literal> は
-
0 を初期推定値として<literal>dassl</literal>で推定されたもの)
-
または行列<literal>[y0 ydot0]</literal>のどちらかです.
-
<literal>g(t,y0,ydot0)</literal>は 0 に等しい必要があります.
-
<literal>ydot0</literal>の推定値のみが既知の場合,
-
<literal>info(7)=1</literal> と設定してください.
-
</para>
-
<variablelist>
-
<varlistentry>
-
<term>y0</term>
-
<listitem>
-
<para>初期条件の実数列ベクトル.</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>ydot0</term>
-
<listitem>
-
<para>
-
<literal>t0</literal>における<literal>y</literal>の時間微分の
-
実数列ベクトル (推定値の場合もある).
-
</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>t0</term>
-
<listitem>
-
<para>初期時間を表す実数.</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>t</term>
-
<listitem>
-
<para>実数のスカラーおよびベクトル. 解を求める時間を指定します.
-
<literal>info(2)=1</literal>を設定することにより,
-
dasslの各ステップ点での解を得ることができることに注意してください.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>atol, rtol</term>
-
<listitem>
-
<para>
-
実数スカラーまたは<literal>y</literal>と同じ大きさの列ベクトル.
-
<literal>atol, rtol</literal> はそれぞれ絶対および相対許容誤差を指定します.
-
ベクトルの場合,<literal>y</literal>の各要素毎に許容誤差が指定されます.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>res</term>
-
<listitem>
-
<para>
-
<link linkend="external">外部</link> (関数またはリストまたは文字列).
-
<literal>g(t,y,ydot)</literal>の値を計算. 以下のようにします :
-
</para>
-
<itemizedlist>
-
<listitem>
-
<para>Scilab関数.</para>
-
<para>呼び出し手順は
-
<literal>[r,ires]=res(t,y,ydot)</literal> とする必要があり,
-
<literal>res</literal> は残差
-
<literal>r=g(t,y,ydot)</literal> と誤差フラグ<literal>ires</literal>を
-
返す必要があります.
-
<literal>res</literal>が<literal>r</literal>の計算に
-
成功した場合には<literal>ires = 0</literal>,
-
<literal>(t,x,xdot)</literal>のローカルな残差が定義されない
-
場合には <literal>=-1</literal>,
-
パラメータが許容範囲外の場合は <literal>=-2</literal> となります.
-
</para>
-
</listitem>
-
<listitem>
-
<para>リスト.</para>
-
<para>外部ルーチンのこの形式は
-
関数にパラメータを指定する際に使用されます.
-
以下のような形式とします:
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(res,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>
-
ただし,ここで関数<literal>res</literal>の呼び出し手順は
-
以下のようになります
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
-r=res(t,y,ydot,x1,x2,...)
- ]]></programlisting>
-
+ <screen><![CDATA[
+r = res(t,y,ydot,x1,x2,...)
+]]></screen>
<para>
-
この場合も<literal>res</literal> は
-
<literal>(t,y,ydot,x1,x2,...)</literal>の関数として<literal>r=g(t,y,ydot)</literal>
-
を返します.
-
</para>
-
</listitem>
-
<listitem>
-
<para>文字列.</para>
-
<para>
-
ScilabにリンクされるCまたはFortranサブルーチンの名前を指します.
-
</para>
-
<para>Cの場合,呼び出し手順は次のようになります:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void res(double *t, double y[], double yd[], double r[],
int *ires, double rpar[], int ipar[])
- ]]></programlisting>
-
+]]></screen>
<para>Fortranの場合,呼び出し手順は次のようになります:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine res(t,y,yd,r,ires,rpar,ipar)
double precision t, y(*),yd(*),r(*),rpar(*)
integer ires,ipar(*)
- ]]></programlisting>
-
+]]></screen>
<para>
-
<literal>rpar</literal> と <literal>ipar</literal> 配列は存在しますが,
-
使用できません.
-
</para>
-
</listitem>
-
</itemizedlist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>jac</term>
-
<listitem>
-
<para>外部 (関数またはリストまたは文字列).
-
指定したパラメータ<literal>cj</literal>の値を用いて
-
<literal>dg/dy+cj*dg/dydot</literal> の値を計算します.
-
</para>
-
<itemizedlist>
-
<listitem>
-
<para>Scilab関数.</para>
-
<para>その呼び出し手順は
-
<literal>r=jac(t,y,ydot,cj)</literal> とする必要があり,
-
<literal>jac</literal>関数は
-
<literal>r=dg(t,y,ydot)/dy+cj*dg(t,y,ydot)/dydot</literal>を
-
返す必要があります.
-
ただし, <literal>cj</literal>は実数スカラーです
-
</para>
-
</listitem>
-
<listitem>
-
<para>リスト.</para>
-
<para>以下のようにする必要があります</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
list(jac,x1,x2,...)
- ]]></programlisting>
-
+]]></screen>
<para>
-
ただしこの場合の関数<literal>jac</literal>の呼び出し手順は
-
以下のようになります
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
-r=jac(t,y,ydot,cj,x1,x2,...)
- ]]></programlisting>
-
+ <screen><![CDATA[
+r = jac(t,y,ydot,cj,x1,x2,...)
+]]></screen>
<para>
-
この場合も<literal>jac</literal> は
-
<literal>(t,y,ydot,cj,x1,x2,...)</literal>の関数として
-
<literal>dg/dy+cj*dg/dydot</literal>を返します.
-
</para>
-
</listitem>
-
<listitem>
-
<para>文字列.</para>
-
<para>
-
ScilabにリンクされるFortranサブルーチンの名前を指す
-
必要があります
-
</para>
-
<para>Cの場合,呼び出し手順は次のようになります:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
void jac(double *t, double y[], double yd[], double pd[],
double *cj, double rpar[], int ipar[])
- ]]></programlisting>
-
+]]></screen>
<para>Fortranの場合, 呼び出し手順は次のようになります:</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
subroutine jac(t,y,yd,pd,cj,rpar,ipar)
double precision t, y(*),yd(*),pd(*),cj,rpar(*)
integer ipar(*)
- ]]></programlisting>
-
+]]></screen>
</listitem>
-
</itemizedlist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info</term>
-
<listitem>
-
<para>
-
<literal>7</literal> 個の要素を含むオプションのリスト.
-
デフォルト値は list([],0,[],[],[],0,0);
-
</para>
-
<variablelist>
-
<varlistentry>
-
<term>info(1)</term>
-
<listitem>
-
<para>
-
<literal>g</literal>を評価できる時間の最大値
-
を指定する実数スカラー.
-
または,時間に制限を課さない場合には,空の行列
-
<literal>[]</literal>
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(2)</term>
-
<listitem>
-
<para>
-
<literal>dassl</literal> が
-
その中間的な計算値を返す(<literal>flag=1</literal>),
-
もしくはユーザが指定した時間においてのみ値を返す
-
(<literal>flag=0</literal>)
-
かどうかを指定するフラグ.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(3)</term>
-
<listitem>
-
<para>
-
<literal>2</literal> 個の要素のベクトルで,
-
<literal>jac</literal>; <literal>r(i - j + ml + mu + 1,j) =
-
"dg(i)/dy(j)+cj*dg(i)/dydot(j)"
-
</literal>
-
により計算された
-
帯行列の定義<literal>[ml,mu]</literal>を指定します.
-
<literal>jac</literal> が完全な行列を返す場合,
-
<literal>info(3)=[]</literal> を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(4)</term>
-
<listitem>
-
<para>実数スカラーで,
-
ステップの大きさの最大値を指定します.
-
制限をしない場合には,<literal>info(4)=[]</literal>
-
を指定してください.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(5)</term>
-
<listitem>
-
<para>実数スカラーで,ステップの大きさの初期値を指定します.
-
指定しない場合には,<literal>info(5)=[]</literal>
-
を指定してください.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(6)</term>
-
<listitem>
-
<para>
-
解が非負であることが既知の場合に<literal>info(6)=1</literal>に設定,
-
そうでない場合は <literal>info(6)=0</literal> に設定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>info(7)</term>
-
<listitem>
-
<para>
-
<literal>ydot0</literal>が単なる推定値の場合は
-
<literal>info(7)=1</literal>,
-
<literal>g(t0,y0,ydot0)=0</literal>の場合は
-
<literal>info(7)=0</literal> を指定します.
-
</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>hd</term>
-
<listitem>
-
<para>
-
実数ベクトルで,<literal>dassl</literal>のコンテキストを保存し,
-
積分を再開することを可能にします.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>r</term>
-
<listitem>
-
<para>実数行列. 各列はベクトル[t;x(t);xdot(t)]です.
-
ただし t は解を計算した時間です.
-
</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</refsection>
-
<refsection>
-
<title>説明</title>
-
<para>
-
<literal>dassl</literal>関数は代数微分方程式を積分し,
-
指定した時間における <literal>y</literal> の
-
値を返します.
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
-g(t,y,ydot)=0
+ <screen><![CDATA[
+g(t,y,ydot) = 0
y(t0)=y0 and ydot(t0)=ydot0
- ]]></programlisting>
-
+]]></screen>
</refsection>
-
<refsection>
-
<title>例</title>
-
<programlisting role="example"><![CDATA[
function [r,ires]=chemres(t,y,yd)
r=[-0.04*y(1)+1d4*y(2)*y(3)-yd(1)
y(1)+y(2)+y(3)-1];
ires=0
endfunction
+
function pd=chemjac(x,y,yd,cj)
pd=[-0.04-cj , 1d4*y(3) , 1d4*y(2);
0.04 ,-1d4*y(3)-2*3d7*y(2)-cj ,-1d4*y(2);
1 , 1 , 1 ]
endfunction
-y0=[1;0;0];
-yd0=[-0.04;0.04;0];
-t=[1.d-5:0.02:.4,0.41:.1:4,40,400,4000,40000,4d5,4d6,4d7,4d8,4d9,4d10];
-y=dassl([y0,yd0],0,t,chemres);
-info=list([],0,[],[],[],0,0);
-info(2)=1;
-y1=dassl([y0,yd0],0,4d10,chemres,info);
-y2=dassl([y0,yd0],0,4d10,chemres,chemjac,info);
+
+y0 = [1;0;0];
+yd0 = [-0.04;0.04;0];
+t = [1.d-5:0.02:.4,0.41:.1:4,40,400,4000,40000,4d5,4d6,4d7,4d8,4d9,4d10];
+y = dassl([y0,yd0],0,t,chemres);
+
+info = list([],0,[],[],[],0,0);
+info(2) = 1;
+y1 = dassl([y0,yd0],0,4d10,chemres,info);
+y2 = dassl([y0,yd0],0,4d10,chemres,chemjac,info);
+
//パラメータとしてオプション引数を使用
//-----------------------------------
function [r,ires]=chemres(t,y,yd ,a,b,c)
a ,-b*y(3)-2*c*y(2)-cj ,-b*y(2);
1 , 1 , 1 ]
endfunction
-y3=dassl([y0,yd0],0,t,list(chemres,0.04,1d4,3d7),list(chemjac,0.04,1d4,3d7));
+y3 = dassl([y0,yd0],0,t,list(chemres,0.04,1d4,3d7),list(chemjac,0.04,1d4,3d7));
+
// Cコードを使用
//------------
// - Cコードを作成
cd TMPDIR
-rescode=['void chemres(double *t, double y[], double yd[], double r[], int *ires, double rpar[], int ipar[])'
- ' {'
- ' r[0] = -0.04*y[0]+1.0e4*y[1]*y[2] -yd[0];'
- ' r[1] = 0.04*y[0]-1.0e4*y[1]*y[2]-3.0e7*y[1]*y[1]-yd[1];'
- ' r[2] = y[0]+y[1]+y[2]-1;'
- ' *ires = 0;'
- ' }'];
-jaccode=['void chemjac(double *t, double y[], double yd[], double pd[], double *cj, double rpar[], int ipar[])'
- ' {'
- ' /* first column*/'
- ' pd[0] = -0.04-*cj;'
- ' pd[1] = 0.04;'
- ' pd[2] = 1.0;'
- ' /* second column*/'
- ' pd[3] = 1.0e4*y[2];'
- ' pd[4] = -1.0e4*y[2]-2*3.0e7*y[1]-*cj;'
- ' pd[5] = 1.0;'
- ' /* third column*/'
- ' pd[6] = 1.0e4*y[1];'
- ' pd[7] = -1.0e4*y[1];'
- ' pd[8] = 1.0;'
- ' }'];
+rescode = ['void chemres(double *t, double y[], double yd[], double r[], int *ires, double rpar[], int ipar[])'
+ ' {'
+ ' r[0] = -0.04*y[0]+1.0e4*y[1]*y[2] -yd[0];'
+ ' r[1] = 0.04*y[0]-1.0e4*y[1]*y[2]-3.0e7*y[1]*y[1]-yd[1];'
+ ' r[2] = y[0]+y[1]+y[2]-1;'
+ ' *ires = 0;'
+ ' }'];
+jaccode = ['void chemjac(double *t, double y[], double yd[], double pd[], double *cj, double rpar[], int ipar[])'
+ ' {'
+ ' /* first column*/'
+ ' pd[0] = -0.04-*cj;'
+ ' pd[1] = 0.04;'
+ ' pd[2] = 1.0;'
+ ' /* second column*/'
+ ' pd[3] = 1.0e4*y[2];'
+ ' pd[4] = -1.0e4*y[2]-2*3.0e7*y[1]-*cj;'
+ ' pd[5] = 1.0;'
+ ' /* third column*/'
+ ' pd[6] = 1.0e4*y[1];'
+ ' pd[7] = -1.0e4*y[1];'
+ ' pd[8] = 1.0;'
+ ' }'];
mputl([rescode;jaccode],'mycode.c') //Cファイルを作成
// - コンパイル
ilib_for_link(['chemres','chemjac'],'mycode.c',[],'c','','loader.sce');//compile
// - Scilabにリンク
exec(loader.sce') //incremental linking
// - dasslをコール
-y4=dassl([y0,yd0],0,t,'chemres','chemjac');
+y4 = dassl([y0,yd0],0,t,'chemres','chemjac');
]]></programlisting>
-
</refsection>
-
<refsection role="see also">
-
<title>参照</title>
-
<simplelist type="inline">
-
<member>
-
<link linkend="ode">ode</link>
-
</member>
-
<member>
-
<link linkend="dasrt">dasrt</link>
-
</member>
-
<member>
-
<link linkend="impl">impl</link>
-
</member>
-
<member>
-
<link linkend="call">call</link>
-
</member>
-
<member>
-
<link linkend="link">link</link>
-
</member>
-
<member>
-
<link linkend="external">external</link>
-
</member>
-
</simplelist>
-
</refsection>
-
-</refentry>
-
+</refentry>
\ No newline at end of file
<?xml version="1.0" encoding="UTF-8"?>
-
<!--
* Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
* Copyright (C) 2008 - INRIA
* along with this program.
*
-->
-
<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:id="diff" xml:lang="ja">
-
<refnamediv>
-
<refname>diff</refname>
-
<refpurpose>差分および離散微分</refpurpose>
-
</refnamediv>
-
<refsynopsisdiv>
-
<title>呼び出し手順</title>
-
<synopsis>
-
- y=diff(x)
-
- y=diff(x [,n [,dim]])
-
+ y = diff(x)
+ y = diff(x, n)
+ y = diff(x, n, dim)
</synopsis>
-
</refsynopsisdiv>
-
<refsection>
-
<title>引数</title>
-
<variablelist>
-
<varlistentry>
-
<term>x</term>
-
<listitem>
-
<para>ベクトルまたは行列 (実数, 複素数, 疎行列 または多項式)</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>n</term>
-
<listitem>
-
<para>整数,差分の次数</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>dim</term>
-
<listitem>
-
<para>整数または"r","c" および "*"を値とする文字列</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>y</term>
-
<listitem>
-
<para>スカラーまたはベクトル</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</refsection>
-
<refsection>
-
<title>説明</title>
-
<para>
-
<literal>y=diff(x)</literal> は差分関数
-
y=x(2:$)-x(1:$-1)を計算します.
-
</para>
-
<para>
-
<literal>diff(x,n,dim)</literal> は,
-
次元<literal>dim</literal>方向の<literal>n</literal> 次
-
差分関数です.
-
</para>
-
<para>
-
<literal>diff(x,n,"*")</literal> は
-
diff(x(:),n)と等価です.
-
</para>
-
<para>
-
<literal>n</literal>のデフォルト値は 1です.
-
<literal>dim</literal> のデフォルト値は <literal>"*"</literal>です.
-
</para>
-
<para>
-
<literal>dim='r'</literal> は <literal>dim=1</literal>と
-
等価で,
-
<literal>dim='c'</literal> は
-
<literal>dim=2</literal>と等価です.
-
</para>
-
</refsection>
-
<refsection>
-
<title>例</title>
-
<programlisting role="example"><![CDATA[
-v=(1:8)^3;
+v = (1:8)^3;
diff(v)
diff(v,3)
-A=[(1:8)^2
- (1:8)^3
- (1:8)^4];
+A = [(1:8)^2
+ (1:8)^3
+ (1:8)^4];
diff(A,3,2)
+
//approximate differentiation
-step=0.001
-t=0:step:10;
-y=sin(t);
-dy=diff(sin(t))/step; //approximate differentiation of sine function
+step = 0.001
+t = 0:step:10;
+y = sin(t);
+dy = diff(sin(t))/step; //approximate differentiation of sine function
norm(dy-cos(t(1:$-1)),%inf)
]]></programlisting>
-
</refsection>
-
<refsection role="see also">
-
<title>参照</title>
-
<simplelist type="inline">
-
<member>
-
<link linkend="interp">interp</link>
-
</member>
-
<member>
-
<link linkend="interp2d">interp2d</link>
-
</member>
-
<member>
-
<link linkend="splin">splin</link>
-
</member>
-
<member>
-
<link linkend="eval_cshep2d">eval_cshep2d</link>
-
</member>
-
<member>
-
<link linkend="numderivative">numderivative</link>
-
</member>
-
</simplelist>
-
</refsection>
-
-</refentry>
-
+</refentry>
\ No newline at end of file
<?xml version="1.0" encoding="UTF-8"?>
-
<!--
* Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
* Copyright (C) 2008 - INRIA
* along with this program.
*
-->
-
<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:id="impl" xml:lang="ja">
-
<refnamediv>
-
<refname>impl</refname>
-
<refpurpose>微分代数方程式</refpurpose>
-
</refnamediv>
-
<refsynopsisdiv>
-
<title>呼出し手順</title>
-
- <synopsis>y=impl([type],y0,ydot0,t0,t [,atol, [rtol]],res,adda [,jac])</synopsis>
-
+ <synopsis>y = impl([type], y0, ydot0, t0, t [,atol, [rtol]], res, adda [,jac])</synopsis>
</refsynopsisdiv>
-
<refsection>
-
<title>引数</title>
-
<variablelist>
-
<varlistentry>
-
<term>y0,ydot0</term>
-
<listitem>
-
<para>実数のベクトルまたは行列 (初期条件).</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>t0</term>
-
<listitem>
-
<para>実数のスカラー (初期時間).</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>t</term>
-
<listitem>
-
<para>実数ベクトル (解を計算する時刻).</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>res,adda</term>
-
<listitem>
-
<para>外部 (関数または文字列またはリスト).</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>type</term>
-
<listitem>
-
<para>
-
文字列 <literal>'adams'</literal> または
-
<literal>'stiff'</literal>
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>atol,rtol</term>
-
<listitem>
-
<para>
-
<literal>y</literal>と同じ大きさの実数スカラーまたは実数ベクトル.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>jac</term>
-
<listitem>
-
<para>外部 (関数または文字列またはリスト).</para>
-
</listitem>
-
</varlistentry>
-
</variablelist>
-
</refsection>
-
<refsection>
-
<title>説明</title>
-
<para>陰的な線形微分方程式の解</para>
-
<para>A(t,y) dy/dt=g(t,y), y(t0)=y0</para>
-
<para>
-
<literal>t0</literal> は初期時刻, <literal>y0</literal>
-
は初期条件,
-
<literal>t0</literal>における<literal>y</literal>の時間微分
-
ベクトル <literal>ydot0</literal> も指定する必要があります.
-
入力 <literal>res</literal>は,
-
<link linkend="external">外部ルーチン</link> ,すなわち,
-
指定された構文を有する関数,また規定の呼び出し手順を有する
-
FortranサブルーチンまたはC関数の名前(文字列),
-
またはリストです.
-
</para>
-
<para>
-
<literal>res</literal> が関数の場合, その構文は
-
以下のようにする必要があります:
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r = res(t,y,ydot)
- ]]></programlisting>
-
+]]></screen>
<para>
-
ただし,<literal>t</literal>は,実数ベクトル (時間) そして
-
<literal>y</literal> および <literal>ydot</literal> は実数ベクトル
-
(状態量および状態量の微分)です. この関数は
-
<literal>r=g(t,y)-A(t,y)*ydot</literal>を返す必要があります.
-
</para>
-
<para>
-
<literal>res</literal>が文字列の場合,
-
FortranサブルーチンまたはC関数の名前を指します.
-
この例については,
-
<literal>SCI/modules/differential_equations/sci_gateway/fortran/Ex-impl.f</literal>を参照ください.
-
</para>
-
<para>
-
<literal>res</literal> をリストとすることもできます:
-
<link linkend="ode">ode</link>のヘルプを参照ください.
-
</para>
-
<para>
-
入力 <literal>adda</literal> も外部ルーチンですi.
-
</para>
-
<para>
-
<literal>adda</literal> が関数の場合, その構文は以下のようにする必要があります:
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
r = adda(t,y,p)
- ]]></programlisting>
-
+]]></screen>
<para>
-
これは<literal>r=A(t,y)+p</literal>を返します. ただし,
-
<literal>p</literal> は
-
<literal>A(t,y)</literal>に加算される行列です.
-
</para>
-
<para>
-
<literal>adda</literal> が文字列の場合,
-
FortranサブルーチンまたはC関数の名前を指します.
-
この例については,
-
<literal>SCI/modules/differential_equations/sci_gateway/fortran/Ex-impl.f</literal>を参照ください.
-
</para>
-
<para>
-
<literal>adda</literal> をリストとすることもできます:
-
<literal>ode</literal>のヘルプを参照ください.
-
</para>
-
<para>
-
入力 <literal>jac</literal> も外部ルーチンです.
-
</para>
-
<para>
-
<literal>jac</literal> が関数の場合, その構文は以下のようにする必要があります:
-
</para>
-
- <programlisting role="no-scilab-exec"><![CDATA[
+ <screen><![CDATA[
j = jac(t,y,ydot)
- ]]></programlisting>
-
+]]></screen>
<para>これは
-
<literal>r=g(t,y)-A(t,y)*ydot</literal>の
-
<literal>y</literal>に関するヤコビアンを返します.
-
</para>
-
<para>
-
<literal>jac</literal> が文字列の場合,
-
FortranサブルーチンまたはC関数の名前を指します.
-
この例については,
-
<literal>SCI/modules/differential_equations/sci_gateway/fortran/Ex-impl.f</literal>を参照ください.
-
</para>
-
<para>
-
<literal>jac</literal> をリストとすることも可能です:
-
<link linkend="ode">ode</link>のヘルプを参照ください.
-
</para>
-
</refsection>
-
<refsection>
-
<title>例</title>
-
<programlisting role="example"><![CDATA[
-y=impl([1;0;0],[-0.04;0.04;0],0,0.4,'resid','aplusp');
+y = impl([1;0;0], [-0.04;0.04;0], 0, 0.4, 'resid', 'aplusp');
// Using hot restart
-//[x1,w,iw]=impl([1;0;0],[-0.04;0.04;0],0,0.2,'resid','aplusp');
+// [x1, w, iw] = impl([1;0;0], [-0.04;0.04;0], 0, 0.2, 'resid', 'aplusp');
// hot start from previous call
-//[x1]=impl([1;0;0],[-0.04;0.04;0],0.2,0.4,'resid','aplusp',w,iw);
-//maxi(abs(x1-x))
+// x1 = impl([1;0;0], [-0.04;0.04;0], 0.2, 0.4, 'resid', 'aplusp', w, iw);
+// maxi(abs(x1-x))
]]></programlisting>
-
</refsection>
-
<refsection role="see also">
-
<title>参照</title>
-
<simplelist type="inline">
-
<member>
-
<link linkend="dassl">dassl</link>
-
</member>
-
<member>
-
<link linkend="ode">ode</link>
-
</member>
-
<member>
-
<link linkend="external">external</link>
-
</member>
-
</simplelist>
-
</refsection>
-
-</refentry>
-
+</refentry>
\ No newline at end of file
<?xml version="1.0" encoding="UTF-8"?>
-
<!--
* Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
* Copyright (C) 2008 - INRIA
* along with this program.
*
-->
-
<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink"
xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml"
xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook"
xmlns:scilab="http://www.scilab.org" xml:id="int2d" xml:lang="ja">
<refnamediv>
-
<refname>int2d</refname>
-
<refpurpose>求積法および立体求積法により2次元積分を定義</refpurpose>
-
</refnamediv>
-
<refsynopsisdiv>
-
<title>呼び出し手順</title>
-
<synopsis>
[I, err] = int2d(X, Y, f)
[I, err] = int2d(X, Y, f, params)
</synopsis>
-
</refsynopsisdiv>
-
<refsection>
-
<title>引数</title>
-
<variablelist>
-
<varlistentry>
-
<term>X</term>
-
<listitem>
-
<para>
-
N個の三角形の頂点の横座標軸の値を有する3 行 N 列の配列.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>Y</term>
-
<listitem>
-
<para>
-
N個の三角形の頂点の縦座標軸の値を有する3 行 N 列の配列.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>f</term>
-
<listitem>
-
<para>
-
被積分関数 <literal>f(u,v)</literal>を定義する
-
外部 (関数またはリストまたは文字列)
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>params</term>
-
<listitem>
-
<para>
-
実数ベクトル <literal>[tol, iclose, maxtri, mevals,
-
iflag]
-
</literal>
-
.デフォルト値は <literal>[1.d-10, 1, 50, 4000,
-
1]
-
</literal>
-
.
-
</para>
-
<variablelist>
-
<varlistentry>
-
<term>tol</term>
-
<listitem>
-
<para>
-
誤差の境界の指定値.
-
<literal>iflag=0</literal>の場合, <literal>tol</literal> は
-
相対誤差の境界として解釈されます;
-
<literal>iflag=</literal>1の場合,
-
絶対誤差の境界と解釈されます.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>iclose</term>
-
<listitem>
-
<para>
-
LQM0 または LQM1の選択を定義する積分パラメータ.
-
<literal>iclose=1</literal>の場合, LQM1 が使用されます.
-
<literal>iclose</literal>がその他の値の場合,
-
LQM0 が使用されます.
-
LQM0 は三角形の内部の点における値のみ関数の値を使用します.
-
LQM1 は通常 LQM0 より正確ですが, 三角形の境界点を含む
-
より多くの点で被積分関数を評価します.
-
通常,被積分関数が三角形の境界において特異とならない限り,
-
LQM1 を使用するのが良いでしょう.
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>maxtri</term>
-
<listitem>
-
<para>
-
領域の最後の三角形分割における三角形の最大数
-
</para>
-
</listitem>
-
</varlistentry>
-
<varlistentry>
-
<term>mevals</term>
-
<listitem>
-
<para>
-
&nbs
projects / scilab.git / commitdiff