设有一块厚度为2δ的二维无限大平板初始温度为t_0在初始瞬间将它放置于温度为t_∞的流体中如图1所示。设t_∞t_0流体与板面间的表面传热系数h为常数。设平板的厚度为01m导热系数10Wm· K密度200Kgm3比热容50JKg·K总时长1min使用MATLAB软件仿真分别绘制0s10s3060s的温度随厚度变化曲线。

根据热传导定律，可以得到平板内部的温度分布方程：

∂^2T/∂x^2 = (1/(α^2)) * ∂T/∂t

其中，T为温度，x为厚度方向的坐标，t为时间，α为热扩散系数，α = sqrt(κ/(ρc))，κ为导热系数，ρ为密度，c为比热容。

根据边界条件，可得到平板表面的边界条件：

-k * ∂T/∂x = h * (T - t∞)

其中，k为导热系数，h为表面传热系数，t∞为流体温度。

为了简化计算，将厚度向量x分为N个等间距的点，时间t分为M个等间距的点。

根据差分近似，可以得到离散形式的温度分布方程：

T(m,n+1) = T(m,n) + (α^2) * (T(m-1,n) - 2T(m,n) + T(m+1,n)) * Δt/Δx^2

其中，T(m,n)表示第m个厚度点在第n个时间点的温度，Δt和Δx分别为时间和厚度的步长。

根据边界条件，可以得到离散形式的边界条件：

T(1,n+1) = T(2,n+1) - (2hΔx/k) * (T(2,n+1) - t∞)

为了计算方便，定义一个矩阵T，T(i,j)表示第i个厚度点在第j个时间点的温度。

根据上述离散方程和边界条件，可以使用MATLAB编写如下代码进行仿真计算：

% 参数设置
L = 0.1; % 平板厚度
N = 100; % 厚度分段数
Δx = L/N; % 厚度步长

T0 = 25; % 初始温度
T∞ = 100; % 流体温度
h = 20; % 表面传热系数
k = 1; % 导热系数
ρ = 200; % 密度
c = 50; % 比热容
α = sqrt(k/(ρ*c)); % 热扩散系数

t = 60; % 总时长
M = 6001; % 时间分段数
Δt = t/(M-1); % 时间步长

% 初始化温度矩阵
T = zeros(N+2, M);
T(:,1) = T0;

% 进行时间步进计算
for n = 1:M-1
    % 计算内部温度
    T(2:N+1,n+1) = T(2:N+1,n) + (α^2) * (T(1:N,n) - 2*T(2:N+1,n) + T(3:N+2,n)) * Δt/Δx^2;
    
    % 计算表面温度
    T(1,n+1) = T(2,n+1) - (2*h*Δx/k) * (T(2,n+1) - T∞);
end

% 绘制温度随厚度变化曲线
x = linspace(0, L, N+2);
plot(x, T(:,1), 'b-', 'LineWidth', 2);
hold on;
plot(x, T(:,101), 'r-', 'LineWidth', 2);
plot(x, T(:,301), 'g-', 'LineWidth', 2);
plot(x, T(:,601), 'm-', 'LineWidth', 2);
hold off;
xlabel('厚度(m)');
ylabel('温度(℃)');
legend('0s', '10s', '30s', '60s');

运行上述代码即可得到0s、10s、30s和60s时刻的温度随厚度变化曲线