MathGroup Archive 2010

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

Search the Archive

Re: NDSolve - how to bypass safety chceck?

slawek wrote:
> U=BFytkownik "Daniel Lichtblau" <danl at> napisa=B3 w wiadomo=B6ci grup
> dyskusyjnych:i2ok3i$7ni$1 at
>> You neglected to mention what happened when you tried it. Was it
>> something bad (e.g. program hang or frozen computer)?
>> Daniel Lichtblau
>> Wolfram Research
> I solve sets of DDE equations.  The simplest example is:
> {x1'[t]   ==   x1[t]   -   x1[t - 2 d] * x2[t - d]^2  -  2 * x1[x - d] *
> x2[x - 2 d] * x2[x - d],
>  x2'[t]   ==   x2[t]   -   x2[t - 2 d] * x1[t - d]^2  -  2 * x2[x - d] *
> x1[x - 2 d] * x1[x - d]}
> Above equations are related to integro-differential equations, which are can
> not be solved non-numerically and which are nonlinear plasma physics
> equations. The history function is a main trick to solve DDE, because it
> extended domain of time t. (This "trick" and DDE have been never applied in
> published works on the topic.)
> The ODE-INT "old method" take hours to obtain a solution (Fortran), the DDE
> "new method" gives a solution after about half a second (Mathematica 7, the
> same PC).  Wow!
> There are - and it is ok - some subdomains of the time shift d, where
> solutions are chaotic. I use an "module", to generate plots and animations
> (or use Slider/Dynamic):
> plot[d_] := Module[{...}, NDSolve[...]; ...; Plot[...]]
> plots = Table[plot[d],{d,0.05,5.0}]
> ListAnimate[plots]
> For d = 5.0 Mathematica refuses give a solutions and warnig is generate d. It
> is ok, because the solution is chaotic. Nevertheless I would like have even
> completly noisy "wrong" solution instead a blank plot. Thus the question how
> to weak "quality checking" in NDSolve and to force NDSolve to always produce
> an output, even if the covergence is poor.
> NDSolve return when it finish computations, thus - I suppose - Infinity as
> the number of steps would generate infinite computation time and no results.
> I do not look for more accurate results. I look for an option to switch
> NDSolve to a "dumb mode".
> slawek

It is pretty difficult to give a useful response in the absence of 
working code. Below are modifications of your partial code that seem to 
"work", but do not readily exhibit the chaos or early termination of the 

You will notice that I changed "x" to "t" in the DDE independent 
variables, so that it would not be a partial delayed DE; this is in 
accord with your notation, where derivatives are only with respect to t. 
In order to convince NDSolve to actually do anything, I also had to set 
up initial conditions. And sacrifice a small reptile.

eq = {x1'[t] ==
     x1[t] - x1[t - 2 d]*x2[t - d]^2 -
      2*x1[t - d]*x2[t - 2 d]*x2[t - d],
    x2'[t] ==
     x2[t] - x2[t - 2 d]*x1[t - d]^2 -
      2*x2[t - d]*x1[t - 2 d]*x1[t - d], x1[t /; t <= 0] == -2,
    x2[t /; t <= 0] == 1};

plot[p_] :=
  Module[{sol}, sol = NDSolve[eq /. d -> p, {x1[t], x2[t]}, {t, 0, 1}];
    Plot[{x1[t], x2[t]} /. First[sol], {t, 0, 1}]]


Okay, so what to do when you actually run into an early termination? 
Here are a few possibilities.

(1) Use "EventHandler" as a method, and code a restart.

(2) Jettison the DDE approach, and use recurrences to approximate the 
behavior. This is a way of emulating a blind fixed step approach (that 
is to say, there will be no error assessment, so the noise is allowed to 

(3) Possibly lowering Accuracy?PrecisionGoal might help.

(4) Possibly you can use Min and Max to restrict the sizes of 
derivatives. But I'm not sure that is what matters for the situation you 
describe (again, without your working example, it is quite difficult to 
test anything).

Daniel Lichtblau
Wolfram Research

  • Prev by Date: Re: Which inside Module causes problems with ReplaceAll
  • Next by Date: Re: different results for Standard vs Prefix forms
  • Previous by thread: Re: NDSolve - how to bypass safety chceck?
  • Next by thread: Re: NDSolve - how to bypass safety chceck?