XIAN-FEM-2026June/三维matlab代码/2023-2-端口激励问题(四面体网格)/MakeData.m

83 lines
2.7 KiB
Matlab

clc
clear all
%description: utilize COMSOL to produce the mesh data
%read comsol project and get mesh data
model=mphload('Fem5xA2.mph');
[~, m2]=mphmeshstats(model);
%processing mesh data
Mesh.Nodes=m2.vertex';%vertex position
Elements=m2.elem{2}'+1;
Mesh.Elements = sort(Elements,2);%mesh connectivity n*3
Mesh.Domains=m2.elementity{2};%domain index n*1
Faces=m2.elem{3}'+1;%boundary face
Mesh.Faces = sort(Faces,2);
Mesh.FacesIndex=m2.elementity{3};%boundary face index
%材料参数
epsilon=[1 4 11.9 11.9 11.9 11.9 11.9 11.9 1]';
mu=[1.000021 1 1 1 1 1 1 1 1]';
sigma=[38000000 0 0 10 124603 4437.4 80290 1583.02 0]';
Material.epsilonr=ones(24,1);
Material.mur=ones(24,1);
Material.sigma=ones(24,1);
%材料赋值
material=cell(9,1);
material{1}=[5 9 19 20 24];%Ai
material{2}=[2 4 6 10 11 17 21 23];%SiO2
material{3}=[3 14 22];%Si
material{4}=1;%Si_sub
material{5}=18;%Si_N_1e20
material{6}=[15 16];%Si_N_1e18
material{7}=8;%Si_P_1e20
material{8}=[12 13];%Si_P_5e17
material{9}=7;%air
for i=1:9
tempIndex=material{i};
Material.epsilonr(tempIndex(1))=epsilon(i);
Material.sigma(tempIndex(1))=sigma(i);
Material.mur(tempIndex(1))=mu(i);
for j=1:length(tempIndex)
Material.epsilonr(tempIndex(j))=Material.epsilonr(tempIndex(1));
Material.mur(tempIndex(j))=Material.mur(tempIndex(1));
Material.sigma(tempIndex(j))=Material.sigma(tempIndex(1));
end
end
Material.epsilonrF=ones(174,1);
Material.murF=ones(174,1);
Material.sigmaF=ones(174,1);
material=cell(9,1);
material{1}=[19 37 87 92 112,...
20 38 88 93 113];%Ai
material{2}=[7 15 23 42 46 77 97 107,...
8 16 24 43 47 78 98 108];%SiO2
material{3}=[11 63 102,...
12 64 103];%Si
material{4}=[3 4];%Si_sub
material{5}=[82 83];%Si_N_1e20
material{6}=[68 73 69 74];%Si_N_1e18
material{7}=[32 33];%Si_P_1e20
material{8}=[52 57 53 58];%Si_P_5e17
material{9}=[27 28];%air
for i=1:9
tempIndex=material{i};
Material.epsilonrF(tempIndex(1))=epsilon(i);
Material.murF(tempIndex(1))=mu(i);
Material.sigmaF(tempIndex(1))=sigma(i);
for j=1:length(tempIndex)
Material.epsilonrF(tempIndex(j))=Material.epsilonrF(tempIndex(1));
Material.murF(tempIndex(j))=Material.murF(tempIndex(1));
Material.sigmaF(tempIndex(j))=Material.sigmaF(tempIndex(1));
end
end
Physic.Port1=[3, 7, 11, 15, 19, 23, 27, 32, 37, 42, 46, 52, 57, 63, 68, 73, 77, 82, 87, 92, 97, 102, 107, 112];%端口1 入射端
Physic.Port2=[4, 8, 12, 16, 20, 24, 28, 33, 38, 43, 47, 53, 58, 64, 69, 74, 78, 83, 88, 93, 98, 103, 108, 113];%端口2 出射端
Physic.PEC=[1, 2, 5, 9, 13, 17, 21, 25, 29, 115, 116, 117, 118, 119, 120, 121];%PEC
%save mesh data as a file
save('MeshDataxxxx','Mesh','Physic','Material');