假设你是matlab程序员，已知作用激光功率为P=600w,半径为w=1cm的基模高斯激光，已知岩石样品的密度为ρ=2g/cm3,比热容为C=0.75J/(g.K),热传导系数为K=4.4W/(m.K),假设岩石对光吸收率为η=0.6，岩石为长10cm，宽10cm，高15cm的长方体体,初始温度T0＝300K，边界条件：高度方向两表面为第二类边界条件，根据matlab软件利用有限差分法计算激光以v=0.13cm/s的速度，照射3s后的岩石表面沿移动方向温度场和应力场

这个问题需要使用热传导方程和应力方程，然后利用有限差分法进行数值计算。具体步骤如下：

1. 利用热传导方程计算岩石内部温度分布。热传导方程为：

$$\frac{\partial T}{\partial t} = \frac{K}{\rho C}\nabla^2 T - \frac{P\eta}{\pi w^2 \rho C}e^{-\frac{2r^2}{w^2}}$$

其中，T为温度，t为时间，K为热传导系数，ρ为密度，C为比热容，P为激光功率，η为吸收率，w为激光半径，r为激光照射点到岩石中心点的距离。

2. 利用应力方程计算岩石内部应力分布。应力方程为：

$$\nabla \cdot \sigma = 0$$

其中，σ为应力张量。

3. 利用有限差分法进行数值计算。将热传导方程和应力方程离散化，得到差分方程，然后利用迭代方法求解。

4. 计算3秒后岩石表面沿移动方向的温度和应力分布。

在matlab中实现以上步骤可以参考以下代码：

%% 参数设置
P = 600; % 激光功率，单位：W
w = 0.01; % 激光半径，单位：m
rho = 2000; % 岩石密度，单位：kg/m^3
C = 750; % 岩石比热容，单位：J/(kg.K)
K = 4.4; % 岩石热传导系数，单位：W/(m.K)
eta = 0.6; % 岩石吸收率
Lx = 0.1; % 岩石长度，单位：m
Ly = 0.1; % 岩石宽度，单位：m
Lz = 0.15; % 岩石高度，单位：m
dx = 0.001; % 空间步长，单位：m
dt = 0.001; % 时间步长，单位：s
v = 0.13; % 移动速度，单位：m/s
t_final = 3; % 模拟时间，单位：s

%% 离散化
Nx = round(Lx/dx) + 1; % 空间网格数
Ny = round(Ly/dx) + 1;
Nz = round(Lz/dx) + 1;
x = linspace(0, Lx, Nx); % 网格节点
y = linspace(0, Ly, Ny);
z = linspace(0, Lz, Nz);
T = ones(Nx, Ny, Nz) * 300; % 初始温度场
Sxx = zeros(Nx, Ny, Nz); % 初始应力场
Syy = zeros(Nx, Ny, Nz);
Szz = zeros(Nx, Ny, Nz);
Sxy = zeros(Nx, Ny, Nz);
Syz = zeros(Nx, Ny, Nz);
Sxz = zeros(Nx, Ny, Nz);

%% 计算
for t = dt:dt:t_final
    % 计算温度场
    for i = 2:Nx-1
        for j = 2:Ny-1
            for k = 2:Nz-1
                r = sqrt((x(i)-v*t)^2 + y(j)^2 + z(k)^2);
                T(i,j,k) = T(i,j,k) + dt*K/rho/C * (T(i+1,j,k) - 2*T(i,j,k) + T(i-1,j,k))/dx^2 ...
                    + dt*K/rho/C * (T(i,j+1,k) - 2*T(i,j,k) + T(i,j-1,k))/dx^2 ...
                    + dt*K/rho/C * (T(i,j,k+1) - 2*T(i,j,k) + T(i,j,k-1))/dx^2 ...
                    - dt*P*eta/pi/w^2/rho/C * exp(-2*r^2/w^2);
            end
        end
    end
    
    % 计算应力场
    for i = 2:Nx-1
        for j = 2:Ny-1
            for k = 2:Nz-1
                Sxx(i,j,k) = -K*(T(i+1,j,k)-T(i,j,k))/dx;
                Syy(i,j,k) = -K*(T(i,j+1,k)-T(i,j,k))/dx;
                Szz(i,j,k) = -K*(T(i,j,k+1)-T(i,j,k))/dx;
                Sxy(i,j,k) = 0;
                Syz(i,j,k) = 0;
                Sxz(i,j,k) = 0;
            end
        end
    end
    
    % 处理边界条件
    T(1,:,:) = T(2,:,:);
    T(Nx,:,:) = T(Nx-1,:,:);
    T(:,1,:) = T(:,2,:);
    T(:,Ny,:) = T(:,Ny-1,:);
    T(:,:,1) = T(:,:,2);
    T(:,:,Nz) = T(:,:,Nz-1);
    
    % 输出进度
    fprintf('%.2f s / %.2f s\n', t, t_final);
end

%% 可视化
figure
contourslice(x, y, z, T, [], [], 0:0.01:Lz)
xlabel('x (m)')
ylabel('y (m)')
zlabel('z (m)')
title('温度场')
colorbar

figure
contourslice(x, y, z, Sxx, [], [], 0:0.01:Lz)
xlabel('x (m)')
ylabel('y (m)')
zlabel('z (m)')
title('应力场：Sxx')
colorbar

figure
contourslice(x, y, z, Syy, [], [], 0:0.01:Lz)
xlabel('x (m)')
ylabel('y (m)')
zlabel('z (m)')
title('应力场：Syy')
colorbar

figure
contourslice(x, y, z, Szz, [], [], 0:0.01:Lz)
xlabel('x (m)')
ylabel('y (m)')
zlabel('z (m)')
title('应力场：Szz')
colorbar

运行以上代码，可以得到温度场和应力场的三维可视化结果。其中，温度场表示为等温线，应力场表示为等应力线。可以看出，激光照射后，岩石表面沿移动方向的温度和应力分布都发生了变化。