Calculate the temperature distr

Welcome !

I am attaching my MATLAB code for calculating a temperature

Please advise and appropriate modifications to the solution.

My regards

clear all
close all
clc
%% Input the domain dimension and the required grid.
m = 5; % Number of nodes in x direction in the grid
n = 5; % Number of nodes in y direction in the grid
L = 40; % Length of the plate (cm)
H = 40; % Height of the plate (cm)
dx = L/(m-1); % Step size in x direction (cm)
dy = H/(n-1); % Step size in y direction (cm)
k = 0.49; % coefficient of thermal conductivity (Cal/(s.cm.C))
lamda = 1.5; % Over Relaxation Weighting factor
epsilon = 1e-6;% Pre-specified error as stopping criterion (%)
%% Input temperature distribution around the plate boundary and its 
properties
T_lower = 0; % Lower edge temperature (C)
T_left = 75; % Left edge temperature (C)
T_upper = 100; % Upper edge temperature (C)
T_right = 50; % Right edge temperature (C)
%% Initialize the temoperature distributions
T = zeros(n,m); 
T_old = zeros(n,m);
%% Set up boundary conditions
T(end,2:end-1) = T_lower; 
T(2:end-1,1) = T_left;
T(1,2:end-1) = T_upper; 
T(2:end-1,end) = T_right;
error = 1; % Initial error counter
iteration = 0; % Iteration counter
while error >= epsilon
 for j = m-1:-1:2
 for i = 2:n-1
 T_old(i,j) = T(i,j); 
 T(i,j) = (T(i+1,j)+T(i-1,j)+T(i,j+1)+T(i,j-1))/4;
 T(i,j) = lamda*T(i,j) + (1-lamda)*T_old(i,j);
 error_T(i,j) = 100*abs((T(i,j)-T_old(i,j))/T(i,j));
 end
 end
 error = max(error_T(:)); % Calculate maximum error 
 iteration = iteration + 1; % Calculate number of iterations
end

[Image Analyst]

User edited away the question. Basically he was asking for comments on how to improve his code. That code is pretty much in @VBBV answer below: https://www.mathworks.com/matlabcentral/answers/1947483-is-it-matlab-code-edite#answer_1216573

You might want to reconsider answering any of @mr.Robert's questions in the future.

[Rik]

I recovered the removed content from the Google cache (something which anyone can do). Editing away your question is very rude. Someone spent time reading your question, understanding your issue, figuring out the solution, and writing an answer. Now you repay that kindness by ensuring that the next person with a similar question can't benefit from this answer.

[VBBV]

@Rik, you're right. Wasting others time, by sending private email addresses as well. Definitely so rude and aweful

[Rik]

Why did you edit your question again? Now I have to go back and revert your edit. This is tiresome. The past has also suggested I'm more stubborn. You evidently don't care about other people's time, but you are also wasting your own time. If you want to make an edit, post a comment with what you want to change and why.

[Rik]

For reference, the original question:

Calculate the temperature distribution of the L-shaped plate

Welcome !

I am attaching my MATLAB code to calculate the temperature distribution of the L-shaped plate as in the figure below.

Please give advice and appropriate modifications for the solution.

My regards

-------------------------------------------------------------------------------------------------------------------------------------------------------------------

clear all 
close all 
clc

m =9; 
n =7;
epsilon = 1e-6;
% Initialize the temoperature distributions
T = zeros(n,m);
T_old = zeros(n,m);

T(:,1) = 120:-20:0 

