74 lines
2.6 KiB
Matlab
74 lines
2.6 KiB
Matlab
clc
|
||
clear all
|
||
% Ridge waveguide port scattering — Fem4 baseline (MeshData2x.mat)
|
||
% 不修改 MainFem2.m / MainFemx.m；本脚本专用于 MeshData2x（Fem4, L=4um）
|
||
|
||
rootDir = fileparts(mfilename('fullpath'));
|
||
meshFile = fullfile(rootDir, 'MeshData2x.mat');
|
||
outDir = fullfile(rootDir, 'OutFile_fem4');
|
||
if ~exist(outDir, 'dir')
|
||
mkdir(outDir);
|
||
end
|
||
|
||
tic
|
||
S = load(meshFile);
|
||
Mesh = S.Mesh;
|
||
if isfield(S, 'Physic')
|
||
Physic = S.Physic;
|
||
Material = S.Material;
|
||
[~, ~, Port1Solver, Port2Solver] = configPhysicMaterialFem4();
|
||
else
|
||
[Physic, Material, Port1Solver, Port2Solver] = configPhysicMaterialFem4();
|
||
end
|
||
if ~isfield(Mesh, 'NbrNodes')
|
||
Mesh.NbrNodes = length(Mesh.Nodes);
|
||
Mesh.NbrElements = length(Mesh.Elements);
|
||
end
|
||
disp('read mesh (MeshData2x / Fem4):');
|
||
toc
|
||
|
||
% Port eigenmode solve
|
||
tic
|
||
Port1Mesh = GetBMesh(Mesh.Nodes, Mesh.Faces, Mesh.FacesIndex, Physic.Port1);
|
||
[Port1Solver, Port1Mesh] = BoundaryEigenMode(Port1Mesh, Physic, Material, Port1Solver);
|
||
Port1Solver = GetPowerCoef(Port1Solver, Port1Mesh, Material.murF, 1);
|
||
Plot2Dfield(Port1Mesh, Port1Solver.modeNum, ...
|
||
Port1Solver.normE*Port1Solver.powerCoef, Port1Solver.neff, 1, 0.5);
|
||
Port1Mesh.normal = [0; 0; -1];
|
||
|
||
Port2Mesh = GetBMesh(Mesh.Nodes, Mesh.Faces, Mesh.FacesIndex, Physic.Port2);
|
||
[Port2Solver, Port2Mesh] = BoundaryEigenMode(Port2Mesh, Physic, Material, Port2Solver);
|
||
Port2Mesh.normal = [0; 0; 1];
|
||
disp('port eigenmode solve:');
|
||
toc
|
||
|
||
% Global FEM
|
||
[Solver, Mesh] = FemMatrixAssembly(Physic, Material, Mesh, ...
|
||
Port1Mesh, Port1Solver, Port2Mesh, Port2Solver);
|
||
|
||
AllFaces = GetBMesh(Mesh.Nodes, Mesh.Faces, Mesh.FacesIndex, ...
|
||
[Physic.Port1, Physic.Port2, Physic.PEC]);
|
||
|
||
tic
|
||
[~, ~, ~, normE_elem] = GetExEyEz2(1, Mesh.Nodes, Mesh.Elements, ...
|
||
Mesh.EdgesOfElements, Solver.Et);
|
||
plotE(Mesh.Nodes, AllFaces.OldFaces, Solver.normE, 3, 0.5);
|
||
In = plotE2(Mesh.Nodes, Port1Mesh.OldFaces, Mesh.AllFaces, ...
|
||
Mesh.FacesOfElements, normE_elem, 4, 0.5);
|
||
Out = plotE2(Mesh.Nodes, Port2Mesh.OldFaces, Mesh.AllFaces, ...
|
||
Mesh.FacesOfElements, normE_elem, 5, 0.5);
|
||
Scale = Out/In;
|
||
S11 = Solver.e1/abs(Port1Solver.powerCoef);
|
||
S21 = Solver.e2/abs(Port1Solver.powerCoef);
|
||
disp('post-process / plot:');
|
||
toc
|
||
|
||
fprintf('\n=== Fem4 port results (MeshData2x) ===\n');
|
||
fprintf('Port1 neff = %g + %gi\n', real(Port1Solver.neff(1)), imag(Port1Solver.neff(1)));
|
||
fprintf('powerCoef = %g\n', Port1Solver.powerCoef);
|
||
fprintf('S11 = %g + %gi\n', real(S11), imag(S11));
|
||
fprintf('S21 = %g + %gi\n', real(S21), imag(S21));
|
||
fprintf('|S11| = %g, |S21| = %g, Out/In = %g\n', abs(S11), abs(S21), Scale);
|
||
|
||
savePortResults(outDir, Mesh, Solver, S11, S21, Port1Solver, Port2Solver);
|