Speeding up a numerical integrator
- To: mathgroup at christensen.cybernetics.net
- Subject: [mg1891] Speeding up a numerical integrator
- From: Scott.A.Hill at williams.edu (Lancelot)
- Date: Sat, 12 Aug 1995 22:51:53 -0400
- Organization: Williams College, Williamstown MA
Hello. I've found it necessary to write a simple numerical integrator in Mathematica, as NIntegrate is not cooperating with me for some reason. Anyway, here is the "Integrator" Integrator::usage="Integrator[f,{x,xmin,xmax},{y,ymin,ymax},{z,zmin,zmax},div] does a numerical integration over the variables x, y, and z, to the limits specified, with div divisions per dimension."; Integrator[f_,{x_,xmin_,xmax_},{y_,ymin_,ymax_},{z_,zmin_,zmax_},div_]:= ( answer=0.; Do[( part=(f/.{x->i,y->j,z->k})(xmax-xmin)(ymax-ymin)(zmax-zmin)/div^3; answer=N[answer+If[NumberQ[part],part,0]] ), {i,xmin,xmax,N[(xmax-xmin)/div]}, {j,ymin,ymax,N[(ymax-ymin)/div]}, {k,zmin,zmax,N[(zmax-zmin)/div]} ]; answer); I have the "If" statement there in case the function blows up to infinity at some point. I figured it was faster than using a Limit command on everything. However, the thing runs rather slowly. For the command Integrator[x*y*z,{x,0,1},{y,0,1},{z,0,1},20] it took approximately 28 seconds. Replacing the 20 by 40 upped the time to 205 seconds. Note that even at 40 the answer is 0.1346, which is still rather far off from the correct answer of .125. At this rate the thing will take 4000 seconds = over an hour to evaluate Integrator[x*y*z,{x,0,1},{y,0,1},{z,0,1},100], and my function is much more complex than x*y*z. Thanks in advance. / :@-) Scott \