find the temperature of the wall as a function of position when dealing with transient two-dimensional heat conduction in matlab in matlab

main.m
% Define parameters
L = 1; % Length of plate (m)
W = 0.5; % Width of plate (m)
T1 = 20; % Initial temperature (C)
T2 = 100; % Temperature of left edge (C)
T3 = 0; % Temperature of right edge (C)
k = 0.05; % Thermal conductivity (W/m*C)
rho = 2500; % Density (kg/m3)
cp = 1000; % Specific heat (J/kg*C)
tEnd = 100; % Time to simulate (s)
dx = 0.01; % Spacing in x-direction (m)
dy = 0.01; % Spacing in y-direction (m)
dt = 0.1; % Time step (s)

% Create meshgrid
x = 0:dx:L;
y = 0:dy:W;
[X,Y] = meshgrid(x,y);

% Initialize temperature matrix
T = T1*ones(size(X));

% Set boundary conditions
T(:,1) = T2;
T(:,end) = T3;

% Calculate heat transfer coefficient
h = 10; % Convective heat transfer coefficient (W/m2*C)
alpha = k/(rho*cp);
hx = h*dx/k;
hy = h*dy/k;

% Solve system of equations iteratively
t = 0;
while t < tEnd
    T_old = T;
    for i = 2:size(X,1)-1
        for j = 2:size(X,2)-1
            Tx = (T_old(i+1,j) - 2*T_old(i,j) + T_old(i-1,j))/dx^2;
            Ty = (T_old(i,j+1) - 2*T_old(i,j) + T_old(i,j-1))/dy^2;
            T(i,j) = T_old(i,j) + alpha*dt*(Tx + Ty) - (hx*(T_old(i,j)-T2) + hy*(T_old(i,j)-T3))*dt/(rho*cp);
        end
    end
    t = t + dt;
end

% Plot temperature distribution
surf(X,Y,T)
colorbar
xlabel('x (m)')
ylabel('y (m)')
zlabel('Temperature (C)')
1292 chars
53 lines