       Re: Any idea how to parallelize this small code

• To: mathgroup at smc.vnet.net
• Subject: [mg110884] Re: Any idea how to parallelize this small code
• From: Patrick Scheibe <pscheibe at trm.uni-leipzig.de>
• Date: Sat, 10 Jul 2010 04:00:10 -0400 (EDT)

Hi,

just think your calculation of the new P through: The algorithm relies
on the order of the Do since it uses not only oldP for the calculation
of the "new P". So for calculating the value P[[i]] it is required that
the values of P for smaller i are already calculated to evaluate the
term

Sum[A[[i, j]] P[[j]], {j, 1, i - 1}]

on the right side of your iteration scheme.
This is bad for data-parallelism where you rely on the fact that all
values P[[i]] can be calculated separately from each other. Remember
that the order of calculation in a parallel algorith is often not
determined.

What you should do is to go on a popular "scholar" web search an find
"parallel Successive Over Relaxation" publications.

Then come back to Mathematica and read more about DistributeDefinitions
to make data (in your case the A and B and omega) available to all
started parallel kernels.

Hope this helps a bit.

Cheers
Patrick

On Thu, 2010-07-08 at 03:14 -0400, pratip wrote:
> Dear Group,
>
> Here is a simple code to solve linear systems. It is motivated from
> Fullerton University numerics example. Now if someone can suggest how
> to parallelize this code.
>
> SORmethod[A0_, B0_, P0_, omega_, max_] :=
>  Module[{A = N[A0], B = N[B0], i, j, k = 0, n = Length[P0], P = P0,
>    oldP = P0},
>   While[ k < max,
>    Do[P[[i]] = (1 - omega) oldP[[i]] +
>       omega/A[[i,
>         i]] (B[[i]] - Sum[A[[i, j]] P[[j]], {j, 1, i - 1}] -
>          Sum[A[[i, j]] oldP[[j]], {j, 1 + i, n}]) , {i, 1, n}];
>    oldP = P;
>    k = k + 1; ];
>   Return[P] ]
>
> To test the solver.
>
> n = 10;
> A = DiagonalMatrix[Table[2., {n}]];
> For[i = 1, i <= n - 1, i++,
>   A[[i, i + 1]] = A[[i + 1, i]] = 1.;];
> B = Table[5. - Abs[i - 5], {i, 1, n}];
> P = Table[1., {i, n}];
>  X3 = SORmethod[A, B, P, 1.5, 200] // Chop
>
> Check the above result with Mathematica LinearSolve
>
> LinearSolve[A, B] // N
>
> Straightforward ParallelDo or ParallelSum is not working. One can see
> for large n (abt 400-4000) we will significantly get help if we can
> parallize the code.
>
> Best regards,
>
> Pratip
>

• Prev by Date: Re: replacement x->y except in Exp[x]
• Next by Date: Re: Simple Question
• Previous by thread: I am posting the problem once more
• Next by thread: why the beep at startup