#COPYRIGHT NOTICE: Copyright (c) 1995 by Sasha Cyganowski. Permission is #granted to anyone to use, modify, or redistribute this software freely, #subject to the following restrictions: 1. The author accepts no #responsibility for any consequences of this software and makes no guarantee #that this software is free of defects. 2. The origin of this software must #not be misrepresented, either by explicit claim or by omission. #3. This notice and the copyright must be included with all copies or #modified versions of this software. 4. This software may not be included or #redistributed as part of any package to be sold for profit unlesss the author #has given explicit written permission to do so. # Sasha Cyganowski, School of Computing and Mathematics, Deakin University, #Waurn Ponds, VIC, Australia. Internet: sash@deakin.edu.au # Telephone: 61 52 271382. # This file contains the function definitions and help files for the package #"stochastic", developed by Sasha Cyganowski. The package contains functions #which can be used to find explicit solutions of Stochastic Differential #Equations (SDEs), construct numerical schemes up to strong order 2 and weak #order 3, check for commutative noise of the first and second kind, and #convert SDEs into their coloured noise form. Other useful routines are also #included. # The general schemes, solutions, and conditions from which the code # was developed can be found in Kloeden, P.E, Platen, E.: Numerical # Solution of Stochastic Differential Equations, (Springer-Verlag, 1992). # Note: This software is known to work for MapleV version 3 (for Windows) # and MapleV version 3 (for Unix). `help/text/stochastic` := TEXT( ` `, `HELP FOR: Introduction to the stochastic package`, ` `, `CALLING SEQUENCE:`, ` (args)`, ` stochastic[](args)`, ` `, `SYNOPSIS:`, `- The stochastic package contains routines useful for finding explicit`, ` solutions of Stochastic Differential Equations (SDEs), and routines`, ` useful for constructing numerical schemes up to strong order 2 and`, ` weak order 3.`, ` Other useful features include a routine which converts SDEs with white`, ` noise into coloured noise form, and routines which check whether`, ` an SDE has commutative noise of the first and/or second kind. Other`, ` useful routines are also available. In particular, the operator`, ` procedures L0 and LJ are included so that users can easily construct`, ` numerical schemes other than those already available.`, ` For more information on a particular function type the command`, ` ?.`, ` The author of the stochastic package is Sasha Cyganowski, from Deakin`, ` University, Waurn Ponds, Australia, 3217. Comments are welcomed via`, ` e-mail, which can be sent to the folllowing address: sash@deakin.edu.au`, ` The general schemes, solutions, and conditions used in the coding of the`, ` stochastic package can be found in Kloeden, P.E, Platen, E.: Numerical`, ` Solution of Stochastic Differential Equations, (Springer-Verlag, 1992).`, ` The package is also described in Cyganowski, S.O., A MAPLE Package for`, ` Stochastic Differential Equations, Computational Techniques and`, ` Applications: CTAC95, editors A. Easton, R. May, World Scientific,`, ` (to appear).`, ` `, `- To use a stochastic function, either define that function alone using the`, ` command with(stochastic, ), or define all stochastic functions`, ` using the command with(stochastic). Alternatively, invoke the function`, ` using the long form stochastic[]. This long form notation is`, ` necessary whenever there is a conflict between a package function name`, ` and another function used in the same session.`, ` `, `- The functions available are:`, ` `, ` linear MLJ explicit L0 Euler Milstein`, ` Taylor1hlf SLO correct Taylor2 wkeuler wktay2`, ` reducible wktay3 colour comm1 comm2 LJ`, ` `, `- As an example, to find the explicit solution of an SDE with drift`, ` coefficient 1/2*a^2*x and diffusion coefficient a*x, use`, ` `, ` with(stochastic,explicit); explicit(1/2*a^2*x,a*x);`, ` `, `- Note the stochastic numerical routines require the drift and diffusion`, ` coefficients of an SDE to be entered in the variables x[N] and t, where`, ` x[N] are the state variables of the N-dimensional SDE and t denotes time.`, ` The routines linear, reducible, and explicit require the drift and`, ` duffusion coefficients to be entered in the variables x and t.`, ` `, `SEE ALSO: with` ): stochastic[linear]:=proc(a:algebraic,b:algebraic) local temp1,alpha,beta,gamma,delta,fundsoln,fundsoln2,soln,default1,default2, default3; if diff(a,x,x) <> 0 or diff(b,x,x) <> 0 then ERROR(`SDE not linear, try a reducible procedure`) else alpha := diff(a,x); alpha := subs(t = s,alpha); beta := diff(b,x); beta := subs(t = s,beta); if diff(beta,s) = 0 then temp1 := beta*W; else temp1:=Int(beta,W = 0 .. t); fi; gamma := coeff(a,x,0); gamma := subs(t = s,gamma); delta := coeff(b,x,0); delta := subs(t = s,delta); fundsoln := exp(int(alpha-1/2*beta^2,s = 0 .. t)+temp1); fundsoln2 := subs(t = s,fundsoln); if beta = 0 then soln := fundsoln*(X[0]+int(1/fundsoln2*(gamma-beta*delta), s = 0 .. t)+Int(1/fundsoln2*delta,W = 0 .. t)) else soln := fundsoln*(X[0]+Int(1/fundsoln2*(gamma-beta*delta), s = 0 .. t)+Int(1/fundsoln2*delta,W = 0 .. t)) fi; default1 := Int(0,W = 0 .. t) = 0; default2 := Int(0,W = 0 .. s) = 0; default3 := Int(0,s = 0 .. t) = 0; soln := X[t] = subs(default1,default2,default3,soln); fi; end: `help/text/linear` := TEXT( `FUNCTION: stochastic[linear] - returns the explicit solution of a linear`, `Stochastic Differential Equation (SDE).`, ` `, `CALLING SEQUENCE: linear(a,b);`, ` `, `PARAMETERS: a - algebraic, given in the variables x and t.`, ` b - algebraic, given in the variables x and t.`, ` `, `SYNOPSIS: `, `- The call linear(a,b); returns the explicit solution for a linear`, ` Stochastic Differential Equation (SDE) with drift coefficient a, and`, ` diffusion coefficient b. The output consists of the variables X[t], X[0],`, ` and W. X[t] denotes the explicit solution, X[0] denotes the initial value`, ` of the solution, W denotes a standard Wiener process, and t denotes time.`, ` If the SDE is not of linear form a suitable error message is returned.`, ` `, `- This routine is used by the routine stochastic[explicit] and, in general,`, ` is not intended to be used on its own.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[linear].`, ` `, `EXAMPLES: `, ` `, `linear(-x,2);`, ` `, ` / t \`, ` | / |`, ` | | 2 |`, ` X[t] = exp(- t) |X[0] + | -------- dW|`, ` | | exp(- s) |`, ` | / |`, ` \ 0 /`, ` `, `SEE ALSO: stochastic, stochastic[reducible], stochastic[explicit] ` ): stochastic[reducible]:=proc(a:algebraic,b:algebraic) local beta,temp1,h,temp3,alpha,soln,soln1; h := int(1/b,x); temp1 := alpha*b*h+1/2*b*simplify(diff(b,x)); temp1 = a; alpha := simplify(solve(",alpha)); beta := alpha*h; if diff(alpha,x) = 0 then if alpha=0 then soln:=h=subs(x=X[0],h)+W; X[t]=simplify(solve(soln,x)); else soln1 := h = exp(alpha*t)*subs(x = X[0],h)+exp(alpha*t)*Int(exp(-alpha *s),W = 0 .. t); X[t] = solve(soln1,x); fi elif diff(beta,x) = 0 then X[t]=simplify(solve(h = beta*t+W+subs(x = X[0],h),x)); else ERROR(`non-linear SDE not reducible`) fi end: `help/text/reducible` := TEXT( `FUNCTION: stochastic[reducible] - returns the explicit solution of a`, `reducible Stochastic Differential Equation (SDE).`, ` `, `CALLING SEQUENCE: reducible(a,b);`, ` `, `PARAMETERS: a - algebraic, given in the variables x and t.`, ` b - algebraic, given in the variables x and t.`, ` `, `SYNOPSIS: `, `- The call reducible(a,b); returns the explicit solution for a reducible`, ` Stochastic Differential Equation (SDE) with drift coefficient a, and`, ` diffusion coefficient b. The output consists of the variables X[t], X[0],`, ` and W. X[t] denotes the explicit solution, X[0] denotes the initial`, ` value of the solution, W denotes a standard Wiener process, and t denotes`, ` time. If the SDE is not of reducible form a suitable error message is`, ` returned.`, ` `, `- This routine is used by the routine stochastic[explicit] and, in general,`, ` is not intended to be used on its own.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be called via the call`, ` stochastic[reducible].`, ` `, `EXAMPLES:`, ` `, `>reducible(1,2*sqrt(x));`, ` `, ` 1/2 2`, ` X[t] = 2 X[0] W + W + X[0]`, ` `, `SEE ALSO: stochastic, stochastic[linear],stochastic[explicit] ` ): stochastic[explicit]:= proc(a:algebraic,b:algebraic) if diff(a,x,x) = 0 and diff(b,x,x) = 0 then linear (a,b) else reducible(a,b) fi end: `help/text/explicit` := TEXT( `FUNCTION: stochastic[explicit] - returns the explicit solution of a`, ` Stochastic Differential Equation (SDE).`, ` `, `CALLING SEQUENCE: explicit(a,b);`, ` `, `PARAMETERS: a - algebraic, given in the variables x and t.`, ` b - algebraic, given in the variables x and t.`, ` `, `SYNOPSIS: `, `- The call explicit(a,b); returns the explicit solution for a Stochastic`, ` Differential Equation (SDE) with drift coefficient a, and diffusion`, ` coefficient b. The output consists of the variables X[t], X[0], and W.`, ` X[t] denotes the explicit solution, X[0] denotes the initial value of the`, ` solution, and W denotes a standard Wiener process. The user is returned`, ` with a suitable error message if no known explicit solution exists.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[explicit].`, ` `, `EXAMPLES: `, ` `, `>explicit(1/2*a^2*x,a*x);`, ` `, ` X[t]=exp(aW)X[0]`, ` `, `SEE ALSO: stochastic, stochastis[linear], stochastic[reducible]` ): stochastic[LJ]:= proc(X:algebraic,b:list(list(algebraic)),j:integer) sum('op(j,op(k,b))* diff(X,x[k])','k' = 1 .. nops(b)) end: `help/text/LJ` := TEXT( `FUNCTION: stochastic[LJ] - applies the partial differential operator LJ, as`, ` described in Kloeden, P.E., Platen, E.: Numerical Solution of Stochastic`, ` Differential Equations, pg 339, (Springer-Verlag 1992).`, ` `, `CALLING SEQUENCE: LJ(X,[[b11,..,b1M],..,[bN1,..,bNM]],j);`, ` `, `PARAMETERS: X - algebraic, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` the variables x[N] and t.`, ` j - integer.`, ` `, `SYNOPSIS: `, `- The call LJ(X,[[b11,..,b1M],..,[bN1,..,bNM]],j); computes the application`, ` of the partial differential operator LJ to X, as described in`, ` Kloeden, P.E., Platen, E: Numerical Solution of Stochastic Differential`, ` Equations, pg 339, (Springer-Verlag 1992).`, ` [b11,..,b1M],..,[bN1,..,bNM] denotes the`, ` diffusion matrix of dimension N*M, where M is the dimension of the`, ` Wiener process and N is the dimension of the SDE. j=1,..,M denotes the`, ` "current" dimension of the Wiener process.`, ` The output variables are consistent with the variables used as input.`, ` `, `- This routine is used by the following routines: stochastic[Euler],`, ` stochastic[Milstein], stochastic[Taylor1.5], stochastic[Taylor2],`, ` stochastic[wkeuler], and stochastic[wktay2].`, ` In general, it is not intended for use on its own.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[LJ]`, ` `, `EXAMPLES: `, ` `, `>LJ(x[2],[[r,0],[0,r]],2);`, ` `, ` r `, ` `, `>LJ(x[2],[[r,0],[0,r]],1);`, ` `, ` 0 `, ` `, `SEE ALSO: stochastic[L0], stochastic[MLJ], stochastic[SL0],`, ` stochastic[Euler], stochastic[Milstein], stochastic[Taylor1hlf]`, ` stochastic[Taylor2], stochastic[wkeuler], stochastic[wktay2],`, ` stochastic.` ): stochastic[L0]:=proc(X:algebraic,a:list(algebraic),b:list(list(algebraic))) local part1,part2,part3; part1 := diff(X,t); part2 := sum('a[k]*diff(X,x[k])','k' = 1 .. nops(a)); part3 := 1/2*sum( 'sum('sum('op(j,op(k,b))*op(j,op(l,b))*diff(X,x[k],x[l])', 'j' = 1 .. nops(op(1,b)))','k' = 1 .. nops(a))','l' = 1 .. nops(a) ); part1+part2+part3; end: `help/text/L0` := TEXT( `FUNCTION: stochastic[L0] - applies the partial differential operator L0, as`, ` described in Kloeden, P.E., Platen, E.: Numerical Solution of Stochastic`, ` Differential Equations, pg 339, (Springer-Verlag 1992).`, ` `, `CALLING SEQUENCE: L0(X,[a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]);`, ` `, `PARAMETERS: X - algebraic, given in the variables x[N] and t.`, ` a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` the variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call L0(X,[a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]); computes the`, ` application of the partial differential operator L0 to X, as described`, ` in Kloeden, P.E., Platen, E: Numerical Solution of Stochastic`, ` Differential Equations, pg 339, (Springer-Verlag 1992). a1,..,aN`, ` denotes the drift coefficients of the N-dimensional SDE.`, ` [b11,..,b1M],..,[bN1,..,bNM] denotes the diffusion matrix of dimension`, ` N*M.`, ` The output variables are consistent with the variables used as input.`, ` `, `- This routine is used by the following routines: stochastic[Euler],`, ` stochastic[Milstein], stochastic[Taylor1.5], stochastic[wkeuler]`, ` and stochastic[wktay2].`, ` In general, it is not intended for use on its own.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[L0]`, ` `, `EXAMPLES: `, ` `, `>L0(x[2],[x[1],x[2]],[[r,0],[0,r]]);`, ` `, ` x[2] `, ` `, `SEE ALSO: stochastic[LJ], stochastic[MLJ], stochastic[SL0],`, ` stochastic[Euler], stochastic[Milstein], stochastic[Taylor1hlf]`, ` stochastic[wkeuler], stochastic[wktay2], stochastic[wktay3],`, ` stochastic.` ): stochastic[Euler]:=proc(a:list(algebraic),b:list(list(algebraic))) local i,u,soln; for i to nops(a) do soln[i] := Y.i[n+1] = Y.i[n]+L0(x[i],a,b)*Delta[n]+sum('LJ(x[i],b,j)* Delta*W.j[n]','j' = 1 .. nops(op(1,b))); for u to nops(a) do soln[i] := subs(x[u] = Y.u[n],soln[i]) od od; RETURN(eval(soln)) end: `help/text/Euler` := TEXT( `FUNCTION: stochastic[Euler] - constructs stochastic Taylor schemes of`, `strong order 0.5, otherwise known as Euler schemes.`, ` `, `CALLING SEQUENCE: Euler([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]);`, ` `, `PARAMETERS: a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call Euler([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]); returns the`, ` Euler scheme for an N-dimensional Stochastic Differential Equation (SDE)`, ` with M-dimensional noise which has as its drift coefficients a1,..,aN and`, ` [b11,..,b1M],..,[bN1,..,bNM] as its diffusion coefficients.`, ` The output consists of the variables YN[n], DeltaWM[n], and Delta[n].`, ` YN[n] denotes the strong order 0.5 stochastic Taylor approximation`, ` to x[N] at the n-th step. DeltaWM[n] denotes the change in the`, ` M-dimensional Wiener process at the n-th step. Delta[n] denotes the step`, ` size at the n-th step.`, ` `, ` This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[Euler].`, ` `, `EXAMPLES: `, ` `, `>Euler([x[2],x[1]],[[r,s],[s,r]]);`, ` `, `table([`, ` 1 = (Y1[n + 1] = Y1[n] + Y2[n] Delta[n] + r Delta W1[n] + s Delta W2[n])`, ` 2 = (Y2[n + 1] = Y2[n] + Y1[n] Delta[n] + s Delta W1[n] + r Delta W2[n])`, ` ])`, ` `, `SEE ALSO: stochastic, stochastic[L0], stochastic[LJ], `, ` stochastic[Milstein], stochastic[Taylor1hlf], stochastic[Taylor2]` ): stochastic[Milstein]:=proc(a:list(algebraic),b:list(list(algebraic))) local u,i,soln; for i to nops(a) do soln[i] := Y.i[n+1] = Y.i[n]+L0(x[i],a,b)*Delta[n]+sum('LJ(x[i],b,j)* Delta*W.j[n]','j' = 1 .. nops(op(1,b)))+ sum('sum('LJ(op(j2,op(i,b)),b,j1)*I[j1,j2]', 'j1' = 1 .. nops(op(1,b)))','j2' = 1 .. nops(op(1,b))); for u to nops(a) do soln[i] := subs(x[u] = Y.u[n],soln[i]) od od; RETURN(eval(soln)) end: `help/text/Milstein` := TEXT( `FUNCTION: stochastic[Milstein] - constructs stochastic Taylor schemes of`, `strong order 1, otherwise known as Milstein schemes.`, ` `, `CALLING SEQUENCE: Milstein([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]);`, ` `, `PARAMETERS: a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call Milstein([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]); returns the`, ` Milstein scheme for an N-dimensional Stochastic Differential Equation`, ` (SDE) with M-dimensional noise which has as its drift coefficients`, ` a1,..,aN and [b11,..,b1M],..,[bN1,..,bNM] as its diffusion coefficients.`, ` The output consists of the variables YN[n], DeltaWM[n], Delta[n],`, ` and I[(j1,j2)]. `, ` YN[n] denotes the strong order 1 stochastic Taylor approximation`, ` to x[N] at the n-th step. DeltaWM[n] denotes the change in the`, ` M-dimensional Wiener process at the n-th step. Delta[n] denotes the step`, ` size at the n-th step. I[(j1,j2)] denotes multiple Ito integrals`, ` (refer to Kloeden, P.E., Platen, E.: Numerical Solution of Stochastic`, ` Differential Equations, (Springer-Verlag, 1992)).`, ` `, ` This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[Milstein].`, ` `, `EXAMPLES: `, ` `, `>Milstein([x[2],x[1]],[[r,s],[s,r]]);`, ` `, `table([`, ` 1 = (Y1[n + 1] = Y1[n] + Y2[n] Delta[n] + r Delta W1[n] + s Delta W2[n])`, ` 2 = (Y2[n + 1] = Y2[n] + Y1[n] Delta[n] + s Delta W1[n] + r Delta W2[n])`, ` ])`, ` `, `SEE ALSO: stochastic, stochastic[L0], stochastic[LJ],`, ` stochastic[Euler], stochastic[Taylor1hlf], stochastic[Taylor2]` ): stochastic[Taylor1hlf]:=proc(a:list(algebraic),b:list(list(algebraic))) local u,i,soln; for i to nops(a) do soln[i] := Y.i[n+1] = Y.i[n]+a[i]*Delta[n]+1/2*L0(a[i],a,b)*Delta[n]^2+ sum('op(j,op(i,b))*Delta*W.j[n]+L0(op(j,op(i,b)),a,b)*I[0,j]+ LJ(a[i],b,j)*I[j,0]','j' = 1 .. nops(op(1,b)))+ sum('sum('LJ(op(j2,op(i,b)),b,j1)*I[j1,j2]', 'j1' = 1 .. nops(op(1,b)))','j2' = 1 .. nops(op(1,b)))+sum( 'sum('sum('LJ(LJ(op(p3,op(i,b)),b,p2),b,p1)*I[p1,p2,p3]', 'p1' = 1 .. nops(op(1,b)))','p2' = 1 .. nops(op(1,b)))', 'p3' = 1 .. nops(op(1,b))); for u to nops(a) do soln[i] := subs(x[u] = Y.u[n],soln[i]) od od; RETURN(eval(soln)) end: `help/text/Taylor1hlf` := TEXT( `FUNCTION: stochastic[Taylor1hlf] - constructs stochastic Taylor schemes of`, `strong order 1.5.`, ` `, `CALLING SEQUENCE: Taylor1hlf([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]);`, ` `, `PARAMETERS: a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call Taylor1hlf([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]); returns`, ` the strong order 1.5 approximation for an N-dimensional Stochastic`, ` Differential Equation (SDE) with M-dimensional noise which has as its`, ` drift coefficients a1,..,aN and [b11,..,b1M],..,[bN1,..,bNM] as its`, ` diffusion coefficients.`, ` The output consists of the variables YN[n], DeltaWM[n], Delta[n],`, ` I[(j1,j2)], and I[(j1,j2,j3)]. `, ` YN[n] denotes the strong order 1.5 stochastic Taylor approximation`, ` to x[N] at the n-th step. DeltaWM[n] denotes the change in the`, ` M-dimensional Wiener process at the n-th step. Delta[n] denotes the step`, ` size at the n-th step. I[(j1,j2)] and I[(j1,j2,j3)] denote multiple Ito`, ` integrals (refer to Kloeden, P.E., Platen, E.: Numerical Solution of`, ` Stochastic Differential Equations, (Springer-Verlag, 1992)).`, ` `, ` This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[Taylor1hlf].`, ` `, `EXAMPLES: `, ` `, `>Taylor1hlf([x[2],x[1]],[[r,s],[s,r]]);`, ` `, `table([`, ` 2`, ` 1 = (Y1[n + 1] = Y1[n] + Y2[n] Delta[n] + 1/2 Y1[n] Delta[n] `, ` + r Delta W1[n] `, ` + s I[1, 0] + s Delta W2[n] + r I[2, 0])`, ` 2`, ` 2 = (Y2[n + 1] = Y2[n] + Y1[n] Delta[n] + 1/2 Y2[n] Delta[n] `, ` + s Delta W1[n] `, ` + r I[1, 0] + r Delta W2[n] + s I[2, 0])`, ` ])`, ` `, `SEE ALSO: stochastic, stochastic[L0], stochastic[LJ],`, ` stochastic[Euler], stochastic[Milstein],`, ` stochastic[Taylor2]` ): stochastic[SL0]:=proc(X:algebraic,a:list(algebraic),b:list(list(algebraic))) local part1,part2; part1 := diff(X,t); part2 := sum('a[k]*diff(X,x[k])','k' = 1 .. nops(a)); part1+part2; end: `help/text/SL0` := TEXT( `FUNCTION: stochastic[SL0] - applies the partial differential operator (1.2),`, ` as described in Kloeden, P.E., Platen, E: Numerical Solution of Stochastic`, ` Differential Equations, pg 339, (Springer-Verlag 1992). This is the`, ` Stratonovich version of the operator L0.`, ` `, `CALLING SEQUENCE: SL0(X,[a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]);`, ` `, `PARAMETERS: X - algebraic, given in the variables x[N] and t.`, ` a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` the variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call SL0(X,[a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]); computes the`, ` application of the Stratonovich version of the operator L0 to X, as`, ` described in Kloeden, P.E., Platen, E: Numerical Solution of Stochastic`, ` Differential Equations, pg 339, (Springer-Verlag 1992).`, ` a1,..,aN denotes the drift coefficients of the N-dimensional SDE.`, ` [b11,..,b1M],..,[bN1,..,bNM] denotes the diffusion matrix of dimension`, ` N*M, where M is the dimension of the Wiener process.`, ` The output variables are consistent with the variables used as input.`, ` `, `- This routine is used by the routine stochastic[Taylor2].`, ` In general, it is not intended for use on its own.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[SL0]`, ` `, `EXAMPLES: `, ` `, `>SL0(x[2],[x[1],x[2]],[[r,0],[0,r]]);`, ` `, ` x[2] `, ` `, `SEE ALSO: stochastic, stochastic[L0], stochastic[LJ],`, ` stochastic[MLJ], stochastic[Taylor2]` ): stochastic[correct]:=proc(a:list(algebraic),b:list(list(algebraic)),i) a[i]- 1/2*sum('LJ(op(j,op(i,b)),b,j)','j' = 1 .. nops(op(1,b))); end: `help/text/correct` := TEXT( `FUNCTION: stochastic[correct] - applies the drift-correction formula for`, `converting Ito Stochastic Differential Equations (SDEs) to Stratonovich`, `SDEs.`, ` `, `CALLING SEQUENCE: correct([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]],i);`, ` `, `PARAMETERS: a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` the variables x[N] and t.`, ` i - integer.`, ` `, `SYNOPSIS: `, `- The call correct([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]],i); converts`, ` the drift coefficient of an Ito SDE to its Stratonovich form.`, ` The drift-correction formula can be found in Kloeden, P.E.,`, ` Platen, E.: Numerical Solution of Stochastic Differential Equations,`, ` pg 339, (Springer-Verlag, 1992).`, ` a1,..,aN denotes the drift coefficients of the N-dimensional SDE.`, ` [b11,..,b1M],..,[bN1,..,bNM] denotes the diffusion matrix of dimension`, ` N*M, where M is the dimension of the Wiener process. i = 1..N denotes`, ` the "current" dimension of the SDE.`, ` The output variables are consistent with the variables used as input.`, ` `, `- This routine is used by the routine stochastic[Taylor2].`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[correct].`, ` `, `EXAMPLES: `, ` `, `>correct([x[1],x[2]],[[r,0],[0,r]],2);`, ` `, ` x[2] `, ` `, `>correct([x[1],x[2]],[[x[1],0],[0,r]],1);`, ` `, ` 1/2 x[1]`, ` `, `SEE ALSO: stochastic, stochastic[Taylor2]` ): stochastic[Taylor2]:= proc(a:list(algebraic),b:list(list(algebraic))) local u,i,soln; for i to nops(a) do soln[i] := Y.i[n+1] = Y.i[n]+correct(a,b,i)*Delta[n]+ 1/2*SL0(correct(a,b,i),a,b)*Delta[n]^2+ sum('op(j,op(i,b))*Delta*W.j[n]+SL0(op(j,op(i,b)),a,b)*J[0,j]+ LJ(correct(a,b,i),b,j)*J[j,0]','j' = 1 .. nops(op(1,b))) +sum('sum('LJ(op(j2,op(i,b)),b,j1)*J[j1,j2]+ SL0(LJ(op(j2,op(i,b)),b,j1),a,b)*J[0,j1,j2]+ LJ(SL0(op(j2,op(i,b)),a,b),b,j1)*J[j1,0,j2]+ LJ(LJ(correct(a,b,i),b,j2),b,j1)*J[j1,j2,0]', 'j1' = 1 .. nops(op(1,b)))', 'j2' = 1 .. nops(op(1,b)))+sum( 'sum('sum('LJ(LJ(op(p3,op(i,b)),b,p2),b,p1)*J[p1,p2,p3]', 'p1' = 1 .. nops(op(1,b)))','p2' = 1 .. nops(op(1,b)))', 'p3' = 1 .. nops(op(1,b)))+sum('sum('sum( 'sum('LJ(LJ(LJ(op(m4,op(i,b)),b,m3),b,m2),b,m1)*J[m1,m2,m3,m4]', 'm1' = 1 .. nops(op(1,b)))','m2' = 1 .. nops(op(1,b))) ','m3' = 1 .. nops(op(1,b)))','m4' = 1 .. nops(op(1,b))); for u to nops(a) do soln[i] := subs(x[u] = Y.u[n],soln[i]) od od; RETURN(eval(soln)) end: `help/text/Taylor2` := TEXT( `FUNCTION: stochastic[Taylor2] - constructs stochastic Taylor schemes of`, `strong order 2.`, ` `, `CALLING SEQUENCE: Taylor2([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]);`, ` `, `PARAMETERS: a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call Taylor2([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]); returns the`, ` strong order 2 stochastic Taylor approximation for an N-dimensional`, ` Stochastic Differential Equation (SDE) with M-dimensional noise which has`, ` as its drift coefficients a1,..,aN and [b11,..,b1M],..,[bN1,..,bNM] as its`, ` diffusion coefficients.`, ` The output consists of the variables YN[n], DeltaWM[n], Delta[n],`, ` J[(j1,j2)], J[(j1,j2,j3)], and J[(j1,j2,j3,j4)].`, ` YN[n] denotes the strong order 2 stochastic Taylor approximation`, ` to x[N] at the n-th step. DeltaWM[n] denotes the change in the`, ` M-dimensional Wiener process at the n-th step. Delta[n] denotes the step`, ` size at the n-th step. J[(j1,j2)], J[(j1,j2,j3)], and J[(j1,j2,j3,j4)]`, ` denote multiple Stratonovich integrals (refer Kloeden, P.E., Platen, E.:`, ` Numerical Solution of Stochastic Differential Equations,`, ` (Springer-Verlag, 1992). `, ` `, ` This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[Taylor2].`, ` `, `EXAMPLES: `, ` `, `>Taylor2([x[2],x[1]],[[r,s],[s,r]]);`, ` `, `table([`, ` 2`, ` 1 = (Y1[n + 1] = Y1[n] + Y2[n] Delta[n] + 1/2 Y1[n] Delta[n] `, ` + r Delta W1[n] `, ` + s J[1, 0] + s Delta W2[n] + r J[2, 0])`, ` 2`, ` 2 = (Y2[n + 1] = Y2[n] + Y1[n] Delta[n] + 1/2 Y2[n] Delta[n] `, ` + s Delta W1[n] `, ` + r J[1, 0] + r Delta W2[n] + s J[2, 0])`, ` ])`, ` `, `SEE ALSO: stochastic, stochastic[LJ], stochastic[SL0],`, ` stochastic[Euler], stochastic[Milstein], stochastic[Taylor1hlf],`, ` stochastic[correct].` ): stochastic[wkeuler]:= proc(a:list(algebraic),b:list(list(algebraic))) local u,i,soln; for i to nops(a) do soln[i] := Y.i[n+1] = Y.i[n]+L0(x[i],a,b)*Delta[n]+ sum('LJ(x[i],b,j)*Delta*Ws.j[n]','j' = 1 .. nops(op(1,b))); for u to nops(a) do soln[i] := subs(x[u] = Y.u[n],soln[i]) od od; RETURN(eval(soln)) end: `help/text/wkeuler` := TEXT( `FUNCTION: stochastic[wkeuler] - constructs stochastic Taylor schemes of`, `weak order 1, otherwise known as weak Euler schemes.`, ` `, `CALLING SEQUENCE: wkeuler([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]);`, ` `, `PARAMETERS: a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call wkeuler([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]); returns the`, ` simplified weak euler scheme for an N-dimensional Stochastic Differential`, ` Equation (SDE) with M-dimensional noise which has as its drift`, ` coefficients a1,..,aN and [b11,..,b1M],..,[bN1,..,bNM] as its diffusion`, ` coefficients.`, ` The simplified scheme can be found in Kloeden, P.E. and Platen, E.:`, ` Numerical Solution of Stochastic Differential Equations, pg. 458,`, ` (Springer-Verlag 1992). The output consists of the variables YN[n],`, ` DeltaWsM[n], and Delta[n]. YN[n] denotes the 1st order simplified`, ` weak approximation to x[N] at the n-th step. DeltaWsM[n] denotes the`, ` change in the M-dimensional noise process at the n-th step. Note here`, ` that WsM[n] does not denote standard Wiener processes, but instead are`, ` independent random variables (refer above text). Delta[n] denotes the`, ` step size at the n-th step.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[wkeuler].`, ` `, `EXAMPLES: `, ` `, `>wkeuler([x[1],x[2]],[[r,0],[0,r]]);`, ` `, `table([`, ` 1 = (Y1[n + 1] = Y1[n] + Y1[n] Delta[n] + r Delta Ws1[n])`, ` 2 = (Y2[n + 1] = Y2[n] + Y2[n] Delta[n] + r Delta Ws2[n])`, ` ])`, ` `, `SEE ALSO: stochastic, stochastic[LO], stochastic[LJ], stochastic[wktay2]` ): stochastic[wktay2]:= proc(a:list(algebraic),b:list(list(algebraic))) local u,i,soln; for i to nops(a) do soln[i] := Y.i[n+1] = Y.i[n]+a[i]*Delta[n]+1/2*L0(a[i],a,b)*Delta[n]^2+ sum('(op(j,op(i,b))+1/2*Delta[n]*(L0(op(j,op(i,b)),a,b)+ LJ(a[i],b,j)))*Delta*Ws.j[n]','j' = 1 .. nops(op(1,b)))+1/2* sum('sum('LJ(op(j2,op(i,b)),b,j1)*(Delta^2*Ws.j1[n]*Ws.j2[n]+ V[j1,j2])','j1' = 1 .. nops(op(1,b)))', 'j2' = 1 .. nops(op(1,b))); for u to nops(a) do soln[i] := subs(x[u] = Y.u[n],soln[i]) od od; RETURN(eval(soln)) end: `help/text/wktay2` := TEXT( `FUNCTION: stochastic[wktay2] - constructs stochastic Taylor schemes of`, `weak order 2.`, ` `, `CALLING SEQUENCE: wktay2([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]);`, ` `, `PARAMETERS: a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call wktay2([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]); returns the`, ` simplified weak order 2 stochastic Taylor scheme for an N-dimensional`, ` Stochastic Differential Equation (SDE) with M-dimensional noise which has`, ` as its drift coefficients a1,..,aN and [b11,..,b1M],..,[bN1,..,bNM] as its`, ` diffusion coefficients. The simplified scheme can be found in`, ` Kloeden, P.E. and Platen, E.: Numerical Solution of Stochastic`, ` Differential Equations, pg. 466, (Springer-Verlag 1992).`, ` The output consists of the variables YN[n],`, ` DeltaWsM[n], V[(j1,j2)], and Delta[n]. YN[n] denotes the 2nd order`, ` simplified weak approximation to x[N] at the n-th step. DeltaWsM[n]`, ` denotes the change in the M-dimensional noise process at the n-th step.`, ` Note here that WsM[n] does not denote standard Wiener processes, but`, ` instead are independent random variables (refer above text). V[(j1,j2)]`, ` denotes independent two-point random variables (refer above text).`, ` Delta[n] denotes the step size at the n-th step.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[wktay2].`, ` `, `EXAMPLES: `, ` `, `>wktay2([x[1],x[2]],[[r,0],[0,r]]);`, ` `, `table([`, ` 2`, ` 1 = (Y1[n + 1] = Y1[n] + Y1[n] Delta[n] + 1/2 Y1[n] Delta[n]`, ` + (r + 1/2 Delta[n] r) Delta Ws1[n])`, ` 2`, ` 2 = (Y2[n + 1] = Y2[n] + Y2[n] Delta[n] + 1/2 Y2[n] Delta[n]`, ` + (r + 1/2 Delta[n] r) Delta Ws2[n])`, ` ])`, ` `, `SEE ALSO: stochastic, stochastic[L0], stochastic[LJ], stochastic[wkeuler]` ): stochastic[MLJ]:= proc(X:algebraic,a:list(algebraic),b:list(list(algebraic)),j:integer) local flag; flag := 0; if j = 0 then flag := L0(X,a,b) fi; if flag = 0 then flag := sum('op(j,op(k,b))*diff(X,x[k])', 'k' = 1 .. nops(b)) fi; RETURN(flag) end: `help/text/MLJ` := TEXT( `FUNCTION: stochastic[MLJ] - applies one of the partial differential`, `operators, either L0 or LJ, to X. These operators can be found in Kloeden,`, `P.E., Platen, E.: Numerical Solution of Stochastic Differential Equations,`, `(Springer-Verlag, 1992).`, ` `, `CALLING SEQUENCE: MLJ(X,[a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]],j);`, ` `, `PARAMETERS: X - algebraic, given in the variables x[N] and t.`, ` a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` variables x[N] and t.`, ` j - integer.`, ` `, `SYNOPSIS: `, `- The call MLJ(X,[a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]],j); computes the`, ` application of either LJ or L0 to X. LJ and L0 are partial differential`, ` operators, as described in Kloeden, P.E., Platen, E.: Numerical`, ` Solution of Stochastic Differential Equations, pg 339,`, ` (Springer-Verlag, 1992).`, ` a1,..,aN are the drift coefficients of the N-dimensional Stochastic`, ` Differential Equation (SDE).`, ` [b11,..,b1M],..,[bN1,..,bNM] denotes the diffusion matrix of dimension`, ` N*M, where M is the dimension of the Wiener process. j=1,..,M refers to`, ` the "current" dimension of the Wiener process.`, ` The output variables are consistent with the variables used as input.`, ` `, `- This routine is used by the routine stochastic[wktay3].`, ` In general, it is not intended for use on its own.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[MLJ].`, ` `, `EXAMPLES: `, ` `, `>MLJ(x[2],[x[2],x[1]],[[r,s],[s,r]],2);`, ` `, ` r `, ` `, `>MLJ(x[2],[x[2],x[1]],[[r,s],[s,r]],0);`, ` `, ` x[1] `, ` `, `SEE ALSO: stochastic, stochastic[L0], stochastic[LJ],`, ` stochastic[SL0], stochastic[wktay3].` ): stochastic[wktay3]:=proc(a:list(algebraic),b:list(list(algebraic))) local u,i,soln; for i to nops(a) do soln[i] := Y.i[n+1] = Y.i[n]+a[i]*Delta[n]+ sum('op(j,op(i,b))*Delta*W.j[n]','j' = 1 .. nops(op(1,b)))+ sum('MLJ(a[i],a,b,j0)*I[j0,0]','j0' = 0 .. nops(op(1,b)))+ sum('sum('MLJ(op(j2,op(i,b)),a,b,j1)*I[j1,j2]', 'j2' = 1 .. nops(op(1,b)))','j1' = 0 .. nops(op(1,b)))+ sum('sum('MLJ(MLJ(a[i],a,b,k2),a,b,k1)*I[k1,k2,0]', 'k1' = 0 .. nops(op(1,b)))','k2' = 0 .. nops(op(1,b)))+sum( 'sum('sum('MLJ(MLJ(op(m3,op(i,b)),a,b,m2),a,b,m1)*I[m1,m2,m3]', 'm3' = 1 .. nops(op(1,b)))','m2' = 0 .. nops(op(1,b)))', 'm1' = 0 .. nops(op(1,b))); for u to nops(a) do soln[i] := subs(x[u] = Y.u[n],soln[i]) od od; RETURN(eval(soln)) end: `help/text/wktay3` := TEXT( `FUNCTION: stochastic[wktay3] - constructs stochastic Taylor schemes of`, `weak order 3.`, ` `, `CALLING SEQUENCE: wktay3([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]);`, ` `, `PARAMETERS: a1,..,aN - algebraics, given in the variables x[N] and t.`, ` [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call wktay3([a1,..,aN],[[b11,..,b1M],..,[bN1,..,bNM]]); returns the`, ` weak order 3 stochastic Taylor scheme for an N-dimensional`, ` Stochastic Differential Equation (SDE) with M-dimensional noise which has`, ` as its drift coefficients a1,..,aN and [b11,..,b1M],..,[bN1,..,bNM] as its`, ` diffusion coefficients. The weak order 3 scheme can be found in`, ` Kloeden, P.E., Platen, E.: Numerical Solution of Stochastic Differential`, ` Equations, pg. 468, (Springer-Verlag 1992). The output consists of the`, ` variables YN[n], DeltaWM[n], I[(j1,j2)], I[(j1,j2,j3)], and Delta[n].`, ` YN[n] denotes the 2nd order weak approximation to x[N] at the n-th step.`, ` DeltaWM[n] denotes the change in the M-dimensional Wiener process`, ` at the n-th step. I[(j1,j2)] and I[(j1,j2,j3)] denote multiple Ito`, ` integrals (refer above text). Delta[n] denotes the step size at the n-th`, ` step.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[wktay3].`, ` `, `EXAMPLES: `, ` `, ` `, `>wktay3([x[2],x[1]],[[r,s],[s,r]]);`, ` `, `table([`, ` 1 = (Y1[n + 1] = Y1[n] + Y2[n] Delta[n] + r Delta W1[n] + s Delta W2[n]`, ` + Y1[n] I[0, 0] + s I[1, 0] + r I[2, 0] + Y2[n] I[0, 0, 0]`, ` + r I[1, 0, 0] + s I[2, 0, 0])`, ` 2 = (Y2[n + 1] = Y2[n] + Y1[n] Delta[n] + s Delta W1[n] + r Delta W2[n]`, ` + Y2[n] I[0, 0] + r I[1, 0] + s I[2, 0] + Y1[n] I[0, 0, 0]`, ` + s I[1, 0, 0] + r I[2, 0, 0])`, ` ])`, ` `, `SEE ALSO: stochastic, stochastic[L0], stochastic[MLJ],`, ` stochastic[wkeuler], stochastic[wktay2].` ): stochastic[colour]:= proc(a:list(algebraic),b:list(algebraic)) local temp1,i; for i to nops(a) do temp1[i] := dx[i][t] = a[i]*dt+b[i]*z[t]*dt od; temp1[i] := dz[t] = -gamma*z[t]*dt+beta*dW[t]; RETURN(eval(temp1)) end: `help/text/colour` := TEXT( `FUNCTION: stochastic[colour] - converts Stochastic Differential Equations `, `(SDEs) with white noise into their coloured noise form.`, ` `, `CALLING SEQUENCE: colour([a1,..,aN],[b1,..,bN]);`, ` `, `PARAMETERS: a1,..,aN - algebraics, given in the variables x[N] and t.`, ` b1,..,bN - algebraics, given in the variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call colour([a1,..,aN],[b1,..,bN]); converts N-dimensional SDEs`, ` with scalar white noise into their coloured noise form. Details of this`, ` procedure, together with reasoning and numerical examples can be found`, ` in Milshtein, G.N., Tret'yakov, M.V.: Numerical Solution of Differential`, ` Equations with Coloured Noise, J. Stat. Physics, (to be released).`, ` a1,..,aN denotes the drift coefficients of the N-dimensional SDE.`, ` b1,..,bN denotes the diffusion coefficients of the SDE (scalar noise).`, ` The output consists of the variables z, x[N], W, gamma, beta, and t.`, ` z denotes the Ornstein-Uhlenbeck process, or exponentially correlated`, ` coloured noise. x[N] denotes the state variable of the (N+1)-dimensional`, ` SDE. W denotes a standard Wiener process. gamma and beta denote`, ` parameters, obtainable from experimental data. t denotes time.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[colour].`, ` `, `EXAMPLES: `, ` `, `>colour([a(x[1],t)],[b(x[1],t)]);`, ` `, `table([`, ` 1 = (dx[1][t] = a(x[1], t) dt + b(x[1], t) z[t] dt)`, ` 2 = (dz[t] = - gamma z[t] dt + beta dW[t])`, ` ])`, ` `, `>colour([x[2],x[1]*(alpha-x[1]^2)-x[2]],[0,sigma*x[1]]);`, ` `, `table([`, ` 1 = (dx[1][t] = x[2] dt)`, ` 2`, ` 2 = (dx[2][t] = (x[1] (alpha - x[1] ) - x[2]) dt + sigma x[1] z[t] dt)`, ` 3 = (dz[t] = - gamma z[t] dt + beta dW[t])`, ` ])`, ` `, `SEE ALSO: stochastic` ): stochastic[comm1]:=proc() local LJ1,LJ2,k,j1,j2,flag,p; for p to nargs do if type(args[p],list) <> true then ERROR(`Expecting input to be an expression sequence of lists`) fi od; for k to nargs do for j1 to nops(args[1]) do for j2 to nops(args[1]) do LJ1 := sum('op(j1,args[l])*diff(op(j2,args[k]),x[l])', 'l' = 1 .. nargs); LJ2 := sum('op(j2,args[l])*diff(op(j1,args[k]),x[l])', 'l' = 1 .. nargs); if LJ1 <> LJ2 then flag := 1 fi od od od; if flag = 1 then RETURN(`Commutative noise of the first kind doesn't exist for this system`) else RETURN(`This system exhibits commutative noise of the first kind`) fi; end: `help/text/comm1` := TEXT( `FUNCTION: stochastic[comm1] - informs the user if the diffusion matrix of a`, `Stochastic Differential Equation (SDE) has commutative noise of the first`, `kind.`, ` `, `CALLING SEQUENCE: comm1([b11,..,b1M],..,[bN1,..,bNM]);`, ` `, `PARAMETERS: [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` the variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call comm1([b11,..,b1M],..,[bN1,..,bNM]); returns the user with a`, ` statement indicating whether or not the diffusion matrix of the SDE has`, ` commutative noise of the first kind.`, ` The commutativity condition of the first kind can be found in Kloeden,`, ` P.E., Platen, E.: Numerical Solution of Stochastic Differential`, ` Equations, pg 348, (Springer-Verlag, 1992).`, ` [b11,..,b1M],..,[bN1,..,bNM] denotes the diffusion matrix of dimension`, ` N*M, where N is the dimension of the SDE and M is the dimension of the`, ` Wiener process.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[comm1]`, ` `, `EXAMPLES: `, ` `, `>comm1([1,0],[0,x[1]]);`, ` `, ` Commutative noise of the first kind doesn't exist for this system`, ` `, `SEE ALSO: stochastic, stochastic[comm2]` ): stochastic[comm2]:= proc() local LJ1LJ2,LJ2LJ1,k,p,j1,j2,j3,flag; for p to nargs do if type(args[p],list) <> true then ERROR(`Expecting input to be an expression sequence of lists`) fi; od; for k to nargs do for j1 to nops(args[1]) do for j2 to nops(args[1]) do for j3 to nops(args[1]) do LJ1LJ2 := sum('op(j1,args[m])*diff(sum('op(j2,args[l])* diff(op(j3,args[k]),x[l])','l' = 1 .. nargs),x[m])', 'm' = 1 .. nargs); LJ2LJ1 := sum('op(j2,args[m])*diff(sum('op(j1,args[l])* diff(op(j3,args[k]),x[l])','l' = 1 .. nargs),x[m])', 'm' = 1 .. nargs); if LJ1LJ2 <> LJ2LJ1 then flag := 1 fi; od; od; od; od; if flag = 1 then RETURN(`Commutative noise of the second kind doesn't exist for this system`) else RETURN(`This system exhibits commutative noise of the second kind`) fi; end: `help/text/comm2` := TEXT( `FUNCTION: stochastic[comm2] - informs the user if the diffusion matrix of a`, `Stochastic Differential Equation (SDE) has commutative noise of the second`, `kind.`, ` `, `CALLING SEQUENCE: comm2([b11,..,b1M],..,[bN1,..,bNM]);`, ` `, `PARAMETERS: [b11,..,b1M],..,[bN1,..,bNM] - lists of algebraics, given in`, ` the variables x[N] and t.`, ` `, `SYNOPSIS: `, `- The call comm2([b11,..,b1M],..,[bN1,..,bNM]); returns the user with a`, ` statement indicating whether or not the diffusion matrix of the SDE has`, ` commutative noise of the second kind.`, ` The commutativity condition of the second kind can be found in`, ` Kloeden, P.E., Platen, E.: Numerical Solution of Stochastic Differential`, ` Equations, pg 355, (Springer-Verlag, 1992).`, ` [b11,..,b1M],..,[bN1,..,bNM] denotes the diffusion matrix of dimension`, ` N*M, where N is the dimension of the SDE and M is the dimension of the`, ` Wiener process.`, ` `, `- This routine is part of the stochastic package and is usually loaded via`, ` the call with(stochastic). It can also be invoked via the call`, ` stochastic[comm2].`, ` `, `EXAMPLES: `, ` `, `>comm2([1,(x[2])^2*(x[1])^4],[0,(x[1])^2]);`, ` `, ` Commutative noise of the second kind doesn't exist for this system`, ` `, `SEE ALSO: stochastic, stochastic[comm1]` ):