femcreate_new_v1.m

The Schumann resonance frequencies are a set of eigenmodes in the extremely low frequency (ELF) port

function fem = femCreate_new_v1(x,y,lambda,glob_epsR,epsR,sizeSB,sizePML,globMesh,subMesh,SymOrder,solve)
% COMSOL Multiphysics Model M-file
% Generated by COMSOL 3.4 (COMSOL 3.4.0.248, $Date: 2007/10/10 16:07:51 $)
tol=1e-12;

flclear fem

fem = InitFem('InPlaneWaves','rfwh','TM',1);

%*************  Geometrie Buildup
%keep and use this
clear s c p
s.objs={};s.name={};s.tags={};  %solid
c.objs={};c.name={};c.tags={};  %curve
p.objs={};p.name={};p.tags={};  %point

s = AddRect(s,[-sizeSB/2,sizeSB/2],[sizeSB/2,0],'R1','r1');

%pml
s = AddRect(s,[-sizeSB/2+sizePML,sizeSB/2-sizePML],[sizeSB/2-sizePML,0],'Area','area');

s = AddRect(s,[-sizeSB/2,sizeSB/2],[-sizeSB/2+sizePML,sizeSB/2-sizePML],'CornerL','cornL');
s = AddRect(s,[sizeSB/2-sizePML,sizeSB/2-sizePML],[sizeSB/2,sizeSB/2],'CornerR','cornR');

alpha=-pi/SymOrder;
xp=[];
yp=[];
xf=fliplr(x);
yf=-fliplr(y);
x=[x xf(2:end)*cos(alpha*2)-yf(2:end)*sin(alpha*2)];
y=[y xf(2:end)*sin(alpha*2)+yf(2:end)*cos(alpha*2)];

xp=x;
yp=y;
for n=1:(SymOrder/2-1)
    xp=[xp  x(2:end)*cos(alpha*2*n)-y(2:end)*sin(alpha*2*n)];
    yp=[yp  x(2:end)*sin(alpha*2*n)+y(2:end)*cos(alpha*2*n)];
end

s = AddPoly(s,xp,yp,'Particle','particle');

fem.draw=struct('s',s,'c',c);

fem.geom=geomcsg(fem);


pmlSDx=[GetSubDomNr(fem,-sizeSB/2+sizePML/2,sizeSB/4) GetSubDomNr(fem,sizeSB/2-sizePML/2,sizeSB/4)];
pmlSDy=GetSubDomNr(fem,0,sizeSB/2-sizePML/4);
pmlSDxy=[GetSubDomNr(fem,-sizeSB/2+sizePML/2,sizeSB/2-sizePML/2) GetSubDomNr(fem,sizeSB/2-sizePML/2,sizeSB/2-sizePML/2)];

particleSD=GetSubDomNr(fem,0,-x(1)/2);


clear bnd
bnd = InitBoundaryAppl(fem,'H0',{{0;0;0}});
bnd = SetBoundaryCond(bnd,'E0','PEC',GetSubDomainBounds(fem.geom,0));
bnd = SetBoundaryCond(bnd,'H','excit',[1 3],'H0',{{0;0;'1'}});

fem.appl{1}.bnd=bnd;

equ= InitSubDomProp(fem,'epsilonr',{glob_epsR},'Stype',{'none'},'coordOn',{{0;0}});

equ = SetSubDomProp(equ,'MatParticle',particleSD,'epsilonr',{'epsR'});
equ = SetSubDomProp(equ,'pmlx',pmlSDx,'epsilonr',{1},'Stype',{'coord'},'coordOn',{{1;0}});
equ = SetSubDomProp(equ,'pmly',pmlSDy,'epsilonr',{1},'Stype',{'coord'},'coordOn',{{0;1}});
equ = SetSubDomProp(equ,'pmlxy',pmlSDxy,'epsilonr',{1},'Stype',{'coord'},'coordOn',{{1;1}});

fem.appl{1}.equ = equ;

fem.const = {'lambda0',num2str(lambda,'%14.13e'),'epsR',num2str(epsR)};

fem.appl{1}.var = {'lambda0','lambda0'};

% Multiphysics
fem=multiphysics(fem);

fem=InitMeshFem(fem,globMesh,particleSD,subMesh);

if (solve>0)
    % Solve problem
    disp('Solving...');
    fem=SolveFemStatic(fem);
end