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

testsuite: new test: cavity.m (resonance frequency)

pull/1/head
Sebastian Held 2010-03-31 15:39:46 +02:00
parent b5bea21529
commit 09bbf49b1a
2 changed files with 123 additions and 32 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

124
TESTSUITE/combinedtests/cavity.m
View File

@ -7,9 +7,18 @@ ENABLE_PLOTS = 1;
CLEANUP = 0; % if enabled and result is PASS, remove simulation folder
STOP_IF_FAILED = 1; % if enabled and result is FAILED, stop with error
% LIMITS
upper_error = 0.036; % max +3.6%
lower_error = 0; % max -0.0%
% LIMITS - inside
lower_rel_limit = 1.3e-3; % -0.13%
upper_rel_limit = 1.3e-3; % +0.13%
lower_rel_limit_TM = 2.5e-3; % -0.25%
upper_rel_limit_TM = 0; % +0%
min_rel_amplitude = 0.6; % 60%
min_rel_amplitude_TM = 0.27; % 27%
% LIMITS - outside
outer_rel_limit = 0.02;
max_rel_amplitude = 0.17;
% structure
a = 5e-2;
@ -19,7 +28,7 @@ if ~((b<a) && (a<d))
error 'correct the dimensions of the cavity'
end
f_start = 0;
f_start = 1e9;
f_stop = 10e9;
Sim_Path = 'tmp';
@ -35,23 +44,28 @@ FDTD = SetBoundaryCond(FDTD,BC);
%setup CSXCAD geometry
CSX = InitCSX();
mesh.x = linspace(0,a,35);
mesh.y = linspace(0,b,35);
mesh.z = linspace(0,d,35);
grid_res = 2e-3;
mesh.x = 0:grid_res:a; %linspace(0,a,25);
mesh.y = 0:grid_res:b; %linspace(0,b,25);
mesh.z = 0:grid_res:d; %linspace(0,d,25);
CSX = DefineRectGrid(CSX, 1,mesh);
% excitation
pos1 = [mesh.x(10) mesh.y(10) mesh.z(10)];
pos2 = [mesh.x(10) mesh.y(11) mesh.z(10)];
CSX = AddExcitation(CSX,'excite',1,[0 1 0]);
CSX = AddBox(CSX, 'excite', 1, pos1, pos2);
CSX = AddExcitation(CSX,'excite1',0,[1 1 1]);
p(1,1) = mesh.x(floor(end*2/3));
p(2,1) = mesh.y(floor(end*2/3));
p(3,1) = mesh.z(floor(end*2/3));
p(1,2) = mesh.x(floor(end*2/3)+1);
p(2,2) = mesh.y(floor(end*2/3)+1);
p(3,2) = mesh.z(floor(end*2/3)+1);
CSX = AddCurve( CSX, 'excite1', 0, p );
% %dump
%dump
% CSX = AddDump(CSX,'Et_',0,2);
% pos1 = [mesh.x(1) mesh.y(10) mesh.z(1)];
% pos2 = [mesh.x(end) mesh.y(10) mesh.z(end)];
% pos1 = [mesh.x(1) mesh.y(1) mesh.z(1)];
% pos2 = [mesh.x(end) mesh.y(end) mesh.z(end)];
% CSX = AddBox(CSX,'Et_',0 , pos1,pos2);
%
% %dump
% CSX = AddDump(CSX,'Et2_',0,2);
% pos1 = [mesh.x(1) mesh.y(1) mesh.z(1)];
@ -65,13 +79,19 @@ CSX = AddBox(CSX, 'excite', 1, pos1, pos2);
% CSX = AddBox(CSX,'Et3_',0 , pos1,pos2);
%voltage calc
pos1 = [mesh.x(15) mesh.y(15) mesh.z(15)];
pos2 = [mesh.x(15) mesh.y(16) mesh.z(15)];
CSX = AddProbe(CSX,'ut1x',0);
pos1 = [mesh.x(floor(end/4)) mesh.y(floor(end/2)) mesh.z(floor(end/5))];
pos2 = [mesh.x(floor(end/4)+1) mesh.y(floor(end/2)) mesh.z(floor(end/5))];
CSX = AddBox(CSX,'ut1x', 0 ,pos1,pos2);
CSX = AddProbe(CSX,'ut1y',0);
pos1 = [mesh.x(floor(end/4)) mesh.y(floor(end/2)) mesh.z(floor(end/5))];
pos2 = [mesh.x(floor(end/4)) mesh.y(floor(end/2)+1) mesh.z(floor(end/5))];
CSX = AddBox(CSX,'ut1y', 0 ,pos1,pos2);
pos1 = [mesh.x(floor(end/2)) mesh.y(floor(end/2)) mesh.z(end)];
pos2 = [mesh.x(floor(end/2)) mesh.y(floor(end/2)) mesh.z(end-1)];
CSX = AddProbe(CSX,'ut1z',0);
pos1 = [mesh.x(floor(end/2)) mesh.y(floor(end/2)) mesh.z(floor(end/5))];
pos2 = [mesh.x(floor(end/2)) mesh.y(floor(end/2)) mesh.z(floor(end/5)+1)];
CSX = AddBox(CSX,'ut1z', 0 ,pos1,pos2);
%Write openEMS compatible xml-file
@ -90,6 +110,16 @@ cd(savePath);
% analysis
%
% remove excitation from time series
t_start = 7e-10; % FIXME to be calculated
t_idx_start = interp1( UI.TD{1}.t, 1:numel(UI.TD{1}.t), t_start, 'nearest' );
for n=1:numel(UI.TD)
UI.TD{n}.t = UI.TD{n}.t(t_idx_start:end);
UI.TD{n}.val = UI.TD{n}.val(t_idx_start:end);
[UI.FD{n}.f,UI.FD{n}.val] = FFT_time2freq( UI.TD{n}.t, UI.TD{n}.val );
end
f = UI.FD{1}.f;
ux = UI.FD{1}.val;
uy = UI.FD{2}.val;
@ -106,14 +136,31 @@ uz = uz(f_idx_start:f_idx_stop);
k = @(m,n,l) sqrt( (m*pi/a)^2 + (n*pi/b)^2 + (l*pi/d)^2 );
f_TE101 = c0/(2*pi) * k(1,0,1);
f_TE102 = c0/(2*pi) * k(1,0,2);
f_TE103 = c0/(2*pi) * k(1,0,3);
f_TE201 = c0/(2*pi) * k(2,0,1);
f_TE202 = c0/(2*pi) * k(2,0,2);
f_TM110 = c0/(2*pi) * k(1,1,0);
f_TM210 = c0/(2*pi) * k(2,1,0);
f_TM111 = c0/(2*pi) * k(1,1,1);
f_TE = [f_TE101 f_TE102 f_TE201 f_TE202];
f_TM = [f_TM110 f_TM111];
% calculate frequency limits
temp = [f_start f_TE f_stop];
f_outer1 = [];
f_outer2 = [];
for n=1:numel(temp)-1
f_outer1 = [f_outer1 temp(n) .* (1+outer_rel_limit)];
f_outer2 = [f_outer2 temp(n+1) .* (1-outer_rel_limit)];
end
temp = [f_start f_TM f_stop];
f_outer1_TM = [];
f_outer2_TM = [];
for n=1:numel(temp)-1
f_outer1_TM = [f_outer1_TM temp(n) .* (1+outer_rel_limit)];
f_outer2_TM = [f_outer2_TM temp(n+1) .* (1-outer_rel_limit)];
end
f_TE = [f_TE101 f_TE102 f_TE103 f_TE201 f_TE202];
f_TM = [f_TM110 f_TM210];
if ENABLE_PLOTS
figure
@ -121,25 +168,38 @@ if ENABLE_PLOTS
max1 = max(abs(uy));
hold on
plot( repmat(f_TE,2,1)/1e9, repmat([0; max1],1,numel(f_TE)), 'm-.', 'LineWidth', 2 )
plot( (repmat(f_TE,2,1) .* repmat(1-lower_rel_limit,2,numel(f_TE)))/1e9, repmat([0; max1],1,numel(f_TE)), 'r-', 'LineWidth', 1 )
plot( (repmat(f_TE,2,1) .* repmat(1+upper_rel_limit,2,numel(f_TE)))/1e9, repmat([0; max1],1,numel(f_TE)), 'r-', 'LineWidth', 1 )
plot( (repmat(f_TE,2,1) .* repmat([1-outer_rel_limit;1+outer_rel_limit],1,numel(f_TE)))/1e9, repmat(max1*min_rel_amplitude,2,numel(f_TE)), 'r-', 'LineWidth', 1 ) % freq limits
plot( [f_outer1;f_outer2]/1e9, repmat(max1*max_rel_amplitude,2,numel(f_outer1)), 'g-', 'LineWidth', 1 ) % amplitude limits
xlabel('Frequency (GHz)')
legend( {'u_y','theoretical'} );
legend( {'u_y','theoretical'} )
title( 'TE-modes' )
figure
plot(f/1e9,abs(uz))
max1 = max(abs(uz));
hold on
plot( repmat(f_TM,2,1)/1e9, repmat([0; max1],1,numel(f_TM)), 'm-.', 'LineWidth', 2 )
plot( (repmat(f_TM,2,1) .* repmat(1-lower_rel_limit_TM,2,numel(f_TM)))/1e9, repmat([0; max1],1,numel(f_TM)), 'r-', 'LineWidth', 1 )
plot( (repmat(f_TM,2,1) .* repmat(1+upper_rel_limit_TM,2,numel(f_TM)))/1e9, repmat([0; max1],1,numel(f_TM)), 'r-', 'LineWidth', 1 )
plot( (repmat(f_TM,2,1) .* repmat([1-lower_rel_limit_TM;1+upper_rel_limit_TM],1,numel(f_TM)))/1e9, repmat(max1*min_rel_amplitude_TM,2,numel(f_TM)), 'r-', 'LineWidth', 1 ) % freq limits
plot( [f_outer1_TM;f_outer2_TM]/1e9, repmat(max1*max_rel_amplitude,2,numel(f_outer1_TM)), 'g-', 'LineWidth', 1 ) % amplitude limits
xlabel('Frequency (GHz)')
legend( {'u_z','theoretical'} );
legend( {'u_z','theoretical'} )
title( 'TM-modes' )
end
pass = 1;
% pass = check_limits( Z, upper_limit, lower_limit );
% if pass
% disp( 'combinedtests/Coax.m (characteristic impedance): pass' );
% else
% disp( 'combinedtests/Coax.m (characteristic impedance): * FAILED *' );
% end
pass1 = check_frequency( f, abs(uy), f_TE*(1+upper_rel_limit), f_TE*(1-lower_rel_limit), min_rel_amplitude, 'inside' );
pass2 = check_frequency( f, abs(uz), f_TM*(1+upper_rel_limit_TM), f_TM*(1-lower_rel_limit_TM), min_rel_amplitude_TM, 'inside' );
pass3 = check_frequency( f, abs(uy), f_outer2, f_outer1, max_rel_amplitude, 'outside' );
pass4 = check_frequency( f, abs(uz), f_outer2_TM, f_outer1_TM, max_rel_amplitude, 'outside' );
pass = pass1 && pass2 && pass3 && pass4;
if pass
disp( 'combinedtests/cavity.m (resonance frequency): pass' );
else
disp( 'combinedtests/cavity.m (resonance frequency): * FAILED *' );
end

31
TESTSUITE/helperscripts/check_frequency.m Normal file
View File

@ -0,0 +1,31 @@
function pass = check_frequency( f, val, f_upper, f_lower, rel_amplitude, type )
pass = true;
max1 = max(val);
if numel(f_upper) ~= numel(f_lower)
error 'inconsistant vectors'
end
for n=1:numel(f_upper)
f1 = f_lower(n);
f2 = f_upper(n);
f1_idx = interp1( f, 1:numel(f), f1, 'nearest' );
% if f(f1_idx) < f1, f1_idx = f1_idx + 1; end
f2_idx = interp1( f, 1:numel(f), f2, 'nearest' );
% if f(f2_idx) > f2, f2_idx = f2_idx - 1; end
if strcmp( type, 'inside' )
if max( val(f1_idx:f2_idx) ) < max1 * rel_amplitude
pass = false;
return
end
elseif strcmp( type, 'outside' )
if max( val(f1_idx:f2_idx) ) > max1 * rel_amplitude
pass = false;
return
end
else
error 'unsupported operation'
end
end