error = 1; % Initial error counter
iteration = 0; % Iteration counter
while error >= epsilon
    for j = m-5:-1:2
        for i = 2:n-1
            T_old(i,j) = T(i,j);
            T(i,j) = (T(i+1,j)+T(i-1,j)+T(i,j+1)+T(i,j-1))/4
            
            error_T(i,j) = 100*abs((T(i,j)-T_old(i,j))/T(i,j));
        end
    end
    
    for j = m-1:-1:5
        for i = 5:n-1
            T_old(i,j) = T(i,j);
            T(i,j) = (T(i+1,j)+T(i-1,j)+T(i,j+1)+T(i,j-1))/4
            
            error_T(i,j) = 100*abs((T(i,j)-T_old(i,j))/T(i,j));
            
        end
    end
    
    
    for j = m-1:-1:2
        T_old(7,j) = T(7,j);
        T(7,j) = (2*T(i-1,j)+T(i,j+1)+T(i,j-1))/4
        
        error_T(7,j) = 100*abs((T(7,j)-T_old(7,j))/T(7,j));
        
    end
    
    
    for j = m-4:-1:2
        
        T_old(1,j) = T(1,j);
        T(1,j) = ((2*T(i+1,j))+T(i,j+1)+T(i,j-1))/4
        
        error_T(1,j) = 100*abs((T(1,j)-T_old(1,j))/T(1,j));
        
    end    
    
    for i = n-3:-1:2;
        T_old(i,5) = T(i,5);
        T(i,5) = (T(i+1,j)+T(i-1,j)+(2*T(i,j-1)))/4
        error_T(i,5) = 100*abs((T(i,5)-T_old(i,5))/T(i,5));
    end
    
    for j = m-1:-1:6;
        T_old(4,j) = T(4,j);
        T(4,j) = ((2*T(i+1,j))+T(i,j+1)+T(i,j-1))/4
        error_T(4,j) = 100*abs((T(4,j)-T_old(4,j))/T(4,j));
    end    
    
    
    error = max(error_T(:)); % Calculate maximum error 
    iteration = iteration + 1; % Calculate number of iterations
end

Answer 1

[Les Beckham]

First of all, format your code (type Ctrl-A then Ctrl-I) in the Matlab editor.

Second, add more comments. For example, is m the number of grid squares in the horizontal direction? And explain what each of the multiple for loops is doing.

It also looks like you have at least one mistake in indexing. It seems like the row 4 column 6 square and the row 7 column 1 square are switched.

clear all 
close all 
clc
m = 9; 
n = 7;
epsilon = 1e-6;
% Initialize the temperature distributions
T = zeros(n,m);
T_old = zeros(n,m);
T(:,1) = 120:-20:0;
error = 1; % Initial error counter
iteration = 0; % Iteration counter
while error >= epsilon
    for j = m-5:-1:2
        for i = 2:n-1
            T_old(i,j) = T(i,j);
            T(i,j) = (T(i+1,j)+T(i-1,j)+T(i,j+1)+T(i,j-1))/4;
            error_T(i,j) = 100*abs((T(i,j)-T_old(i,j))/T(i,j));
        end
    end
    for j = m-1:-1:5
        for i = 5:n-1
            T_old(i,j) = T(i,j);
            T(i,j) = (T(i+1,j)+T(i-1,j)+T(i,j+1)+T(i,j-1))/4;
            error_T(i,j) = 100*abs((T(i,j)-T_old(i,j))/T(i,j));
          end
    end
    for j = m-1:-1:2
        T_old(7,j) = T(7,j);
        T(7,j) = (2*T(i-1,j)+T(i,j+1)+T(i,j-1))/4;
        error_T(7,j) = 100*abs((T(7,j)-T_old(7,j))/T(7,j));
     
    end
    for j = m-4:-1:2
        T_old(1,j) = T(1,j);
        T(1,j) = ((2*T(i+1,j))+T(i,j+1)+T(i,j-1))/4;
        error_T(1,j) = 100*abs((T(1,j)-T_old(1,j))/T(1,j));
 
    end
 
    for i = n-3:-1:2;
        T_old(i,5) = T(i,5);
        T(i,5) = (T(i+1,j)+T(i-1,j)+(2*T(i,j-1)))/4;
        error_T(i,5) = 100*abs((T(i,5)-T_old(i,5))/T(i,5));
    end
    for j = m-1:-1:6;
        T_old(4,j) = T(4,j);
        T(4,j) = ((2*T(i+1,j))+T(i,j+1)+T(i,j-1))/4;;
        error_T(4,j) = 100*abs((T(4,j)-T_old(4,j))/T(4,j));
    end    
    error = max(error_T(:)); % Calculate maximum error 
    iteration = iteration + 1; % Calculate number of iterations
