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$
  • $\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 Jan 3 at 17:50
3
$\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$
  • $\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 Jan 3 at 16:30
  • $\begingroup$ @xzczd I do not know what OP wants, but he incorrectly uses parallel computing. $\endgroup$ – Alex Trounev Jan 3 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 Jan 3 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 Jan 3 at 17:31
  • $\begingroup$ @Edmund Are you sure there will be a time gain? $\endgroup$ – Alex Trounev Jan 3 at 18:22

Your Answer

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

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