Calculate the temperature distr

42 ビュー (過去 30 日間)
Ayoub Ahmed
Ayoub Ahmed 2023 年 4 月 14 日
回答済み: Bidyut 2023 年 7 月 8 日
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 件のコメント
Rik
Rik 2023 年 4 月 17 日
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
Rik 2023 年 4 月 18 日
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

サインインしてコメントする。

採用された回答

Les Beckham
Les Beckham 2023 年 4 月 14 日
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

その他の回答 (4 件)

Image Analyst
Image Analyst 2023 年 4 月 14 日
編集済み: Image Analyst 2023 年 4 月 16 日
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.

John D'Errico
John D'Errico 2023 年 4 月 14 日
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 件のコメント
Les Beckham
Les Beckham 2023 年 4 月 14 日
Except, of course it actually was
C Create C matrix here
:-)
Walter Roberson
Walter Roberson 2023 年 4 月 16 日
Could have been
! Create C matrix here
if freeform input was being permitted.

サインインしてコメントする。


VBBV
VBBV 2023 年 4 月 16 日
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

Bidyut
Bidyut 2023 年 7 月 8 日
% 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');

カテゴリ

Help Center および File ExchangeMathematics についてさらに検索

Community Treasure Hunt

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

Start Hunting!

Translated by