3
$\begingroup$

Say I have a system of differential equations that I need to solve for different initial conditions and say the process is sufficiently expensive and I have sufficiently many that I am calculating in parallel. This system contains a WhenEvent.

As minimal example a model for a bouncing ball with different acceleration due to gravity.

state =
  First@NDSolve`ProcessEquations[
      {
       y''[t] == #, y[0] == 5,
       y'[0] == 0,
       WhenEvent[y[t] == 0, y'[t] -> -0.95 y'[t]]
       },
      y, t,
      Method -> {"StiffnessSwitching", "EquationSimplification" -> "Solve"}
      ] & /@ {-9.81, -5.5, -3.3};

This evaluates in serial without issue with the initial height being varied. Code below is setup for parallel processing.

resS = ConstantArray[Null, 3];
Do[
 Module[{state1, r},
  state1 =
   NDSolve`Reinitialize[state[[i]], {y[0] == i 10}] // First;
  NDSolve`Iterate[state1, {0, 10}];
  r = NDSolve`ProcessSolutions[state1];
  resS[[i]] = r;
  ],
 {i, 3}]
Plot[y[t] /. #, {t, 0, 10}] & /@ resS // GraphicsRow

Mathematica graphics

The affect of WhenEvent is seen with each bounce when y[t] == 0.

However in parallel the WhenEvent does not trigger.

resP = ConstantArray[Null, 3];
LaunchKernels[];
SetSharedVariable[resP];

ParallelDo[
 Module[{state1, r},
  state1 =
   NDSolve`Reinitialize[state[[i]], {y[0] == i 10}] // First;
  NDSolve`Iterate[state1, {0, 10}];
  r = NDSolve`ProcessSolutions[state1];
  resP[[i]] = r;
  ],
 {i, 3}]
Plot[y[t] /. #, {t, 0, 10}] & /@ resP // GraphicsRow

Mathematica graphics

WhenEvent never triggers and the ball falls through the floor forever.

Any ideas on how to trigger WhenEvent for NDSolve in parallel without reinventing the wheel?

$\endgroup$
1
  • $\begingroup$ Seems that Parallel* is having trouble with NDSolve`StateData. It's not even related to WhenEvent. A possible boiling down: state = First@ NDSolve`ProcessEquations[{y''[t] == #, y[0] == 5, y'[0] == 0}, y, t] &[-9.81]; resP = ParallelTable[ NDSolve`Iterate[state, {0, 10}]; NDSolve`ProcessSolutions[state], {1}]. DistributedContexts option doesn't seem to help. $\endgroup$
    – xzczd
    Commented Jan 3, 2019 at 17:50

1 Answer 1

4
$\begingroup$

Your code is complex. It's all easier

g = {-9.81, -5.5, -3.3}; y0 = {10, 20, 30};
ParallelTable[{sol = 
   NDSolve[{y''[t] == g[[i]], y[0] == y0[[i]], y'[0] == 0, 
     WhenEvent[y[t] == 0, y'[t] -> -0.95 y'[t]]}, y, {t, 0, 10}]; 
  Plot[y[t] /. sol, {t, 0, 10}]}, {i, 3}]

fig1

$\endgroup$
5
  • 1
    $\begingroup$ I'm afraid this isn't what OP wants. Those NDSolve`ProcessEquations etc. cannot be circumvented. Notice this question is probably a follow-up of this question. $\endgroup$
    – xzczd
    Commented Jan 3, 2019 at 16:30
  • $\begingroup$ @xzczd I do not know what OP wants, but he incorrectly uses parallel computing. $\endgroup$ Commented Jan 3, 2019 at 16:45
  • $\begingroup$ @xzczd You are correct that the above is not in line with what I am seeking. Keep in mind that the example is just a minimal example to demonstrate the issue. The real problem will have multiple NDSolve`Reinitialize calls per system. $\endgroup$
    – Edmund
    Commented Jan 3, 2019 at 17:26
  • $\begingroup$ @AlexTrounev The actual system is sufficiently expensive such that there is a performance gain when utilizing the NDSolve`* functions instead of NDSolve directly. $\endgroup$
    – Edmund
    Commented Jan 3, 2019 at 17:31
  • $\begingroup$ @Edmund Are you sure there will be a time gain? $\endgroup$ Commented Jan 3, 2019 at 18:22

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.