for i=1:m
xx(i,1:n)=0.0025/n*(1:n);
xx(i,n+1:2*n)=0.0025+(0.382-0.0025)/n*(1:n);
xx(i,2*n+1:3*n)=0.382+0.228/5*(1:n);
xx(i,3*n+1:4*n)=0.501+0.228/5*(1:n)
end %确定每个点的x坐标
for j=1:4*n
if xx(1,j)<0.0025
zz(1,j)=sqrt(0.0025^2-(0.0025-xx(1,j))^2);
elseif xx(1,j)<0.382
zz(1,j)=0.0025+(xx(1,j)-0.0025)*tan(4*pi/180);
elseif xx(1,j)<0.610
zz(1,j)=D/2;
else
zz(1,j)=D/2+0.0015;
end
end %确定每一个x处的半径r
for i=2:m
zz(i,:)=zz(1,:)*cos(2*pi/m*(i-1));
yy(i,:)=zz(1,:)*sin(2*pi/m*(i-1));
end %由极坐标求x处圆周上各点的y,z坐标
x=[zeros(m,1) xx];
y=[zeros(m,1) yy];
z=[zeros(m,1) zz];
N=zeros(4*n*m,3);
n1=zeros(4*n*m,3);
m1=zeros(4*n*m,3);
t=zeros(4*n*m,3);
t1=zeros(4*n*m,3);
xb=zeros(4*m*n,1);
yb=zeros(4*m*n,1);
zb=zeros(4*m*n,1);
xyzk=cell(4*m*n,1);
xx(i,1:n)=0.0025/n*(1:n);
xx(i,n+1:2*n)=0.0025+(0.382-0.0025)/n*(1:n);
xx(i,2*n+1:3*n)=0.382+0.228/5*(1:n);
xx(i,3*n+1:4*n)=0.501+0.228/5*(1:n)
end %确定每个点的x坐标
for j=1:4*n
if xx(1,j)<0.0025
zz(1,j)=sqrt(0.0025^2-(0.0025-xx(1,j))^2);
elseif xx(1,j)<0.382
zz(1,j)=0.0025+(xx(1,j)-0.0025)*tan(4*pi/180);
elseif xx(1,j)<0.610
zz(1,j)=D/2;
else
zz(1,j)=D/2+0.0015;
end
end %确定每一个x处的半径r
for i=2:m
zz(i,:)=zz(1,:)*cos(2*pi/m*(i-1));
yy(i,:)=zz(1,:)*sin(2*pi/m*(i-1));
end %由极坐标求x处圆周上各点的y,z坐标
x=[zeros(m,1) xx];
y=[zeros(m,1) yy];
z=[zeros(m,1) zz];
N=zeros(4*n*m,3);
n1=zeros(4*n*m,3);
m1=zeros(4*n*m,3);
t=zeros(4*n*m,3);
t1=zeros(4*n*m,3);
xb=zeros(4*m*n,1);
yb=zeros(4*m*n,1);
zb=zeros(4*m*n,1);
xyzk=cell(4*m*n,1);