NDSolve very very slow

*To*: mathgroup at smc.vnet.net*Subject*: [mg128496] NDSolve very very slow*From*: popov.ghost at gmail.com*Date*: Thu, 25 Oct 2012 01:42:14 -0400 (EDT)*Delivered-to*: l-mathgroup@mail-archive0.wolfram.com*Delivered-to*: l-mathgroup@wolfram.com*Delivered-to*: mathgroup-newout@smc.vnet.net*Delivered-to*: mathgroup-newsend@smc.vnet.net

I'm using NDSolve to solve the ODE: (*Definitions:*) M=1; \[Lambda][l_] = l (l + 1); rinf = 15000; $MinPrecision = 40; wp = $MinPrecision; ac = $MinPrecision - 8; pg = wp/2; eq[\[Omega]_, l_] := \[CapitalPhi]''[r] + (2 (r - M))/( r (r - 2 M)) \[CapitalPhi]'[ r] + ((\[Omega]^2 r^2)/(r - 2 M)^2 - \[Lambda][l]/( r (r - 2 M))) \[CapitalPhi][r] == 0; (*The solution:*) \[CapitalPhi]out[\[Omega]_, l_] := \[CapitalPhi]out[\[Omega], l] = {\[CapitalPhi], \[CapitalPhi]'} /. Block[{$MaxExtraPrecision = 100}, NDSolve[{eq[\[Omega], l], \[CapitalPhi][rinf] == init, \[CapitalPhi]'[rinf] == dinit}, {\[CapitalPhi], \[CapitalPhi]'}, {r, 29, 39}, WorkingPrecision -> wp, AccuracyGoal -> ac, PrecisionGoal -> pg, MaxSteps -> \[Infinity]]][[1]]; (*Run as:*) \[CapitalPhi]out[2, 1]//AbsoluteTiming for example, using the appropriate init, dinit from below. The initial conditions, init and dinit, depend on the parameters \[Omega] and l, and I compute them in a complicated way, but to save on posting all that code I will give them here for just a sample of the troubling values where NDSolve is super slow (namely the large \[Omega] gets the slower it gets because of the solution being very oscillatory): e.g. \[Omega]=2, l=1 init=-0.00003158476708483920624537825468234683185264 + 0.00005870985383714427760467962946161878449483 I; dinit=-0.0001174332599414352606295233010134701480989 - 0.00006318187181746517318131907026619247606411 I; \[Omega]=15, l=5 init=-0.00006498016608463250053612003781430969895230 + 0.00001490041816106928145235191804618494074335 I; dinit=-0.0002235317451490380993971487099054187817855 - 0.0009748334620044207867560932198118869943621 I; \[Omega]=30, l=7 init=0.00006000392711252709381442292707073690333932 - 0.00002905121643121747797219318220828295775626 I; dinit=0.0008716487129159803540979973990019187444996 +0.001800359797589481474882016847763171211252 I; The time that NDSolve takes to solve this system goes from 5minutes in the \[Omega]=2 case, to 3 or 4hrs in the \[Omega]=30 case; I really need to speed this up drastically. I can't help feel that C++ or another system may do it much quicker, but before I take that leap, I would like to ensure Mathematica can't be speeded up for this calculation. 1) Could 'Compile' somehow be used to speed this up? I know supposedley Mathematica compiles automatically to byte code in NDSolve, and there is the Compile->True option (which I have tried without any speed increase), but what about explicitely compiling the functions being fed into NDSolve (I'm thinking of something like what was done here http://forums.wolfram.com/mathgroup/archive/2006/Jan/msg00373.html, could anyone show me how to adapt that to my code?) 2) Specifying somehow that my functions are only numeric and not symbolic with attributes? 3) Different 'Method' in NDSolve? I should mention that for large \[Omega] I do need the result to be reliable to something like 14dp, and that why I am working with high WorkingPrecision etc. Thanks for any help and suggestions.