end
disp(error)

   9.2429e-07

disp(iteration)

    93

imagesc(T)
colorbar

Answer 2

[Image Analyst]

The two pieces of advice I give every programmer:

1. Use A LOT of comments. Your code falls way short here. Every for loop, at least, should have comments saying what that for loop does.
2. Use descriptive variable names. No one likes single letter variable names or the program soon falls into an alphabet soup mess of a code. I see i, j, m, n, etc. i and j are usually reserved for the imaginary variable. So use something like row and column for them. m and n are not descriptive at all. You have to find where in the code they were first instantiated. It looks like, from the way you're using it in zeros that n is rows and m is columns. So just name them rows and columns to clarify that. And since the badly-named i iterates over n (rows), i should be called row, and likewise j should be renamed column or col.

Also, you're using a while loop with no failsafe. What happens if your error never goes below 1e-6 (epsilon)? Answer: you get an infinite loop. Look at this more robust while loop to prevent an infinite loop and warn the user if one would have occurred.

% Demonstration of how to avoid an infinite loop by setting up a failsafe.
% Set up a failsafe
maxIterations = 100; % Way more than you think it would ever need.
loopCounter = 0;
% Now loop until we obtain the required condition: a random number equals exactly 0.5.
% If that never happens, the failsafe will kick us out of the loop so we do not get an infinite loop.
r = nan; % Initialize so we can enter the loop the first time.
while (r ~= 0.5) && loopCounter < maxIterations
    loopCounter = loopCounter + 1;
    fprintf('Iteration #%d.\n', loopCounter)
    r = rand;
end
% Alert user if we exited normally, or if the failsafe kicked us out to avoid an infinite loop.
if loopCounter < maxIterations
    % Then the loop found the condition and exited early, which means normally.
    fprintf('Loop exited normally after %d iterations.\n', loopCounter);
else
    % Then the loop never found the condition and exited when the number of iterations
    % hit the maximum number of iterations allowed, which means abnormally.
    fprintf('Loop exited abnormally after iterating the maximimum number of iterations (%d) without obtaining the exit criteria.\n', maxIterations);
end
fprintf('All done after %d iterations.\n', loopCounter)

Next, don't use built-in function names for your variable names. So don't use "error" which is a built-in function name. Call it "maxError" instead.

Answer 3

[John D'Errico]

Little to say about the code. Does it compute the correct result? If it does, then it is code, code that works. Why do you need tips?

Unless code runs too slowly to be useful, as long as it does what it should do, then it is fine as it is. Spend your time writing the next piece of code. Your coding will improve, and surely you would find ways where your code would improve. But don't worry about that. It will happen.

If you really want coding style tips, I would STRONGLY suggest you learn what the comment lines do.

% You know, those things that look like this!

I saw about 2 or 3 lines in the beginning, but then it looks like you got bored with the idea.

Without comments, your code is just a bunch of unreadable lines of code. Useless to anyone in the future. Impossible to debug if you need to do that next month, or to make changes to it next year. If this were code that someone else might need to use, do you expect them to understand what you are doing? Seriously?

My recommendation is always to imagine that one day, you were to get run down by the crosstown bus. Yeah, I know, that sucks. But now imagine that you are the person who inherited this piece of code. They need to be able to read it. To maintain it. Debug it as needed, make changes.

Honestly, if I were the person to inherit this code, I would toss it in the bit bucket, and write it from scratch, IF I needed to do what it does. If I need to spend any significant amount of time to understand inherited code, then I can write it myelf faster. And MY code will be fully internally documented.

I could tell you the story about (long ago) when I inherited a few hundred lines of Fortran code. The guy who wrote it left the company by surprise. But he left some code behind, that was now mine to deal with. In a few hundred lines, there was only ONE comment line in the entire code. It was (and I recall it perfectly):

% Create C matrix here

And then he went on to fill in the elements of a matrix called C. No hints as to what C did, what it meant. I tossed the code, then rewrote it myself.

[Les Beckham]

Except, of course it actually was

C Create C matrix here

:-)

