【MATLAB】这个代码怎么修改？

#这个代码的问题怎么修改？

对于此运算，数组的大小不兼容。
出错 untitled6 (第 49 行)

clear;
close all;
clc;
k = 0.8; % 热导率常数
c = 960; % 热容量常数
rho = 2400; % 密度常数
h = 6; % 对流传热系数

% 定义域信息
Lx = 0.4; % x方向域长度
Ly = 0.2; % y方向域长度
tmax = 100; % 时间上限

% 定义网格信息
nx = 100; % x方向网格数量
ny = 100; % y方向网格数量
dx = Lx / nx; % x方向网格长度
dy = Ly / ny; % y方向网格长度

% 定义时间步长
dt = 0.1; % 时间步长
nt = ceil(tmax / dt); % 时间步数

% 初始化温度数组
T = 70 * ones(nx, ny); % 设置初始温度为70度

% 设置边界条件
q1 = 500; % 下边界热流密度
q2 = 200; % 上边界热流密度
Q = 20000; % 中心热流密度
tw = 70; % 温度
tf = 25; % 无穷远处温度

% 计算边界温度
T(:,1) = q1 * dx / k + T(:,2); % 左边界温度
T(:,nx) = q2 * dx /k + T(:,nx-1); % 右边界温度
T(ny,:) = Q / (pi * k * Lx * Ly) * ...
    atan(Lx ./ (2 * sqrt((dy * (1:ny)-dy/2).^2 + (Lx/2)^2))) + h * tw / k...
    + T(ny-1,:); % 下边界温度
T(1,:) = tf + (T(2,:) - tf) * exp(-h * dy / (k * 1)); % 上边界温度

% 循环计算温度场
for i = 1:nt
    % 计算二阶导数
    Tx = diff(T,2,1); % 沿x方向
    Ty = diff(T,2,2); % 沿y方向
    
    % 计算二维热扩散方程
    T(2:nx-1,2:ny-1) = T(2:nx-1,2:ny-1) + dt *(k * (Tx + Ty) / (c * rho));
    
    % 设置边界条件
    T(:,1) = q1 * dx / k + T(:,2); % 左边界温度
    T(:,nx) = q2 * dx / k + T(:,nx-1); % 右边界温度
    T(ny,:) = Q / (pi * k * Lx * Ly) * ...
        atan(Lx ./ (2 * sqrt((dy * (1:ny)-dy/2).^2 + (Lx/2)^2))) + h * tw / k ...
        + T(ny-1,:); % 下边界温度
    T(1,:) = tf + (T(2,:) - tf) * exp(-h * dy / (k * 1)); % 上边界温度
    
    % 绘制温度分布图像
    surf(T);
    axis([0 Lx 0 Ly 0 100]);
    drawnow;
end

报错截图发我看看