Continuous evaluation of complex calculations

Question

When having a first look at some system, it's often very useful to just play with the parameters and see how it develops. Therefore, I would like to let a certain evaluation run "until I cancel it".

Take a nonlinear differential equation for example, say the Duffing oscillator, $\ddot{x} + \delta \dot{x} + \alpha x + \beta x^3 = \gamma \cos (\omega t)$, and try to investigate when the solutions start to become different when you alter the initial condition of the one evaluation by 1/1000.

I'd like to NDSolve this system, with the ability to look at the current state while doing so.

Note that this was only an example, the answer does not necessarily have to be constrained to NDSolve. Another example would be a cellular automaton that one observes until a certain pattern emerges in some region.

Answer 1

There have been several good examples of nice interactive visualizations in the answers so far.

Your question was a very general one, so I will mention the intended general solution: Monitor.

For example, this combines very nicely with NDSolve's StepMonitor along with a ProgressIndicator:

In[13]:= Monitor[
 NDSolve[{D[u[t, x], {t, 2}] == D[u[t, x], {x, 2}] + Sin[u[t, x]], 
   u[0, x] == E^-x^2, Derivative[1, 0][u][0, x] == 0, 
   u[t, -10] == u[t, 10]}, u, {t, 0, 100}, {x, -10, 10}, 
  StepMonitor :> (sol = u[t, x]; time = t)], 
 Column[{Plot[sol, {x, -10, 10}, PlotRange -> {0, 8}, 
    Filling -> Axis], ProgressIndicator[time, {0, 100}]}]
 ]

It looks like this while it is evaluating:

(This is a modified version of an example in the documentation for StepMonitor.)

Using Monitor[..., ProgressIndicator[...]] is a very useful idiom. I guess I use it several times per week. (With a bit of extra code you can monitor e.g. ParallelTable this way too: https://stackoverflow.com/questions/7352461/monitoring-progress-of-a-parallel-computation-in-mathematica/7376332#7376332.)

Answer 2

I'll try to address the "I want to look at the state while the calculation is running bit". I will demonstrate one way, using Dynamic functionality.

Suppose I am running something like Nest[f,x0,n] which takes a long time and I'd like to look at the state. I can do something like this. I define

f[{x_, y_}] := (Pause[.1]; {Cos[x^2 + y^2], Sin[x^2 - y^2]});

(arbitrarily selected and artificially slowed), and my goal is to obtain

Nest[
 f,
 {.8, .4},
 30
]

but it's taking too long. I can do this:

(*setup*)
lst = {};
(*fancy visualization*)
Dynamic[
 Show[
  Graphics[
     {Arrowheads[.05], Arrow[#]},
     PlotRange -> {{-.3, 1}, {0, 1}},
     Axes -> True
     ] & /@ Partition[lst, 2, 1]]
 ]

(*and now the calculation itself*)
Nest[
 (AppendTo[lst, f[#]]; f[#]) &,
 {.8, .4},
 80
 ]

Basically, I continuously visualise the state of lst and, at each step, add the new value. It looks like this:

Answer 3

The following is a fairly hacky way of achieving continuous evaluation with history.

de[t0_, y0_, yp0_, dt_, a_, b_, c_, d_, \[Omega]_] := 
 Module[{y, t, sol},
  sol = y /. First@NDSolve[{
       y''[t] + d y'[t]  + a y[t] +  b y[t]^3 == c Cos[\[Omega] t],
       y[t0] == y0,
       y'[t0] == yp0},
      y, {t, t0, t0 + dt}];
  sol]

evalpts[f_, t_, dt_, n_] := Table[{x, f[x]}, {x, t, t + dt, dt/n}]

t0 = 0;
y0 = 1;
yp0 = 0;
z0 = 1;
zp0 = 0;

history = 20;
length = 0;
oldest = 0;
ysols = {};
zsols = {};
n = 10;

Manipulate[
 (* Solve differential equation *)
 ysol = de[t0, y0, yp0, dt, a1, b1, c1, d1, w1];
 zsol = de[t0, z0, zp0, dt, a2, b2, c2, d2, w2];

 (* Get points in this segment *)
 ypts = evalpts[ysol, t0, dt, n];
 zpts = evalpts[zsol, t0, dt, n];

 (* Build up our solution history.
 Replace the oldest part of the solution
 with the newest one *)
 If[length < history,
  ysols = Append[ysols, ypts];
  zsols = Append[zsols, zpts];
  length++,
  ysols = ReplacePart[ysols, oldest -> ypts];
  zsols = ReplacePart[zsols, oldest -> zpts];
  ];
 oldest++;
 If[oldest == history + 1, oldest = 1];

 (* Set our initial conditions for next segment *)
 t0 = t0 + dt;
 {y0, yp0} = {ysol[t0], ysol'[t0]};
 {z0, zp0} = {zsol[t0], zsol'[t0]};

 (* Wait to make the animation a bit slower *)
 Pause[0.05];

 (* Display the solution *)
 ListLinePlot[{Sort@Flatten[ysols, 1], Sort@Flatten[zsols, 1]}, 
  PlotRange -> {ymin, ymax}, Axes -> {False, False}],

 {{a1, 0.1}, -1, +1}, {{a2, 0.1}, -1, +1},
 {{b1, 0}, -1, +1}, {{b2, 0.01}, -1, +1},
 {{c1, 0}, -1, +1}, {{c2, 0}, -1, +1},
 {{d1, 0}, -1, +1}, {{d2, 0}, -1, +1},
 {{w1, 0}, 0, 10}, {{w2, 0}, 0, 10},

 {{dt, 0.5}, 0.01, 1, 0.01},
 {{ymin, -1.2}, -10, -1},
 {{ymax, +1.2}, 1, 10}
 ]

You can adjust the parameters in real time. It is fairly fragile though. It basically works by creating a circular queue of previous solutions and continuously replacing the oldest portion of the solution at each step. I added a small delay so it didn't evaluate too fast.

You can then visualize two different solutions with different initial parameters.

There's LOTS of room for improvement for this.

Answer 4

A solution in some cases is to use reassignment of the type x = f[x] inside the Dynamic[…] structures. This example (Game of Life) will run as long as you wish. You can also pause it. This is the power of interactive front end in Mathematica. But generally you can always invent an infinite loop with many means. For numerical computations, you have to be very careful with error accumulation though; after some long time your solution may become completely wrong.

Manipulate[
 CAC = If[run, 
   CellularAutomaton[{224, {2, {{2, 2, 2}, {2, 1, 2}, {2, 2, 2}}}, {1,
       1}}, CAC], CAC];
 ArrayPlot[CAC, PixelConstrained -> 7, Mesh -> All, 
  Frame -> True], {{run, False, "run"}, {False, True}}, 
 Initialization :> (CAC = RandomInteger[1, {50, 50}]), 
 FrameMargins -> 0]