how can I make the following for loop faster?

Hello everyone, I'm working in a code and I would like to speed up a for loop. I have this piece of code that I would like to make faster. Maybe someone could give me a hint of something I could do. p and th are (Nx+1,Nz+1) matrices. wp, wth are scalar constants.

    while err > tol

      ptemp = p;
      thtemp = th;

      for i = 2:Nx
          for j = 2:Nz

              Thpnm = thn(i,j);
              hpn   = Hn(i,j);

              if (ptemp(i,j) > 0 || thtemp(i,j) >= 1)

                 PE = ((Thpnm*hpn - h(i,j)*th(i,j))/dt - (U*(h(i,j)*th(i,j) - h(i-1,j)*th(i-1,j)))/(2*dx) + (p(i+1,j)*(htemp3(i,j) ...
                 + htemp3(i+1,j))*(PHx(i,j)/4 + PHx(i+1,j)/4) + p(i-1,j)*(htemp3(i,j) + htemp3(i-1,j))*(PHx(i,j)/4 ...
                 + PHx(i-1,j)/4))/(12*dx^2*nu) + (p(i,j+1)*(htemp3(i,j) + htemp3(i,j+1))*(PHz(i,j)/4 + PHz(i,j+1)/4);

                 p(i,j) = wp*PE + (1 - wp)*ptemp(i,j); 

                      if (p(i,j) >= 0)
                      th(i,j) = 1;
                  else
                      p(i,j) = 0;
                  end
              end

              if (p(i,j) <= 0 || th(i,j) < 1)

                THETA = -(((htemp3(i,j) + htemp3(i+1,j))*(p(i,j) - p(i+1,j))*(PHx(i,j)/4 + PHx(i+1,j)/4) + (htemp3(i,j) + htemp3(i-1,j))*(p(i,j)...
                - p(i-1,j))*(PHx(i,j)/4 + PHx(i-1,j)/4))/(12*dx^2*nu) + ((htemp3(i,j) + htemp3(i,j+1))*(p(i,j) ...
                - p(i,j+1))*(PHz(i,j)/4 + PHz(i,j+1)/4) + (htemp3(i,j) + htemp3(i,j-1))*(p(i,j) - p(i,j-1))*(PHz(i,j)/4;
                th(i,j) = wth*THETA + (1 - wth)*thtemp(i,j);
  %                 th(i,j) = 1;

                  if (th(i,j) < 1)
                      p(i,j) = 0;
                  else
                      th(i,j) = 1;
                  end
              end
          end
      end

      %Compute error
      pmax = max(max(p));
      pmaxtemp = max(max(ptemp));
      chpres = (pmax - pmaxtemp)/pmax;
      thmax = max(max(p));
      thmaxtemp = max(max(ptemp));
      thchpres = (thmax - thmaxtemp)/thmax;    

      err = abs(chpres) + abs(thchpres);
  %   errth = abs(thchpres);
      %Break condition
      if k > 50000
          break
      end

      %Loop increment
      k = k + 1;

end

Answer 1

[Shubham]

Hey Gabriel,

I could give you some hints to speed up your code:

• I noticed that some computations are being repeated, for example:

h(i,j)*th(i,j) % This is being computed multiple times for the variable PE

The repeated computations could be stored and reused.

• Another instance of similar repeated computation is

h(i-1,j)*th(i-1,j) % Notice the indices i-1
% if h(i,j)*th(i,j) is already computed then h(i-1,j)*th(i-1,j) would have
% been computed in a prior iteration

• MATLAB is optimized for matrix and vector operations. Try vectorizing the operations instead of computing within nested for loops.
• Simplify and combine mathematical operations where possible to reduce computational overhead.
• Ensure that all matrices (p, th, etc.) are preallocated before the loop to avoid dynamic resizing during each iteration.
• If you have the Parallel Computing Toolbox, you can use parfor to parallelize the loop.
• Use MATLAB’s profile function to identify bottlenecks in your code.

For further information on improving performance, refer to the following documentation: https://www.mathworks.com/help/matlab/matlab_prog/techniques-for-improving-performance.html#buwj1l5

I hope the above hints were useful!