NDSolve - how to bypass safety chceck?
in [Mathematica]
|
From: slawek on 29 Jul 2010 06:44 U=BFytkownik "Daniel Lichtblau" <danl(a)wolfram.com> napisa=B3 w wiadomo=B6ci grup dyskusyjnych:i2ok3i$7ni$1(a)smc.vnet.net... > 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". TIA slawek
From: Daniel Lichtblau on 30 Jul 2010 06:56 slawek wrote: > U=BFytkownik "Daniel Lichtblau" <danl(a)wolfram.com> napisa=B3 w wiadomo=B6ci grup > dyskusyjnych:i2ok3i$7ni$1(a)smc.vnet.net... >> 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". > > TIA > 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 solver. 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}]] plot[.5] 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 stay. (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
From: slawek on 30 Jul 2010 06:55 U=BFytkownik "sean" <sean_incali(a)yahoo.com> napisa=B3 w wiadomo=B6ci grup dyskusyjnych:i2rm1k$r5h$1(a)smc.vnet.net... > One way to circumvent is to increase the number of steps. For above > system, something like 500000 will do it. I found it by trial and > error. maybe it will work for your system. It is a bad idea, because the set of equations have the chaotic behaviour . By the way, a more precise non-aproximate equations are currently computing by the Adams-Bashforth method with spline integration. In pseudo-code: y'[[i]](t) == y[[i]]](t) - Integrate[ g(t) Sum[y[[i]](t-ta) y[[j]](t-ta-tb) y[[k]](t-2 ta-tb)], ta, tb] , and y is a vector depended on time. Nice thing, a set od ODE-INT, but - I think - may be quite well be approximated by DDE. slawek
From: slawek on 31 Jul 2010 02:41 U=BFytkownik "Daniel Lichtblau" <danl(a)wolfram.com> napisa=B3 w wiadomo=B6ci grup dyskusyjnych:i2ub5q$jfj$1(a)smc.vnet.net... > (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 > stay. Yes, yes, yes. Nevertheless, you can see, that there is no "an easy option" in NDSolve to switch off the "safety check", like: dirtysolution = NDSolve[..., BypassAllChecksAnGiveAnyAndPossibleWrongAnswer->True] or evenmorestupid = NDSolve[..., UseDumbFixedEulerMethod->True] Obviously, it would be funny to write an homemade DumbNDSolve[...] myself. The idea is to use Check[NDSolve[...], DumbNDSolve[...]]//Quiet instead of Check[NDSolve[...],"Unstable"]//Quiet slawek
First
|
Prev
|
Pages: 1 2 Prev: A Question About Directive [off-topic] Next: Scoping constructs Block, Module, ModuleBlock |