[Walter Roberson]

Could have been

 ! Create C matrix here

if freeform input was being permitted.

Answer 4

[VBBV]

Without knowing more about the governing equations , its difficult to say about the process involved in your problem definition.

clear all 
close all 
clc

m =9; 
n =7;
epsilon = 1e-6;
% Initialize the temoperature distributions
T = zeros(n,m);
T_old = zeros(n,m);

T(:,1) = 120:-20:0; 

error = 1; % Initial error counter
iteration = 0; % Iteration counter
p = 1;
while error(p) >= epsilon
    for j = m-5:-1:2
        for i = 2:n-1
            T_old(i,j) = T(i,j);
            T(i,j) = (T(i+1,j)+T(i-1,j)+T(i,j+1)+T(i,j-1))/4;
            
            error_T(i,j) = 100*abs((T(i,j)-T_old(i,j))/T(i,j));
        end
    end
    
    for j = m-1:-1:5
        for i = 5:n-1
            T_old(i,j) = T(i,j);
            T(i,j) = (T(i+1,j)+T(i-1,j)+T(i,j+1)+T(i,j-1))/4;
            
            error_T(i,j) = 100*abs((T(i,j)-T_old(i,j))/T(i,j));
            
        end
    end
    
    
    for j = m-1:-1:2
        
        T_old(7,j) = T(7,j);
        T(7,j) = (T(7,j+1)+T(7,j-1)+2*T(7,j))/4;
        
        error_T(7,j) = 100*abs((T(7,j)-T_old(7,j))/T(7,j));
        
    end
    
    
    for j = m-4:-1:2
        
        T_old(1,j) = T(1,j);
        T(1,j) = (T(1,j+1)+T(1,j-1)+2*T(1,j))/4;
        
        error_T(1,j) = 100*abs((T(1,j)-T_old(1,j))/T(1,j));
        
    end    
    
    for i = n-3:-1:2;
        T_old(i,5) = T(i,5);
        T(i,5) = (T(i+1,5)+T(i-1,5)+(2*T(i,5)))/4;
        error_T(i,5) = 100*abs((T(i,5)-T_old(i,5))/T(i,5));
    end
    
    for j = m-1:-1:6;
        T_old(4,j) = T(4,j);
        T(4,j) = (T(4,j+1)+T(4,j-1)+2*T(4,j))/4;
        error_T(4,j) = 100*abs((T(4,j)-T_old(4,j))/T(4,j));
    end    
    
    
    error(p+1) = max(error_T(:)); % Calculate maximum error 
    iteration = iteration + 1; % Calculate number of iterations
    p = p+1;
end

