I'm trying to parallelize a set of NDSolves that share a few variables. Here is a short working example:


n = 20; nAll = 27; m = 0; del = 1/10; theta = 2;
s = RandomReal[10, {100, nAll}]; SetSharedVariable[s];

map0 = Compile[{{y, _Real, 1}, {t, _Real}},
     cur = If[t - .1 < # < t, Exp[(# - t)/0.1], 0.] & /@ y; 
    CompilationTarget -> "C", RuntimeOptions -> "Speed"];
map[t_?NumericQ, j_] := map0[s[[j]] + del, t];

Do[process[i] = Hold[NDSolve`Iterate[z, {m, m + del}];
     NDSolve`ProcessSolutions[z, "Forward"];] /. 
    z -> Symbol["eq" <> ToString[i]], {i, n}];
vars = Table[Symbol["eq" <> ToString[i]], {i, n}];
vars // Evaluate = 
  Table[NDSolve`ProcessEquations[{u'[t] == -u[t] + 
        Table[map[t,j], {j, nAll}].RandomReal[1, nAll], 
      u[0] == 0, 
      WhenEvent[u[t] > theta, {u[t] -> 0, AppendTo[s[[i]], t]}, 
       "LocationMethod" -> "StepEnd"]}, u, t, 
     Method -> {"StiffnessSwitching", 
       "EquationSimplification" -> "Solve"}, MaxSteps -> 10^7] // 
    First, {i, n}];

Now, if I evaluate process[5] once, I get this

Hold[NDSolve`Iterate[eq5, {m, m + del}]; NDSolve`ProcessSolutions[eq5, "Forward"];]

which should iterate NDSolve once I apply ReleaseHold since eq5 will point to NDSolve'StateData["<" 0. ">"]. This works fine in the main kernel:

process[5] // ReleaseHold

(*output*) NDSolve`StateData[<0.,0.1>]

I've been trying to parallelize this using something like:

SetSharedVariable[#] & /@ vars    
ParallelDo[process[i] // ReleaseHold, {i,n}]

but I couldn't find a way to keep the symbol and its values separate. I tried to change eq1 with eq[1] but then NDSolve'Iterate doesn't accept it as an input. Any help would be appreciated. Thank you!

edit: The reason I want to share these variables is that I'll iterate this over different values of m and need the variables to be passed between kernels.

edit2: I have tried using NDSolve'Reinitialize to avoid setting the equations as shared variables. With same definitions:

Do[volt[i] = 0, {i, n}]; SetSharedFunction[volt];
temp = {}; Do[
 eq[i] = NDSolve`ProcessEquations[{u'[t] == -u[t] + 
        Table[map[t,j], {j, nAll}].RandomReal[1, nAll], 
      u[0] == 0, 
      WhenEvent[u[t] > theta, {u[t] -> 0, AppendTo[s[[i]], t]}, 
       "LocationMethod" -> "StepEnd"]}, u, t, 
     Method -> {"StiffnessSwitching", 
       "EquationSimplification" -> "Solve"}] // First, {i, n}];
DistributeDefinitions[eq, temp, map];

 Module[{state0 = eq[i], state1, saved = volt[i], sol}, 
  state1 = NDSolve`Reinitialize[state0, {u[m] == saved}] // First; 
  NDSolve`Iterate[state1, {m, m +del}]; 
  volt[i] = 
   NDSolve`ProcessSolutions[state1, "Forward"] // Values // First;
  If[temp != {}, s[[i]] = Join[s[[i]], temp]; temp = {}]], {i, n}]

which still doesn't work. NDSolve doesn't seem to accept s as a shared variable.

  • $\begingroup$ What do you expect to get? $\endgroup$ Dec 30 '18 at 12:59
  • 1
    $\begingroup$ Your issue is in your call to SetSharedVariable. You can no longer make use of the list vars which does not point to a symbol any more. Try Replace[ToExpression[Table["m"<>i, ...], StandardForm, Hold], Hold[v_]:>SetDharedVar[v], 1]` or some standard pattern like that. Also make sure to read up on Hold and delayed evaluation or whatever the tutorial is called. $\endgroup$
    – b3m2a1
    Dec 30 '18 at 15:49
  • $\begingroup$ @b3m2a1 Please see my update. I got rid of that expression but still can't get NDSolve to work inside ParallelDo, probably because 's' is shared $\endgroup$
    – thedude
    Dec 30 '18 at 21:33

You may operate on lists so you do not need a symbol for every solution.


n = 20; nAll = 27; m = 0; del = 1/10; theta = 2;

and map and map0 as in OP. Then define volt and state as lists.

volt = Range@n;

state =
      u'[t] == -u[t] + Table[map[t, j], {j, nAll}].RandomReal[1, nAll],
      u[0] == 0,
      WhenEvent[u[t] > theta, {u[t] -> 0, AppendTo[s[[i]], t]}, 
       "LocationMethod" -> "StepEnd"]
     u, t,
     Method -> {"StiffnessSwitching", 
       "EquationSimplification" -> "Solve"}] // First
   ] /@ Range@n;

Next a result list, res, is setup. This is the only variable that has a side-effect from ParallelDo so is the only one that needs SetSharedVariable.

res = ConstantArray[Null, n];

Finally ParallelDo is called. Note that since NDSolve`ProcessSolutions can take the longest of the three NDSolve steps its result is not directly assigned to the shared variable, res. Doing so could lock res in the other kernels for extended periods. Instead the result is assigned to a temporary variable, r, for the potentially extended processing and then quickly assigned to the shared variable.

 Module[{state1, r},
  state1 = 
   NDSolve`Reinitialize[state[[i]], {u[m] == volt[[i]]}] // First;
  NDSolve`Iterate[state1, {m, m + del}];
  r = NDSolve`ProcessSolutions[state1, "Forward"] // Values // First;
  res[[i]] = r;
 {i, n}]

ParallelDo automatically DistributeDefinitions on variables in its body. Therefore, the above solution works best when these variables (state, volt, r, m,del) are not memory intensive.

Hope this helps.

  • $\begingroup$ Thank you! For some reason, WhenEvent is not working inside ParallelDo even with this setup. I want to reset u[t] to 0 when it hits the threshold (theta). $\endgroup$
    – thedude
    Dec 31 '18 at 23:45
  • 1
    $\begingroup$ @thedude I agree that WhenEvent is not firing in parallel but is firing in serial. I am reporting to Wolfram to see if this is a known issue with workaround. $\endgroup$
    – Edmund
    Jan 3 '19 at 11:07

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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