% 2D Steady Conduction - Finite Difference clear; clc;% Grid setup nx = 21; ny = 21; % grid points (including boundaries) Lx = 0.1; Ly = 0.1; dx = Lx/(nx-1); dy = Ly/(ny-1);
% Initialize temperature matrix T = zeros(ny, nx);
% Boundary conditions T(1,:) = 100; % top (y=0) T(end,:) = 0; % bottom (y=Ly) T(:,1) = 50; % left (x=0) T(:,end) = 50; % right (x=Lx)
% Gauss-Seidel iteration max_iter = 5000; tolerance = 1e-6; error = 1; iter = 0;
while error > tolerance && iter < max_iter T_old = T; for i = 2:ny-1 for j = 2:nx-1 T(i,j) = (T(i-1,j) + T(i+1,j) + T(i,j-1) + T(i,j+1)) / 4; end end error = max(max(abs(T - T_old))); iter = iter + 1; end % 2D Steady Conduction - Finite Difference clear;
fprintf('Converged in %d iterations\n', iter);
% Plot results x = linspace(0, Lx, nx); y = linspace(0, Ly, ny); [X, Y] = meshgrid(x, y);
figure; contourf(X, Y, T, 20, 'LineColor', 'none'); colorbar; xlabel('x (m)'); ylabel('y (m)'); title('2D Temperature Distribution (°C)'); colormap(jet); axis equal;
% Temperature along vertical centerline mid_x_idx = ceil(nx/2); figure; plot(T(:,mid_x_idx), y, 'k-', 'LineWidth', 2); ylabel('y (m)'); xlabel('Temperature (°C)'); title('Temperature Profile at Center x = 0.05 m'); grid on;Output: Converged in 1073 iterations
Output:
Converged in 1073 iterations
The specific phrasing of the title provides a history of how the file was distributed:
The old days of “rapidshare added patched” are over – and good riddance. Today, you have: The specific phrasing of the title provides a
Goal: temperature vs time for small Biot number (lumped) and for 1D slab by finite difference.
Key equations:
Example (lumped): Sphere, ρ=7800 kg/m3, c=470 J/kgK, r=0.01 m, h=50 W/m2K, T0=200°C, T_inf=20°C. Compute T at t=10 s.
MATLAB (lumped):
rho=7800; c=470; r=0.01; h=50; T0=200; Tinf=20; t=10;
V=4/3*pi*r^3; A=4*pi*r^2;
T = Tinf + (T0-Tinf)*exp(-h*A/(rho*V*c)*t);
fprintf('T(10s)=%.2f °C\n',T);
Example (1D slab explicit FD): slab thickness L=0.02 m, k=16 W/mK, rho=7800, c=460, initial T0=100°C, boundaries T=20°C, simulate to 50 s.
MATLAB (explicit FD):
L=0.02; nx=51; dx=L/(nx-1);
k=16; rho=7800; c=460; alpha=k/(rho*c);
dt=0.01; nt=5000; % ensure dt <= dx^2/(2*alpha)
x=linspace(0,L,nx);
T = 100*ones(1,nx);
T([1,end])=20;
for n=1:nt
Tn=T;
for i=2:nx-1
T(i)=Tn(i)+alpha*dt/dx^2*(Tn(i+1)-2*Tn(i)+Tn(i-1));
end
end
plot(x,T); xlabel('x'); ylabel('T (°C)');