54 lines
1.4 KiB
Matlab
54 lines
1.4 KiB
Matlab
clc
|
|
clear all
|
|
|
|
model=mphload('SBC+EPD.mph');
|
|
[m1, m2]=mphmeshstats(model);
|
|
|
|
%顶点
|
|
mesh.NbrVertex=max(size(m2.vertex));%顶点数
|
|
mesh.Vertex=m2.vertex'; %顶点坐标 3*n
|
|
|
|
%四面体
|
|
mesh.NbrTet=max(size(m2.elem{2}));%四面体数
|
|
mesh.DomainOfTet=m2.elementity{2};%四面体区域索引 n*1
|
|
mesh.Tet=m2.elem{2}+1;%四面体顶点索引 4*n
|
|
mesh.Tet = sort(mesh.Tet)';
|
|
|
|
%边
|
|
%边排序
|
|
el2no=mesh.Tet';
|
|
n1=el2no([1 1 1 2 2 3],:);
|
|
n2=el2no([2 3 4 3 4 4],:);
|
|
el_ed2no_array=[n1(:) n2(:)];
|
|
%按1-2 1-3 1-4 2-3 2-4 3-4排序
|
|
%Edges 全局边n*2 EdgesOfElements四面体边索引 6N*1 sign四面体边方向 6N*1
|
|
[mesh.Edge,~,mesh.EdgeOfTet]=unique(el_ed2no_array,'rows');
|
|
mesh.NbrEdge=max(size(mesh.Edge));%边数
|
|
mesh.EdgeOfTet=reshape(mesh.EdgeOfTet,6,mesh.NbrTet)';
|
|
|
|
%面
|
|
el2no=mesh.Tet';
|
|
n1=el2no([1 1 1 2],:);
|
|
n2=el2no([2 2 3 3],:);
|
|
n3=el2no([3 4 4 4],:);
|
|
el_ed2no_array=[n1(:) n2(:) n3(:)];
|
|
mesh.Tri=sort(m2.elem{3})'+1;
|
|
mesh.DomainOfTri=m2.elementity{3};
|
|
mesh.NbrTri=length(mesh.Tri);
|
|
mesh.ConnOfTri=zeros(mesh.NbrTri,2);
|
|
|
|
for i=1:mesh.NbrTri
|
|
[is,index]=ismember(mesh.Tri(i,:),el_ed2no_array,"rows");
|
|
mesh.ConnOfTri(i,1)=fix((index-1)/4)+1;
|
|
mesh.ConnOfTri(i,2)=index-(mesh.ConnOfTri(i,1)-1)*4;
|
|
end
|
|
|
|
mesh.NormOfFace=zeros(14,3);
|
|
mesh.NormOfFace(1,:)=[-1,0,0];
|
|
mesh.NormOfFace(2,:)=[0,-1,0];
|
|
mesh.NormOfFace(3,:)=[0,0,-1];
|
|
mesh.NormOfFace(5,:)=[0,1,0];
|
|
mesh.NormOfFace(14,:)=[1,0,0];
|
|
mesh.NormOfFace(4,:)=[0,0,1];
|
|
|
|
save('EPDmesh.mat',"mesh"); |