caiyu/openEMS
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

simplified inf dipole example

Browse Source
This commit is contained in:
Thorsten Liebig 2011-11-28 14:12:16 +01:00
parent 43179fb6e8
commit 5b8dfba6cf
1 changed files with 46 additions and 90 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

136
matlab/examples/NF2FF/infDipol.m
View File

@ -1,4 +1,3 @@
function infDipol
% %
% infinitesimal dipole example % infinitesimal dipole example
% %
@ -6,110 +5,71 @@ function infDipol
close all close all
clear clear
clc clc
drawnow
postprocessing_only = 0;
postprocessing_only = 1;
global g
setup
dipole_orientation = 3; % 1,2,3: x,y,z
CSX = createStructure(dipole_orientation);
if ~postprocessing_only
writeCSX( CSX );
CSXGeomPlot( [g.Sim_Path '/' g.Sim_CSX] )
run;
end
postprocess;
function setup
global g
physical_constants physical_constants
% setup the simulation % setup the simulation
g.drawingunit = 1e-6; % specify everything in um drawingunit = 1e-6; % specify everything in um
g.Sim_Path = 'tmp'; Sim_Path = 'tmp';
g.Sim_CSX = 'tmp.xml'; Sim_CSX = 'tmp.xml';
g.f_max = 1e9; f_max = 1e9;
g.lambda = c0/g.f_max; lambda = c0/f_max;
% setup geometry values % setup geometry values
g.dipole_length = g.lambda/50 /g.drawingunit; dipole_length = lambda/50 /drawingunit;
dipole_orientation = 3; % 1,2,3: x,y,z
function CSX = createStructure(dipole_orientation)
global g
physical_constants
CSX = InitCSX(); CSX = InitCSX();
% create an equidistant mesh % create an equidistant mesh
mesh.x = -g.dipole_length*10:g.dipole_length/2:g.dipole_length*10; mesh.x = -dipole_length*10:dipole_length/2:dipole_length*10;
mesh.y = -g.dipole_length*10:g.dipole_length/2:g.dipole_length*10; mesh.y = -dipole_length*10:dipole_length/2:dipole_length*10;
mesh.z = -g.dipole_length*10:g.dipole_length/2:g.dipole_length*10; mesh.z = -dipole_length*10:dipole_length/2:dipole_length*10;
mesh = AddPML( mesh, [8 8 8 8 8 8] ); % add space for PML
CSX = DefineRectGrid( CSX, g.drawingunit, mesh );
% excitation % excitation
ex_vector = [0 0 0]; ex_vector = [0 0 0];
ex_vector(dipole_orientation) = 1; ex_vector(dipole_orientation) = 1;
start = ex_vector * -g.dipole_length/2; start = ex_vector * -dipole_length/2;
stop = ex_vector * g.dipole_length/2; stop = ex_vector * dipole_length/2;
CSX = AddExcitation( CSX, 'infDipole', 1, ex_vector ); CSX = AddExcitation( CSX, 'infDipole', 1, ex_vector );
CSX = AddBox( CSX, 'infDipole', 1, start, stop ); CSX = AddBox( CSX, 'infDipole', 1, start, stop );
% NFFF contour % NFFF contour
s1 = [-4.5, -4.5, -4.5] * g.dipole_length/2; start = [mesh.x(1) mesh.y(1) mesh.z(1) ]
s2 = [ 4.5, 4.5, 4.5] * g.dipole_length/2; stop = [mesh.x(end) mesh.y(end) mesh.z(end) ]
[CSX g.nf2ff] = CreateNF2FFBox(CSX, 'nf2ff', s1, s2); [CSX nf2ff] = CreateNF2FFBox(CSX, 'nf2ff', start, stop);
% add space for PML
mesh = AddPML( mesh, [8 8 8 8 8 8] );
% define the mesh
CSX = DefineRectGrid( CSX, drawingunit, mesh );
if ~postprocessing_only
% setup FDTD parameters & excitation function
max_timesteps = 2000;
min_decrement = 1e-6;
FDTD = InitFDTD( max_timesteps, min_decrement, 'OverSampling',10 );
FDTD = SetGaussExcite( FDTD, f_max/2, f_max/2 );
BC = {'PML_8' 'PML_8' 'PML_8' 'PML_8' 'PML_8' 'PML_8'};
FDTD = SetBoundaryCond( FDTD, BC );
% Write openEMS compatible xml-file
[~,~,~] = rmdir(Sim_Path,'s');
[~,~,~] = mkdir(Sim_Path);
WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
% define openEMS options and start simulation
openEMS_opts = '';
RunOpenEMS( Sim_Path, Sim_CSX, openEMS_opts );
end
%% post processing
function writeCSX(CSX)
global g
% setup FDTD parameters & excitation function
max_timesteps = 2000;
min_decrement = 1e-6;
FDTD = InitFDTD( max_timesteps, min_decrement, 'OverSampling',10 );
FDTD = SetGaussExcite( FDTD, g.f_max/2, g.f_max/2 );
BC = {'PML_8' 'PML_8' 'PML_8' 'PML_8' 'PML_8' 'PML_8'};
FDTD = SetBoundaryCond( FDTD, BC );
% Write openEMS compatible xml-file
[~,~,~] = rmdir(g.Sim_Path,'s');
[~,~,~] = mkdir(g.Sim_Path);
WriteOpenEMS([g.Sim_Path '/' g.Sim_CSX],FDTD,CSX);
function run
global g
% define openEMS options and start simulation
openEMS_opts = '';
openEMS_opts = [openEMS_opts ' --engine=fastest'];
RunOpenEMS( g.Sim_Path, g.Sim_CSX, openEMS_opts );
function postprocess
global g
disp( ' ' ); disp( ' ' );
disp( ' ********************************************************** ' ); disp( ' ********************************************************** ' );
disp( ' ' ); disp( ' ' );
@ -117,15 +77,13 @@ disp( ' ' );
% calculate the far field at phi=0 degrees and at phi=90 degrees % calculate the far field at phi=0 degrees and at phi=90 degrees
thetaRange = 0:2:359; thetaRange = 0:2:359;
r = 1; % evaluate fields at radius r r = 1; % evaluate fields at radius r
disp( 'calculating far field at phi=[0 90] deg...' ); disp( 'calculating far field at phi=[0 90] deg..' );
[E_far_theta,E_far_phi,Prad,Dmax] = AnalyzeNF2FF( g.Sim_Path, g.nf2ff, g.f_max, thetaRange, [0 90], r ); [E_far_theta,E_far_phi,Prad,Dmax] = AnalyzeNF2FF( Sim_Path, nf2ff, f_max, thetaRange, [0 90], r );
% display power and directivity % display power and directivity
disp( ['radiated power: Prad = ' num2str(Prad)] ); disp( ['radiated power: Prad = ' num2str(Prad)] );
disp( ['directivity: Dmax = ' num2str(Dmax)] ); disp( ['directivity: Dmax = ' num2str(Dmax)] );
% calculate the e-field magnitude for phi = 0 deg % calculate the e-field magnitude for phi = 0 deg
E_phi0_far = zeros(1,numel(thetaRange)); E_phi0_far = zeros(1,numel(thetaRange));
for n=1:numel(thetaRange) for n=1:numel(thetaRange)
@ -139,7 +97,6 @@ ylabel( 'theta / deg' );
title( ['electrical far field (V/m) @r=' num2str(r) ' m phi=0 deg'] ); title( ['electrical far field (V/m) @r=' num2str(r) ' m phi=0 deg'] );
legend( 'e-field magnitude', 'Location', 'BestOutside' ); legend( 'e-field magnitude', 'Location', 'BestOutside' );
% calculate the e-field magnitude for phi = 90 deg % calculate the e-field magnitude for phi = 90 deg
E_phi90_far = zeros(1,numel(thetaRange)); E_phi90_far = zeros(1,numel(thetaRange));
for n=1:numel(thetaRange) for n=1:numel(thetaRange)
@ -153,14 +110,11 @@ ylabel( 'theta / deg' );
title( ['electrical far field (V/m) @r=' num2str(r) ' m phi=90 deg'] ); title( ['electrical far field (V/m) @r=' num2str(r) ' m phi=90 deg'] );
legend( 'e-field magnitude', 'Location', 'BestOutside' ); legend( 'e-field magnitude', 'Location', 'BestOutside' );
% calculate the far field at theta=90 degrees % calculate the far field at theta=90 degrees
phiRange = 0:2:359; phiRange = 0:2:359;
r = 1; % evaluate fields at radius r r = 1; % evaluate fields at radius r
disp( 'calculating far field at theta=90 deg...' ); disp( 'calculating far field at theta=90 deg..' );
[E_far_theta,E_far_phi] = AnalyzeNF2FF( g.Sim_Path, g.nf2ff, g.f_max, 90, phiRange, r ); [E_far_theta,E_far_phi] = AnalyzeNF2FF( Sim_Path, nf2ff, f_max, 90, phiRange, r );
E_theta90_far = zeros(1,numel(phiRange)); E_theta90_far = zeros(1,numel(phiRange));
for n=1:numel(phiRange) for n=1:numel(phiRange)
@ -175,14 +129,12 @@ title( ['electrical far field (V/m) @r=' num2str(r) ' m theta=90 deg'] );
legend( 'e-field magnitude', 'Location', 'BestOutside' ); legend( 'e-field magnitude', 'Location', 'BestOutside' );
% calculate 3D pattern % calculate 3D pattern
phiRange = 0:15:360; phiRange = 0:15:360;
thetaRange = 0:10:180; thetaRange = 0:10:180;
r = 1; % evaluate fields at radius r r = 1; % evaluate fields at radius r
disp( 'calculating 3D far field...' ); disp( 'calculating 3D far field...' );
[E_far_theta,E_far_phi] = AnalyzeNF2FF( g.Sim_Path, g.nf2ff, g.f_max, thetaRange, phiRange, r ); [E_far_theta,E_far_phi] = AnalyzeNF2FF( Sim_Path, nf2ff, f_max, thetaRange, phiRange, r );
E_far = sqrt( abs(E_far_theta).^2 + abs(E_far_phi).^2 ); E_far = sqrt( abs(E_far_theta).^2 + abs(E_far_phi).^2 );
E_far_normalized = E_far / max(E_far(:)); E_far_normalized = E_far / max(E_far(:));
[theta,phi] = ndgrid(thetaRange/180*pi,phiRange/180*pi); [theta,phi] = ndgrid(thetaRange/180*pi,phiRange/180*pi);
@ -195,3 +147,7 @@ axis equal
xlabel( 'x' ); xlabel( 'x' );
ylabel( 'y' ); ylabel( 'y' );
zlabel( 'z' ); zlabel( 'z' );
%
DumpFF2VTK([Sim_Path '/FF_pattern.vtk'],E_far_normalized, thetaRange, phiRange);
disp(['view the farfield pattern "' Sim_Path '/FF_pattern.vtk" using paraview' ]);