ParallelDo works sometimes only

The following ParallelDo works sometimes (If I restart Mathematica, sometimes works). But, it usually does not work. Had not it worked all the time, I can be sure there is a coding error, but as it works sometimes, I am not sure what's going on. And If I replace the ParallelDo with Do, it works well always. What would be the reason?

c = 1000;
u = 2;
Array[T, u];
F0[x_] := x^2;
ParallelDo[n = 1;
If[j == 1,
T[j] = AbsoluteTiming[
Do[h[j, n] =
NDSolveValue[{y''[x] == Cos[y[x]], y'[0] == 2.606,
y[0] == (n - 1)/c*2*Pi - Pi}, y, {x, 0, 1},
"Method" -> {"EquationSimplification" -> {Automatic,
"TimeConstraint" -> Infinity}}, AccuracyGoal -> 3], {n,
c}
];];,
T[j] = AbsoluteTiming[
Do[h[j, n] =
NDSolveValue[{y''[x] == Cos[y[x]] + F0[h[j - 1, 1][x]],
y'[0] == 2.606, y[0] == (n - 1)/c*2*Pi - Pi}, y, {x, 0, 1},
"Method" -> {"EquationSimplification" -> {Automatic,
"TimeConstraint" -> Infinity}}, AccuracyGoal -> 3.],
{n, c}];];];, {j, u}];

NDSolveValue::underdet: There are more dependent variables, {h[1,1][x],y[x]}, than equations, so the system is underdetermined.


And using ParallelDo instead of Do would make the computation significantly faster for this code?

• I have the same issue in MMA 10.2. Parallelize[] set up over Do[] doesn't help either... – Quit007 Oct 23 '17 at 22:37

First and foremost I believe, it is the = command with having the index on both sides of the equation which prevents ParallelDo[..] or Parallelize[Do[..]] from MMa comprehending the code, because in my case I have a totally different code, with the only similarity being the =-command and having the index on both sides. (This is just an idea without any further theory..)

Second, what you are trying to solve is an equation where h[j-1,...] is necessary to know before being able to solve h[j,...]. With this, it is a kind of recursive formula, which is not suited for parallel computation, because (at my first glance) it must be solved sequentially.

However, obviously from all the values of h you compute, you need only the first h[j-1,1] to start the computation. Computing h[j,n] with n going from 1 to c can, on the other hand, be parallelized.

Therefore I would suggest the following modifications for your code:

c = 1000;
cVals = Range[c];
u = 3;
F0[x_] := x^2;

(*this calculates the very first values*)
h1 = ParallelMap[
NDSolveValue[{y''[x] == Cos[y[x]], y'[0] == 2.606,
y[0] == (# - 1)/c*2*Pi - Pi}, y, {x, 0, 1},
"Method" -> {"EquationSimplification" -> {Automatic,
"TimeConstraint" -> Infinity}}, AccuracyGoal -> 3] &, cVals];

(*Now we need the recursion - (Do[..] would of course do the job too, but I don't see how to parallelize it further)*)
j = 2;
allHsolutions = {h1};
While[j <= u,
hNew = ParallelMap[
NDSolveValue[{y''[x] == Cos[y[x]] + F0[allHsolutions[[j - 1, 1]][x]], y'[0] == 2.606,
y[0] == (# - 1)/c*2*Pi - Pi}, y, {x, 0,"Method" -> {"EquationSimplification" -> {Automatic, "TimeConstraint" -> Infinity}}, AccuracyGoal -> 3] &, cVals];
AppendTo[allHsolutions, hNew];
j++;
];
allHsolutions


Here I used allHsolutions as an array which will get filled with all the calculated solutions. There might be still faster and more elegant methods, but this might be closer to the form you wanted to use.

As you can see the parallel part is now done without assigning in-line left-hand-side to right-hand-side equations, while storing the calculated functions in an array and accessing them by their index n and j as in allHsolutions[[j,n]]

• Thank you Quit007 – user16308 Nov 27 '17 at 0:24