# Calculate the temperature distr

15 views (last 30 days)

Show older comments

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

##### 5 Comments

Rik
on 17 Apr 2023

Rik
on 18 Apr 2023

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

### Accepted Answer

Les Beckham
on 14 Apr 2023

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)

disp(iteration)

imagesc(T)

colorbar

##### 0 Comments

### More Answers (4)

Image Analyst
on 14 Apr 2023

Edited: Image Analyst
on 16 Apr 2023

The two pieces of advice I give every programmer:

- 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.
- 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.

##### 0 Comments

John D'Errico
on 14 Apr 2023

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.

##### 2 Comments

Les Beckham
on 14 Apr 2023

Except, of course it actually was

C Create C matrix here

:-)

Walter Roberson
on 16 Apr 2023

Could have been

! Create C matrix here

if freeform input was being permitted.

VBBV
on 16 Apr 2023

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

##### 0 Comments

Bidyut
on 8 Jul 2023

% 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');

##### 0 Comments

### See Also

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!