1
$\begingroup$

I am writing a code which solve a list of equations and then create a Table with resulst of equation solving. In applications where I use this code, I create a bunch of equations and my computations can be very long (more than 3-4 days).

I investigate that Table creation with ReplaceAll function is my narrow point, which calculate much more longer than other parts. Here is an example:

  In[171]:= ClearAll["Global`*"]

In[172]:= t1 = AbsoluteTime[];

In[173]:= m = 150;
p[1, k_] = 0;
p[m, k_] = 1;
p[j_, 1] = 0;
p[j_, m] = 0;
dx = dy = 0.01;

In[179]:= 
eqns = Flatten[
   Table[1/dx^2 (p[j + 1, k] - 2 p[j, k] + p[j - 1, k]) + 
      1/dy^2 (p[j, k + 1] - 2 p[j, k] + p[j, k - 1]) == 0, {j, 2, 
     m - 1}, {k, 2, m - 1}]];

In[180]:= vars = Flatten[Table[p[j, k], {j, 2, m - 1}, {k, 2, m - 1}]];

In[181]:= sol = Solve[eqns, vars][[1]];

In[182]:= AbsoluteTime[] - t1

Out[182]= 0.8736015

In[183]:= 
pTable = Table[{j dx, k dy, p[j, k] /. sol}, {j, 1, m}, {k, 1, m}];

In[184]:= AbsoluteTime[] - t1

Out[184]= 38.4112677

Is there any way to improve this operation? Last time when I just move "ReplaceAll" command near p[j, k] it speed up calculation around 40%.

Maybe there is another hint, which can speed up replacement in this example?

I would be grateful for any help.

$\endgroup$
3
$\begingroup$

What you have laid out is inefficient in a number of ways. But, I will tackle your question, as is, first, then I will lay out how I would approach it.

Since your bottleneck is

Table[{j dx, k dy, p[j, k] /. sol}, {j, 1, m}, {k, 1, m}]

you need to use a more efficient layout for your replacement rules: Dispatch. Replace

sol = Solve[eqns, vars][[1]];

with

sol = Dispatch@Solve[eqns, vars][[1]];

and you will notice an enormous speed up.

As to how I would approach this problem, I would recast it as a matrix problem. The simplest method for doing that is to use CoefficientArrays which returns something of the form {b, m}. Then, this problem becomes

{b, m} = CoefficientArrays[eqns, vars];
solution = LinearSolve[m, b];

where solution is your p[j,k] array flattened out.

$\endgroup$

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.