[X Y] = meshgrid(linspace(0,1,7),linspace(120,0,length(T)));
v2 = 0:4:130;
contourf(X,Y,T.',v2)
colorbar
shading interp

figure 
imagesc(T)

figure
plot(1:iteration+1,error);title('% error'); xlabel('Iteration'); grid
ax= gca;
ax.YTick = 0:10:110;
ylim([0 110])

However, there seems to be an apparent bug in your program which is present in below sections of code. As @Image Analyst rightly mentioned, there can be potential problem if the while loop goes into infinite loop, when the error never goes below 1e-6, In such cases, you would need a condition to break the loop and stop the program running forever.

% May be this 
for j = m-1:-1:2
    
    T_old(7,j) = T(7,j);
    T(7,j) = (T(7,j+1)+T(7,j-1)+2*T(7,j))/4;
    
    error_T(7,j) = 100*abs((T(7,j)-T_old(7,j))/T(7,j));
    
end

for j = m-4:-1:2    
    T_old(1,j) = T(1,j);
    T(1,j) = (T(1,j+1)+T(1,j-1)+2*T(1,j))/4;
    
    error_T(1,j) = 100*abs((T(1,j)-T_old(1,j))/T(1,j));    
end    

for i = n-3:-1:2;
    T_old(i,5) = T(i,5);
    T(i,5) = (T(i+1,5)+T(i-1,5)+(2*T(i,5)))/4;
    error_T(i,5) = 100*abs((T(i,5)-T_old(i,5))/T(i,5));
end

for j = m-1:-1:6;
    T_old(4,j) = T(4,j);
    T(4,j) = (T(4,j+1)+T(4,j-1)+2*T(4,j))/4;
    error_T(4,j) = 100*abs((T(4,j)-T_old(4,j))/T(4,j));
end    
    

Answer 5

[Bidyut]

% Parameters
L = 30; % Length of the domain
R = 300; % Right boundary
B = 300; % Bottom boundary
T = 300; % Top boundary
Nx = 540; % Number of grid points in x-direction
Ny = 640; % Number of grid points in y-direction
tol = 1e-4; % Tolerance for convergence
omega = 1.5; % Relaxation factor

% Initialize variables
dx = L/(Nx-1);
dy = L/(Ny-1);
x = linspace(0, L, Nx);
y = linspace(0, L, Ny);
T_old = zeros(Ny, Nx);
T_new = zeros(Ny, Nx);
error = 1;

% Boundary conditions
T_new(:, 1) = L; % Left boundary
T_new(:, Nx) = R; % Right boundary
T_new(1, :) = B; % Bottom boundary
T_new(Ny, :) = T; % Top boundary

% SOR method
while error > tol
    for i = 2:Nx-1
        for j = 2:Ny-1
            T_new(j, i) = (1-omega)*T_old(j, i) + omega*(T_new(j, i-1) + T_old(j, i+1) + T_new(j-1, i) + T_old(j+1, i))/4;
        end
    end
    error = max(abs(T_new(:) - T_old(:)));
    T_old = T_new;
end

% Plotting the temperature distribution
[X, Y] = meshgrid(x, y);
figure;
contourf(X, Y, T_new);
colorbar;
title('Temperature Distribution');
xlabel('x');
ylabel('y');L = 30; % Length of the domainR = 300; % Right boundaryB = 300; % Bottom boundaryT = 300; % Top boundaryNx = 540; % Number of grid points in x-directionNy = 640; % Number of grid points in y-directiontol = 1e-4; % Tolerance for convergenceomega = 1.5; % Relaxation factor% Initialize variablesdx = L/(Nx-1);dy = L/(Ny-1);x = linspace(0, L, Nx);y = linspace(0, L, Ny);T_old = zeros(Ny, Nx);T_new = zeros(Ny, Nx);error = 1;% Boundary conditionsT_new(:, 1) = L; % Left boundaryT_new(:, Nx) = R; % Right boundaryT_new(1, :) = B; % Bottom boundaryT_new(Ny, :) = T; % Top boundary% SOR methodwhile error > tol  for i = 2:Nx-1      for j = 2:Ny-1          T_new(j, i) = (1-omega)*T_old(j, i) + omega*(T_new(j, i-1) + T_old(j, i+1) + T_new(j-1, i) + T_old(j+1, i))/4;      end  end  error = max(abs(T_new(:) - T_old(:)));  T_old = T_new;end% Plotting the temperature distribution[X, Y] = meshgrid(x, y);figure;contourf(X, Y, T_new);colorbar;title('Temperature Distribution');xlabel('x');ylabel('y');