# NDSolveIterate with shared variables

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

LaunchKernels[];

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

Do[process[i] = Hold[NDSolveIterate[z, {m, m + del}];
NDSolveProcessSolutions[z, "Forward"];] /.
z -> Symbol["eq" <> ToString[i]], {i, n}];
vars = Table[Symbol["eq" <> ToString[i]], {i, n}];
vars // Evaluate =
Table[NDSolveProcessEquations[{u'[t] == -u[t] +
Table[map[t,j], {j, nAll}].RandomReal[1, nAll],
u == 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 once, I get this

Hold[NDSolveIterate[eq5, {m, m + del}]; NDSolveProcessSolutions[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 // ReleaseHold
eq5

(*output*) NDSolveStateData[<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 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] = NDSolveProcessEquations[{u'[t] == -u[t] +
Table[map[t,j], {j, nAll}].RandomReal[1, nAll],
u == 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];

ParallelDo[
Module[{state0 = eq[i], state1, saved = volt[i], sol},
state1 = NDSolveReinitialize[state0, {u[m] == saved}] // First;
NDSolveIterate[state1, {m, m +del}];
volt[i] =
NDSolveProcessSolutions[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.

• What do you expect to get? – Alex Trounev Dec 30 '18 at 12:59
• 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. – b3m2a1 Dec 30 '18 at 15:49
• @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 – thedude Dec 30 '18 at 21:33

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

With

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 =
Function[{i},
NDSolveProcessEquations[
{
u'[t] == -u[t] + Table[map[t, j], {j, nAll}].RandomReal[1, nAll],
u == 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];
LaunchKernels[];
SetSharedVariable[res];


Finally ParallelDo is called. Note that since NDSolveProcessSolutions 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.

ParallelDo[
Module[{state1, r},
state1 =
NDSolveReinitialize[state[[i]], {u[m] == volt[[i]]}] // First;
NDSolveIterate[state1, {m, m + del}];
r = NDSolveProcessSolutions[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.

• 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). – thedude Dec 31 '18 at 23:45
• @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. – Edmund Jan 3 '19 at 11:07