Is it possible to make the program faster?

My problem is that i have a program with many operations that runs for more than 12 hours and i would like to ask if there is any possibility to make it run faster. I am using two for loops where the first one runs the timesteps for one year in steps of 2.5s. So there we already have 12.614.400 iterations and within this for loop we have another for loop that runs over the number of layers of a watertank. I would like to know if i could somehow exchange my for loops to make the program run faster.

clear all
roh_w = 983.2;
cp_w = 4183;
lambda_w = 0.6544;
dt = 2.5;
m_strom = 1;
r = 0.5;
h = 2.5;
Tankschichten = 60;
Zeitschritte = 1440;
T_Start = 80;
T_ein = 50;

A_Tank = pi*r*r;
h_Schicht = h/Tankschichten;
v_w = m_strom/(roh_w*A_Tank);
CFL = (v_w*dt)/h_Schicht;
Fo = (lambda_w*dt)/(roh_w*cp_w*h_Schicht*h_Schicht);

T = zeros(Tankschichten,1);

for x = 1 : Zeitschritte
    T_old = T;
    if x == 1
        for y = 1 : Tankschichten
            T(y) = 80;
        end
    else
        for y = 1 : Tankschichten
            if y == 1
                T(y) = (Fo + CFL)*T_ein + (1 - 2*Fo - CFL)*T_old(y) + Fo*T_old(y+1);
            elseif y == Tankschichten
                T(y) = (Fo + CFL)*T_old(y-1) + (1 - 2*Fo - CFL)*T_old(y) + Fo*T_old(y);
            else
                T(y) = (Fo + CFL)*T_old(y-1) + (1 - 2*Fo - CFL)*T_old(y) + Fo*T_old(y+1);
            end
            
        end
        
    end
    T_store(:,x) = T;
end

figure
plot(T_store(50,:))
hold on
ylabel("Temperatur [°C]")
xlabel("Timestep")

Answer 1

[Rik]

Below you will see my attempt and vectorizing your code. As you may notice, this is not actually faster, which shows you that the loops themselves are not really the issue.

The main things missing from your code are documentation and comments. You should also consider pre-allocating T_store with something like T_store = zeros(Tankschichten,Zeitschritte);.

roh_w = 983.2;
cp_w = 4183;
lambda_w = 0.6544;
dt = 2.5;
m_strom = 1;
r = 0.5;
h = 2.5;
Tankschichten = 60;
Zeitschritte = 1440;
T_Start = 80;
T_ein = 50;

A_Tank = pi*r*r;
h_Schicht = h/Tankschichten;
v_w = m_strom/(roh_w*A_Tank);
CFL = (v_w*dt)/h_Schicht;
Fo = (lambda_w*dt)/(roh_w*cp_w*h_Schicht*h_Schicht);

T_store=vectorized(Zeitschritte,Tankschichten,Fo,CFL,T_ein);

plot(T_store(50,:))
hold on
ylabel("Temperatur [°C]")
xlabel("Timestep")

% Run them all as a warm up.
double_loop(Zeitschritte,Tankschichten,Fo,CFL,T_ein);
double_loop2(Zeitschritte,Tankschichten,Fo,CFL,T_ein);
vectorized(Zeitschritte,Tankschichten,Fo,CFL,T_ein);

fprintf('looped: %.5f\n',timeit(@() double_loop(Zeitschritte,Tankschichten,Fo,CFL,T_ein)))

looped: 0.00335

fprintf('looped: %.5f\n',timeit(@() double_loop2(Zeitschritte,Tankschichten,Fo,CFL,T_ein)))

looped: 0.00302

fprintf('vector: %.5f\n',timeit(@() vectorized(Zeitschritte,Tankschichten,Fo,CFL,T_ein)))

vector: 0.00447

isequal(...
    double_loop(Zeitschritte,Tankschichten,Fo,CFL,T_ein) ,...
    vectorized( Zeitschritte,Tankschichten,Fo,CFL,T_ein) )

ans = logical
   1

function T_store=vectorized(Zeitschritte,Tankschichten,Fo,CFL,T_ein)
T = 80*ones(Tankschichten,1);
T_store = 80*ones(Tankschichten,Zeitschritte);
for x = 2 : Zeitschritte
    T_old = T;
    T(end) = [];
    T = (Fo + CFL)*[T_ein;T] + (1 - 2*Fo - CFL)*T_old + Fo*[T_old(2:end);T_old(end)];
    T_store(:,x) = T;
end
end
function T_store=double_loop(Zeitschritte,Tankschichten,Fo,CFL,T_ein)
T = zeros(Tankschichten,1);
for x = 1 : Zeitschritte
    T_old = T;
    if x == 1
        for y = 1 : Tankschichten
            T(y) = 80;
        end
    else
        for y = 1 : Tankschichten
            if y == 1
                T(y) = (Fo + CFL)*T_ein + (1 - 2*Fo - CFL)*T_old(y) + Fo*T_old(y+1);
            elseif y == Tankschichten
                T(y) = (Fo + CFL)*T_old(y-1) + (1 - 2*Fo - CFL)*T_old(y) + Fo*T_old(y);
            else
                T(y) = (Fo + CFL)*T_old(y-1) + (1 - 2*Fo - CFL)*T_old(y) + Fo*T_old(y+1);
            end
            
        end
        
    end
    T_store(:,x) = T;
end
end
function T_store=double_loop2(Zeitschritte,Tankschichten,Fo,CFL,T_ein)
T = zeros(Tankschichten,1);
T_store = zeros(Tankschichten,Zeitschritte);
for x = 1 : Zeitschritte
    T_old = T;
    if x == 1
        for y = 1 : Tankschichten
            T(y) = 80;
        end
    else
        for y = 1 : Tankschichten
            if y == 1
                T(y) = (Fo + CFL)*T_ein + (1 - 2*Fo - CFL)*T_old(y) + Fo*T_old(y+1);
            elseif y == Tankschichten
                T(y) = (Fo + CFL)*T_old(y-1) + (1 - 2*Fo - CFL)*T_old(y) + Fo*T_old(y);
            else
                T(y) = (Fo + CFL)*T_old(y-1) + (1 - 2*Fo - CFL)*T_old(y) + Fo*T_old(y+1);
            end
            
        end
        
    end
    T_store(:,x) = T;
end
end

[Lars Nevermann]

Thank you for your time.