How can i vectorize this for loop? please help.

t = 0:dt:T;
x = xl:dx:xr;
u = zeros(xn+1,tn+1);
for j = 1:tn;
  tm = j*dt;
  if j == 1; 
    for i = 2:xn;
      T1 = ((gamma(2-bta))/((gamma(j-1+2))*(gamma(1-bta-j+1))))*(u(i+1,1)-2*u(i,1)+u(i-1,1));
      T2 = (2.5*(tm.^1.5)+(gamma(3.5)/gamma(3.5-bta))*(tm.^(2.5-bta)))*sin(x(i));
      u(i,j+1)=u(i,j)+mu.*T1+dt.*(T2);
    end
  end
end

Answer 1

[Jan]

Start with moving all repeated but constant expressions outside the loops:

t = 0:dt:T;
x = xl:dx:xr;
u = zeros(xn+1,tn+1);
c1 = gamma(2-bta);
c4 = gamma(3.5) / gamma(3.5-bta);
for j = 1:tn;
  tm = j*dt;
  c5 = tm.^(2.5-bta);
  c6 = tm.^1.5;
  c7 = 2.5*c6 + c4*c5;
  if j == 1   % ??? Why do you use a loop over j, if only j==1 is considered? 
    c2 = gamma(j-1+2);
    c3 = gamma(1-bta-j+1);
    c8 = (c1 / (c2*c3));
    for i = 2:xn
      T1 = c8 * (u(i+1,1) - 2*u(i,1) + u(i-1,1));
      T2 = c7 * sin(x(i));
      u(i,j+1) = u(i,j) + mu.*T1+dt.*T2;
    end
  end
end

Now the innerloop can be vectorized by moving the index from the for loop insice the calculations: Instead of the block "for i... end"

        T1 = c8 .* (u(3:xn+1,1) - 2 * u(2:xn,1) + u(1:xn-1,1));
        T2 = c7 .* sin(x(2:xn));
        u(2:xn,j+1) = u(2:xn,j) + mu .* T1 + dt .* T2;

I cannot test the code, because some variables are missing.

[Priya Khot]

Thank you so much for the reply. actually i further have "else" condition for that if now m concentrating on that vectorization part.

here, mu = 0.787 and bta = 0.476

t = 0:dt:T;
x = xl:dx:xr;
u = zeros(xn+1,tn+1);
for j = 1:tn;
  tm = j*dt;
  if j == 1; 
    for i = 2:xn;
      T1 = ((gamma(2-bta))/((gamma(j-1+2))*(gamma(1-bta-j+1))))*(u(i+1,1)-2*u(i,1)+u(i-1,1));
      T2 = (2.5*(tm.^1.5)+(gamma(3.5)/gamma(3.5-bta))*(tm.^(2.5-bta)))*sin(x(i));
      u(i,j+1)=u(i,j)+mu.*T1+dt.*(T2);
    end
  else
    for i = 2:xn;
            for k = 1:j;
            T1 = ((gamma(2-bta))/((gamma(j-k+2))*(gamma(1-bta-j+k))))*(u(i+1,k)-2*u(i,k)+u(i-1,k));           
            end
                T2=(2.5*(tm.^1.5)+(gamma(3.5)/gamma(3.5-bta))*(tm.^(2.5-bta)))*sin(x(i));
                u(i,j+1)=u(i,j)+mu.*T1+dt*(T2);
        end
  end
end

now i want to vectorize the "for... loops" after "else". please guide me. Thank you.

[Jan]

This solution suffers again from the repeated calculations of the expensive gamma function for the same constant values. It does not seem that you understand the idea of my suggested code.

The above code overwrites "T1" in the "for k" loop repeatedly. This is not useful, because you can replace the loop by "k=j".