for a1=1:jx for b1=1:jy JSX=a1*gj-r0;%计算点x坐标 JSY=b1*gj;%计算点y坐标 for ii=1:zks %钻孔数 SLX=ZK(ii,1); SLY=ZK(ii,2); rr=sqrt((SLX-JSX)^2+(SLY-JSY)^2); Xt=JSX-SLX; %计算点到钻孔的x距离 for j=1:nj %计算点z坐标 if (j==1) z=1; elseif(j==nj) z=H-1; else z=(j-1)*dz; end for k=1:nj %计算区间 if(k==1) a=0; b=dz/2; elseif(k==nj) a=H-dz/2; b=H; else a=(2*k-3)*dz*0.5; b=(2*k-1)*dz*0.5; end end end end end end向量矩阵化代码

以下是将该代码矩阵化的过程：

1. 首先将三层循环中的最内层循环提取出来，得到如下代码：

for ii=1:zks %钻孔数
    SLX=ZK(ii,1);
    SLY=ZK(ii,2);
    rr=sqrt((SLX-JSX)^2+(SLY-JSY)^2);
    Xt=JSX-SLX;    %计算点到钻孔的x距离
    for j=1:nj  %计算点z坐标
        if (j==1)
            z=1;
        elseif(j==nj)
            z=H-1;
        else
            z=(j-1)*dz;
        end
        for k=1:nj  %计算区间
            if(k==1)
                a=0;
                b=dz/2;
            elseif(k==nj)
                a=H-dz/2;
                b=H;
            else
                a=(2*k-3)*dz*0.5;
                b=(2*k-1)*dz*0.5;
            end
        end
    end
end

2. 将外层两层循环的变量 a1 和 b1 替换为矩阵，得到如下代码：

[a1, b1] = meshgrid(1:jx, 1:jy);
JSX = a1 .* gj - r0; %计算点x坐标
JSY = b1 .* gj; %计算点y坐标

for ii=1:zks %钻孔数
    SLX=ZK(ii,1);
    SLY=ZK(ii,2);
    rr=sqrt((SLX-JSX).^2+(SLY-JSY).^2);
    Xt=JSX-SLX;    %计算点到钻孔的x距离
    for j=1:nj  %计算点z坐标
        if (j==1)
            z=1;
        elseif(j==nj)
            z=H-1;
        else
            z=(j-1)*dz;
        end
        for k=1:nj  %计算区间
            if(k==1)
                a=0;
                b=dz/2;
            elseif(k==nj)
                a=H-dz/2;
                b=H;
            else
                a=(2*k-3)*dz*0.5;
                b=(2*k-1)*dz*0.5;
            end
        end
    end
end

3. 将内部循环中的变量 j 和 k 替换为矩阵，得到如下代码：

[a1, b1] = meshgrid(1:jx, 1:jy);
JSX = a1 .* gj - r0; %计算点x坐标
JSY = b1 .* gj; %计算点y坐标

[j, k] = meshgrid(1:nj, 1:nj);
z = (j-1) .* dz .* (j~=nj) + (nj-1) .* dz .* (j==nj);
a = (2.*k-3) .* dz .* 0.5 .* (k~=nj) + 0 .* (k==1) + (H-dz/2) .* (k==nj);
b = (2.*k-1) .* dz .* 0.5 .* (k~=1) + dz/2 .* (k==1) + H .* (k==nj);

for ii=1:zks %钻孔数
    SLX=ZK(ii,1);
    SLY=ZK(ii,2);
    rr=sqrt((SLX-JSX).^2+(SLY-JSY).^2);
    Xt=JSX-SLX;    %计算点到钻孔的x距离
end

4. 将内部循环中的变量 SLX, SLY, rr 和 Xt 替换为矩阵，得到最终的矩阵化代码：

[a1, b1] = meshgrid(1:jx, 1:jy);
JSX = a1 .* gj - r0; %计算点x坐标
JSY = b1 .* gj; %计算点y坐标

[j, k] = meshgrid(1:nj, 1:nj);
z = (j-1) .* dz .* (j~=nj) + (nj-1) .* dz .* (j==nj);
a = (2.*k-3) .* dz .* 0.5 .* (k~=nj) + 0 .* (k==1) + (H-dz/2) .* (k==nj);
b = (2.*k-1) .* dz .* 0.5 .* (k~=1) + dz/2 .* (k==1) + H .* (k==nj);

SLX = repmat(ZK(:,1), 1, jx*jy);
SLY = repmat(ZK(:,2), 1, jx*jy);
rr = sqrt((SLX-JSX(:)).^2+(SLY-JSY(:)).^2);
Xt = repmat(JSX(:), zks, 1) - SLX;