110 lines
3.8 KiB
Matlab
110 lines
3.8 KiB
Matlab
function exportAbMatrixFem4(outDir)
|
|
%EXPORTABMATRIXFEM4 Export reduced A/b (after PEC, before symrcm) for Fem4 port case.
|
|
%
|
|
% exportAbMatrixFem4() % -> OutFile_fem4_ab/
|
|
% exportAbMatrixFem4('OutFile_fem4_ab')
|
|
%
|
|
% Output format matches C++ Test_OutputMatrix:
|
|
% Ai.txt, Aj.txt (0-based indices), Av.txt (re,im), Bv_real.txt, Bv_imag.txt
|
|
%
|
|
% Run from this folder:
|
|
% cd('三维matlab代码/2023-2-端口激励问题(四面体网格)')
|
|
% exportAbMatrixFem4
|
|
|
|
rootDir = fileparts(mfilename('fullpath'));
|
|
if nargin < 1 || isempty(outDir)
|
|
outDir = fullfile(rootDir, 'OutFile_fem4_ab');
|
|
end
|
|
if ~exist(outDir, 'dir')
|
|
mkdir(outDir);
|
|
end
|
|
|
|
exportDir = fullfile(rootDir, '..', 'matlab 3D一阶散射问题');
|
|
addpath(rootDir, exportDir);
|
|
|
|
meshFile = fullfile(rootDir, 'MeshData2x.mat');
|
|
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
|
|
|
|
fprintf('Port eigenmode solve (same as MainFem2x)...\n');
|
|
Port1Mesh = GetBMesh(Mesh.Nodes, Mesh.Faces, Mesh.FacesIndex, Physic.Port1);
|
|
[Port1Solver, Port1Mesh] = BoundaryEigenMode(Port1Mesh, Physic, Material, Port1Solver);
|
|
Port1Solver = GetPowerCoef(Port1Solver, Port1Mesh, Material.murF, 1);
|
|
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];
|
|
|
|
fprintf('Global matrix assembly (wave + ports + PEC reduce)...\n');
|
|
epsilonr = Material.epsilonr - 1i * Material.sigma * Material.temp;
|
|
[A, Mesh.Edges, Mesh.EdgesOfElements, Dof] = PhysicMatrixAssembly( ...
|
|
Mesh.Nodes, Mesh.Elements, Mesh.Domains, Physic.lam0, epsilonr, Material.mur);
|
|
b = complex(zeros(Dof, 1));
|
|
|
|
el2no = Mesh.Elements';
|
|
f1 = el2no([1 1 1 2], :);
|
|
f2 = el2no([2 2 3 3], :);
|
|
f3 = el2no([3 4 4 4], :);
|
|
f_array = [f1(:) f2(:) f3(:)];
|
|
[Faces, ~, FacesOfElements] = unique(f_array, 'rows');
|
|
Mesh.AllFaces = Faces;
|
|
Mesh.FacesOfElements = FacesOfElements;
|
|
addDof = 0;
|
|
|
|
Port1Mesh.FacesConn = GetBFacesConn(Port1Mesh.OldFaces, Faces, FacesOfElements);
|
|
[A, b, addDof] = PortBCMatrixAssembly(A, b, Mesh.Nodes, Mesh.Elements, Mesh.EdgesOfElements, ...
|
|
Mesh.Domains, Material.mur, Port1Mesh, Port1Solver, Dof, addDof);
|
|
|
|
Port2Mesh.FacesConn = GetBFacesConn(Port2Mesh.OldFaces, Faces, FacesOfElements);
|
|
[A, b, addDof] = PortBCMatrixAssembly(A, b, Mesh.Nodes, Mesh.Elements, Mesh.EdgesOfElements, ...
|
|
Mesh.Domains, Material.mur, Port2Mesh, Port2Solver, Dof, addDof);
|
|
addDof(1) = [];
|
|
|
|
[PECFaces, ~] = GetBFaces(Mesh.Faces, Mesh.FacesIndex, Physic.PEC);
|
|
fl2no = PECFaces';
|
|
n1 = fl2no([1 1 2], :);
|
|
n2 = fl2no([2 3 3], :);
|
|
fl_fd2no_array = [n1(:) n2(:)];
|
|
PECEdges = unique(fl_fd2no_array, 'rows');
|
|
PECInd = find(ismember(Mesh.Edges, PECEdges, 'rows'));
|
|
Free_ind = 1:(Dof + sum(addDof));
|
|
Free_ind(PECInd) = [];
|
|
A = A(Free_ind, Free_ind);
|
|
b = b(Free_ind);
|
|
|
|
fprintf('Reduced system: size=%d, nnz=%d, |b|=%g\n', ...
|
|
size(A, 1), nnz(A), norm(b));
|
|
|
|
export_ab_matrix(A, b, outDir);
|
|
|
|
metaFile = fullfile(outDir, 'meta.txt');
|
|
fid = fopen(metaFile, 'w');
|
|
fprintf(fid, 'source=MATLAB exportAbMatrixFem4\n');
|
|
fprintf(fid, 'mesh=MeshData2x.mat\n');
|
|
fprintf(fid, 'stage=after_PEC_before_symrcm\n');
|
|
fprintf(fid, 'n=%d\n', size(A, 1));
|
|
fprintf(fid, 'nnz=%d\n', nnz(A));
|
|
fprintf(fid, 'nPEC=%d\n', numel(PECInd));
|
|
fprintf(fid, 'Dof=%d\n', Dof);
|
|
fprintf(fid, 'addDof=%d\n', sum(addDof));
|
|
fprintf(fid, 'powerCoef=%.15g\n', Port1Solver.powerCoef);
|
|
fclose(fid);
|
|
|
|
fprintf('Done. Compare with C++ OutFile/Ai.txt, Aj.txt, Av.txt, Bv_*.txt\n');
|
|
fprintf(' MATLAB: %s\n', outDir);
|
|
fprintf(' C++: port/Release/OutFile/\n');
|
|
|
|
end
|