clc
clear
clf
%区域几何尺寸及网格划分
H=2;
L=16;
Nx=8;Ny=8;
E=3/4*2*Nx*Ny;%总单元数
N=(Ny+1)*(Nx/2+1)+(Nx/2)*(Ny/2+1);%总节点数
Dx=L/Nx;Dy=H/Ny;%间距
%节点分布拓补
AAA1=zeros(Ny/2,Nx/2+1);AAA2=zeros(Ny/2+1,Nx+1);
AAA1(1,:)=1:Nx/2+1;
for i=2:Ny/2+1
AAA1(i,:)=AAA1(i-1,:)+(Nx/2+1)*ones(1,Nx/2+1);
end
AAA2(1,:)=AAA1(Ny/2,Nx/2+1)+1:AAA1(Ny/2,Nx/2+1)+Nx+1;
for i=2:Ny/2+1
AAA2(i,:)=AAA2(i-1,:)+(Nx+1)*ones(1,Nx+1);
end
%生成三角形单元
for i=1:Ny/2
for j=1:Nx/2+1
JXY(AAA1(i,j),1)=Dx*(j-1);
JXY(AAA1(i,j),2)=Dy*(i-1);
end
end
for i=1:Ny/2+1
for j=1:Nx+1
JXY(AAA2(i,j),1)=Dx*(j-1);
JXY(AAA2(i,j),2)=Dy*(i-1)+Dy*(Ny/2);
end
end
%生成三角形单元JM
k=0;
for i=1:Ny/2
for j=1:Nx/2
k=k+1;
JM(k,:)=[AAA1(i,j),AAA1(i,j+1),AAA1(i+1,j)];
k=k+1;
JM(k,:)=[AAA1(i+1,j+1),AAA1(i+1,j),AAA1(i,j+1)];
end
end
for i=1:Ny/2
for j=1:Nx
k=k+1;
JM(k,:)=[AAA2(i,j),AAA2(i,j+1),AAA2(i+1,j)];
k=k+1;
JM(k,:)=[AAA2(i+1,j+1),AAA2(i+1,j),AAA2(i,j+1)];
end
end
%储存边界节点序号BP1~BP6,用于第一边界条件
BP1=AAA1(1,:)';
BP2=AAA1(:,3);
BP3=AAA2(1,Nx/2+1:Nx+1)';
BP4=AAA2(:,Nx+1);
BP5=AAA2(Ny/2+1,:)';
BP6=[AAA1(:,1);AAA2(2:(Ny/2+1),1)];
%储存第二边界条件的BE,有BE[边界单元序号，边界边序号，边界方向余弦]
%储存第二边界条件的BE,有BE[边界单元序号，边界边序号，边界方向余弦]
tx1=pi/2;ty1=pi;
tx2=0;ty2=pi/2;
tx3=pi/2;ty3=pi;
tx4=0;ty4=pi/2;
tx5=pi/2;ty5=2*pi;
tx6=pi;ty6=pi/2;
%BE1
BE1=[[1:2:Nx]',ones(size([1:Nx/2]')),ones(Nx/2,1)*cos(tx1),ones(Nx/2,1)*cos(ty1)];
%BE2
BE2=[[Nx:Nx:E/3]',3*ones(size([1:Ny/2]')),ones(Ny/2,1)*cos(tx2),ones(Ny/2,1)*cos(ty2)];
%BE3
BE3=[[((Ny/2+1)*Nx+1):2:(E/3+2*Nx-1)]',ones(size([1:Nx/2]')),ones(Nx/2,1)*cos(tx3),ones(Nx/2,1)*cos(ty3)];
%BE4
BE4=[[(Nx*Ny/2+2*Nx):2*Nx:E]',3*ones(size([1:Ny/2]')),ones(Ny/2,1)*cos(tx4),ones(Ny/2,1)*cos(ty4)];
%BE5
BE5=[[(E-2*Nx+2):2:E]',ones(size([1:Nx]')),ones(Nx,1)*cos(tx5),ones(Nx,1)*cos(ty5)];
%BE6
BBB=[[1:Nx:(Nx*(Ny/2-1)+1)],[Nx*Ny/2+1:2*Nx:E-2*Nx+1]]';
BE6=[BBB,3*ones(Ny,1),ones(Ny,1)*cos(tx6),ones(Ny,1)*cos(ty6)];
%调用函数绘制网格
triangle_grid(JM,JXY);
save msh
评论1