MathGroup Archive 2006

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: ndsolve & compile

  • To: mathgroup at smc.vnet.net
  • Subject: [mg64932] Re: ndsolve & compile
  • From: Christopher Klausmeier <klausme1 at msu.edu>
  • Date: Wed, 8 Mar 2006 00:59:45 -0500 (EST)
  • Sender: owner-wri-mathgroup at wolfram.com

Hi Rob --

Thanks for the help.  The vector form of the equations does speed  
things up nicely on these uncoupled equations.  What I really want to  
do is solve coupled equations (like a spatially-discretized version  
of a reaction-diffusion model), which I haven't had as much luck  
with.  In this case, the vector form you suggested does horribly  
compared to the naive form with n[x]'[t]s in it.  Any suggestions?   
[I know Mathematica can do PDE's directly, but the real problem I  
have in mind contains integrals over the spatial variable, so I  
thought it would be best to discretize in space manually.]

In[2]:=
tmax=100;

r=1;
k=2;
kp=1;
c=1;
m=0.1;

d=0.001;

xmax=256;

dmat =NDSolve`FiniteDifferenceDerivative[2, Range[0,1,1.0/xmax],
         DifferenceOrder\[Rule]2,PeriodicInterpolation\[Rule]True]@
       DifferentiationMatrix[];

eqns=Flatten[Table[{

           n[x]'[t]\[Equal]
             r*n[x][t]*(1-n[x][t]/k)-c*n[x][t]*p[x][t]/(n[x][t]+kp)+
               d*(dmat[[Mod[x-1,xmax,1],x]]*n[Mod[x-1,xmax,1]][t]+
                     dmat[[x,x]]*n[x][t]+
                     dmat[[Mod[x+1,xmax,1],x]]*n[Mod[x+1,xmax,1]][t]),

           p[x]'[t]\[Equal]
             c*n[x][t]*p[x][t]/(n[x][t]+kp)-m*p[x][t]+
               d*(dmat[[Mod[x-1,xmax,1],x]]*p[Mod[x-1,xmax,1]][t]+
                     dmat[[x,x]]*p[x][t]+
                     dmat[[Mod[x+1,xmax,1],x]]*p[Mod[x+1,xmax,1]][t])
           },{x,1,xmax}]];

ics=Flatten[
       Table[{n[x][0]\[Equal]If[0.4*xmax<x<0.6*xmax,1,0],
           p[x][0]\[Equal]0.3},{x,1,xmax}]];

unks=Flatten[Table[{n[x],p[x]},{x,1,xmax}]];

First[Timing[
     NDSolve[Flatten[Join[eqns,ics]],unks,{t,0,tmax},MaxSteps\[Rule] 
Infinity]
     ]]

Out[14]=
13.3953 Second

In[15]:=
First[Timing[
     NDSolve[{n'[t]\[Equal]r*n[t]*(1-n[t]/k)-c*n[t]*p[t]/(n[t]+kp) 
+d*dmat.n[t],
         p'[t]\[Equal]c*n[t]*p[t]/(n[t]+kp)-m*p[t]+d*dmat.p[t],
         n[0]\[Equal]Table[If[0.4*xmax<x<0.6*xmax,1,0],{x,1,xmax}],
         p[0]\[Equal]Table[0.3,{x,1,xmax}]},{n,p},{t,0,tmax},
       MaxSteps\[Rule]Infinity]
     ]]

Out[15]=
65.7573 Second

When diffusion is turned off (d=0), the first method takes 4.5  
seconds but the second takes only 1.2 seconds, so I'd like to harness  
that speed if possible.

thanks again -- Chris
--
Christopher Klausmeier
Kellogg Biological Station
Michigan State University
Hickory Corners MI 49060

Web: http://www.kbs.msu.edu/Faculty/Klausmeier/Index.htm
Email: klausme1 at msu.edu
Phone: (269) 671-4987





  • Prev by Date: Re: Plot derivative
  • Next by Date: Re: ndsolve & compile
  • Previous by thread: Re: ndsolve & compile
  • Next by thread: Re: ndsolve & compile