MathGroup Archive 2005

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

Search the Archive

Re: NDSolve and differential equation system

  • To: mathgroup at smc.vnet.net
  • Subject: [mg54074] Re: NDSolve and differential equation system
  • From: Maxim <ab_def at prontomail.com>
  • Date: Thu, 10 Feb 2005 02:46:04 -0500 (EST)
  • References: <cu1vjh$ltj$1@smc.vnet.net> <cud7nd$307$1@smc.vnet.net>
  • Sender: owner-wri-mathgroup at wolfram.com

On Wed, 9 Feb 2005 14:46:05 +0000 (UTC), Christian Moosmann <some@Adress>  
wrote:

> Hi again
>
> Thanks for your suggestions, using Erf[] instead of UnitStep[] gains a
> factor of about 2, so this is a step to the right direction.
>
>
>  > You can a) use the wonderful advanced interface of NDSolve[]
>  > to use a constant step size
>
> Sadly the IDA method does not support fixed step sizes. And - correct me
> if I am wrong - only the IDA method can solve differential algebraic
> equations.
>
>  > However, if Ainv.K1 and Ainv.K2 commute
>
> Thanks for your effort here, but sadly they don't :-(
>
>
>
>
> Thanks again very much for your suggestions. Do you have some more?  ;-)
>
> Christian
>

Here's another suggestion: 1) rewrite the equations as a scalar system; 2)  
define v[t] as v[t_?NumericQ]:=...; 3) include the corner points t=1/5 and  
t=9/20 in the iterator of NDSolve; 4) decrease PrecisionGoal; 5) set  
SolveDelayed to False.

In[1]:=
SeedRandom[1]
Dim = 80;

K1 = Table[Random[], {Dim}, {Dim}];
K2 = Table[Random[], {Dim}, {Dim}];
Ainv = Inverse[Table[Random[], {Dim}, {Dim}]];
B = Table[Random[], {Dim}];

v[t_?NumericQ] := UnitStep[t - 2/10]*4*(t - 2/10) -
   UnitStep[t - 45/100]*4*(t - 45/100);

Lsol = Module[
   {Lf = Array[x, {Dim}], x},
   x = Evaluate[Through[Lf[#]]]&;
   NDSolve[Thread /@
     {x'[t] + Ainv.(K1 + v[t]*K2).x[t] == B,
      x[0] == Table[0, {Dim}]},
     Lf, {t, 0, 1/5, 9/20, 1},
     PrecisionGoal -> 6]
]; // Timing

Out[1]=
{75.625 Second, Null}

In[2]:=
MaxMemoryUsed[]

Out[2]=
394147368

In principle Dim=100 is possible too, the memory requirements will be  
below 2Gb. It is interesting that if we compare the performance of NDSolve  
with the definition for v[t] given as v[t_]=... and with the same  
definition where UnitStep is replaced with (Sqrt[#^2]/# + 1)/2&, then  
CPU/memory usage will be much higher in the second case. This suggests  
that the preprocessing of the equations may be the bottleneck here, not  
necessarily related to piecewise functions, hence the idea to use  
v[t_?NumericQ].

A rather strange issue is that repeated evaluations of the same NDSolve  
command take more and more time:

In[3]:=
Table[
   Module[
     {Lf = Array[x, {Dim}], x},
     x = Evaluate[Through[Lf[#]]]&;
     NDSolve[Thread /@
       {x'[t] + Ainv.(K1 + v[t]*K2).x[t] == B,
        x[0] == Table[0, {Dim}]},
       Lf, {t, 0, 1/5, 9/20, 1},
       PrecisionGoal -> 6]
   ]; // Timing,
   {3}][[All, 1]]

Out[3]=
{118.125 Second, 129.968 Second, 145.344 Second}

This may be a problem with internal NDSolve code or something related to  
memory management; no swapping is taking place but we can't be sure about  
garbage collection.

Maxim Rytin
m.r at inbox.ru


  • Prev by Date: Re: Controlling Print[].
  • Next by Date: Re: bugs in Mathematica 5.1
  • Previous by thread: Re: NDSolve and differential equation system
  • Next by thread: Few notes