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

Merge branch 'master' into sse 

Conflicts:
	TESTSUITE/helperscripts/invoke_openEMS.m
pull/1/head
Sebastian Held 2010-04-30 15:29:21 +02:00
parent 6575e89252 d0427ea140
commit 209c066732
23 changed files with 380 additions and 200 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

13
FDTD/engine.cpp
View File

@ -118,14 +118,14 @@ void Engine::ApplyVoltageExcite()
exc_pos = (int)numTS - (int)Op->E_Exc_delay[n]; exc_pos = (int)numTS - (int)Op->E_Exc_delay[n];
exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->ExciteLength); exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->ExciteLength);
// if (n==0) cerr << numTS << " => " << Op->ExciteSignal[exc_pos] << endl; // if (n==0) cerr << numTS << " => " << Op->ExciteSignal[exc_pos] << endl;
volt[Op->E_Exc_dir[n]][Op->E_Exc_index[0][n]][Op->E_Exc_index[1][n]][Op->E_Exc_index[2][n]] += Op->E_Exc_amp[n]*Op->ExciteSignal[exc_pos]; volt[Op->E_Exc_dir[n]][Op->E_Exc_index[0][n]][Op->E_Exc_index[1][n]][Op->E_Exc_index[2][n]] += Op->E_Exc_amp[n]*Op->ExciteSignal_volt[exc_pos];
} }
// write the first excitation into the file "et1" // write the first excitation into the file "et1"
if (Op->E_Exc_Count >= 1) { if (Op->E_Exc_Count >= 1) {
exc_pos = (int)numTS - (int)Op->E_Exc_delay[0]; exc_pos = (int)numTS - (int)Op->E_Exc_delay[0];
exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->ExciteLength); exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->ExciteLength);
file_et1 << numTS * Op->GetTimestep() << "\t" << Op->E_Exc_amp[0]*Op->ExciteSignal[exc_pos] << "\n"; // do not use std::endl here, because it will do an implicit flush file_et1 << numTS * Op->GetTimestep() << "\t" << Op->E_Exc_amp[0]*Op->ExciteSignal_volt[exc_pos] << "\n"; // do not use std::endl here, because it will do an implicit flush
} }
} }
@ -157,6 +157,15 @@ void Engine::UpdateCurrents()
void Engine::ApplyCurrentExcite() void Engine::ApplyCurrentExcite()
{ {
int exc_pos;
//soft current excitation here (H-field excite)
for (unsigned int n=0;n<Op->Curr_Exc_Count;++n)
{
exc_pos = (int)numTS - (int)Op->Curr_Exc_delay[n];
exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->ExciteLength);
// if (n==0) cerr << numTS << " => " << Op->ExciteSignal[exc_pos] << endl;
curr[Op->Curr_Exc_dir[n]][Op->Curr_Exc_index[0][n]][Op->Curr_Exc_index[1][n]][Op->Curr_Exc_index[2][n]] += Op->Curr_Exc_amp[n]*Op->ExciteSignal_curr[exc_pos];
}
} }
bool Engine::IterateTS(unsigned int iterTS) bool Engine::IterateTS(unsigned int iterTS)

34
FDTD/engine_ext_mur_abc.cpp
View File

@ -31,10 +31,28 @@ Engine_Ext_Mur_ABC::Engine_Ext_Mur_ABC(Operator_Ext_Mur_ABC* op_ext) : Engine_Ex
m_LineNr = m_Op_mur->m_LineNr; m_LineNr = m_Op_mur->m_LineNr;
m_LineNr_Shift = m_Op_mur->m_LineNr_Shift; m_LineNr_Shift = m_Op_mur->m_LineNr_Shift;
m_Mur_Coeff = m_Op_mur->m_Mur_Coeff; m_Mur_Coeff_nyP = m_Op_mur->m_Mur_Coeff_nyP;
m_Mur_Coeff_nyPP = m_Op_mur->m_Mur_Coeff_nyPP;
m_volt_nyP = Create2DArray(m_numLines); m_volt_nyP = Create2DArray(m_numLines);
m_volt_nyPP = Create2DArray(m_numLines); m_volt_nyPP = Create2DArray(m_numLines);
//find if some excitation is on this mur-abc and find the max length of this excite, so that the abc can start after the excitation is done...
int maxDelay=-1;
for (unsigned int n=0;n<m_Op_mur->m_Op->E_Exc_Count;++n)
{
if ( ((m_Op_mur->m_Op->E_Exc_dir[n]==m_nyP) || (m_Op_mur->m_Op->E_Exc_dir[n]==m_nyPP)) && (m_Op_mur->m_Op->E_Exc_index[m_ny][n]==m_LineNr) )
{
if ((int)m_Op_mur->m_Op->E_Exc_delay[n]>maxDelay)
maxDelay = (int)m_Op_mur->m_Op->E_Exc_delay[n];
}
}
m_start_TS = 0;
if (maxDelay>=0)
{
m_start_TS = maxDelay + m_Op_mur->m_Op->ExciteLength + 10; //give it some extra timesteps, for the excitation to travel at least one cell away
cerr << "Engine_Ext_Mur_ABC::Engine_Ext_Mur_ABC: Warning: Excitation inside the Mur-ABC #" << m_ny << "-" << (int)(m_LineNr>0) << " found!!!! Mur-ABC will be switched on after excitation is done at " << m_start_TS << " timesteps!!! " << endl;
}
} }
Engine_Ext_Mur_ABC::~Engine_Ext_Mur_ABC() Engine_Ext_Mur_ABC::~Engine_Ext_Mur_ABC()
@ -47,8 +65,8 @@ Engine_Ext_Mur_ABC::~Engine_Ext_Mur_ABC()
void Engine_Ext_Mur_ABC::DoPreVoltageUpdates() void Engine_Ext_Mur_ABC::DoPreVoltageUpdates()
{ {
if (IsActive()==false) return;
if (m_Eng==NULL) return; if (m_Eng==NULL) return;
if (m_Mur_Coeff==0) return;
unsigned int pos[] = {0,0,0}; unsigned int pos[] = {0,0,0};
unsigned int pos_shift[] = {0,0,0}; unsigned int pos_shift[] = {0,0,0};
pos[m_ny] = m_LineNr; pos[m_ny] = m_LineNr;
@ -60,16 +78,16 @@ void Engine_Ext_Mur_ABC::DoPreVoltageUpdates()
for (pos[m_nyPP]=0;pos[m_nyPP]<m_numLines[1];++pos[m_nyPP]) for (pos[m_nyPP]=0;pos[m_nyPP]<m_numLines[1];++pos[m_nyPP])
{ {
pos_shift[m_nyPP] = pos[m_nyPP]; pos_shift[m_nyPP] = pos[m_nyPP];
m_volt_nyP[pos[m_nyP]][pos[m_nyPP]] = m_Eng->GetVolt(m_nyP,pos_shift) - m_Mur_Coeff * m_Eng->GetVolt(m_nyP,pos); m_volt_nyP[pos[m_nyP]][pos[m_nyPP]] = m_Eng->GetVolt(m_nyP,pos_shift) - m_Op_mur->m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->GetVolt(m_nyP,pos);
m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]] = m_Eng->GetVolt(m_nyPP,pos_shift) - m_Mur_Coeff * m_Eng->GetVolt(m_nyPP,pos); m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]] = m_Eng->GetVolt(m_nyPP,pos_shift) - m_Op_mur->m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->GetVolt(m_nyPP,pos);
} }
} }
} }
void Engine_Ext_Mur_ABC::DoPostVoltageUpdates() void Engine_Ext_Mur_ABC::DoPostVoltageUpdates()
{ {
if (IsActive()==false) return;
if (m_Eng==NULL) return; if (m_Eng==NULL) return;
if (m_Mur_Coeff==0) return;
unsigned int pos[] = {0,0,0}; unsigned int pos[] = {0,0,0};
unsigned int pos_shift[] = {0,0,0}; unsigned int pos_shift[] = {0,0,0};
pos[m_ny] = m_LineNr; pos[m_ny] = m_LineNr;
@ -81,16 +99,16 @@ void Engine_Ext_Mur_ABC::DoPostVoltageUpdates()
for (pos[m_nyPP]=0;pos[m_nyPP]<m_numLines[1];++pos[m_nyPP]) for (pos[m_nyPP]=0;pos[m_nyPP]<m_numLines[1];++pos[m_nyPP])
{ {
pos_shift[m_nyPP] = pos[m_nyPP]; pos_shift[m_nyPP] = pos[m_nyPP];
m_volt_nyP[pos[m_nyP]][pos[m_nyPP]] += m_Mur_Coeff * m_Eng->GetVolt(m_nyP,pos_shift); m_volt_nyP[pos[m_nyP]][pos[m_nyPP]] += m_Op_mur->m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->GetVolt(m_nyP,pos_shift);
m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]] += m_Mur_Coeff * m_Eng->GetVolt(m_nyPP,pos_shift); m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]] += m_Op_mur->m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->GetVolt(m_nyPP,pos_shift);
} }
} }
} }
void Engine_Ext_Mur_ABC::Apply2Voltages() void Engine_Ext_Mur_ABC::Apply2Voltages()
{ {
if (IsActive()==false) return;
if (m_Eng==NULL) return; if (m_Eng==NULL) return;
if (m_Mur_Coeff==0) return;
unsigned int pos[] = {0,0,0}; unsigned int pos[] = {0,0,0};
pos[m_ny] = m_LineNr; pos[m_ny] = m_LineNr;

7
FDTD/engine_ext_mur_abc.h
View File

@ -19,6 +19,7 @@
#define ENGINE_EXT_MUR_ABC_H #define ENGINE_EXT_MUR_ABC_H
#include "engine_extension.h" #include "engine_extension.h"
#include "engine.h"
#include "operator.h" #include "operator.h"
class Operator_Ext_Mur_ABC; class Operator_Ext_Mur_ABC;
@ -36,13 +37,17 @@ public:
protected: protected:
Operator_Ext_Mur_ABC* m_Op_mur; Operator_Ext_Mur_ABC* m_Op_mur;
inline bool IsActive() {if (m_Eng->GetNumberOfTimesteps()<m_start_TS) return false; return true;}
unsigned int m_start_TS;
int m_ny; int m_ny;
int m_nyP,m_nyPP; int m_nyP,m_nyPP;
unsigned int m_LineNr; unsigned int m_LineNr;
int m_LineNr_Shift; int m_LineNr_Shift;
unsigned int m_numLines[2]; unsigned int m_numLines[2];
FDTD_FLOAT m_Mur_Coeff; FDTD_FLOAT** m_Mur_Coeff_nyP;
FDTD_FLOAT** m_Mur_Coeff_nyPP;
FDTD_FLOAT** m_volt_nyP; //n+1 direction FDTD_FLOAT** m_volt_nyP; //n+1 direction
FDTD_FLOAT** m_volt_nyPP; //n+2 direction FDTD_FLOAT** m_volt_nyPP; //n+2 direction
}; };

249
FDTD/operator.cpp
View File

@ -45,18 +45,23 @@ void Operator::Init()
{ {
CSX = NULL; CSX = NULL;
ExciteSignal = NULL; ExciteSignal_volt = NULL;
ExciteSignal_curr = NULL;
E_Exc_delay = NULL; E_Exc_delay = NULL;
E_Exc_amp=NULL; E_Exc_amp=NULL;
E_Exc_dir=NULL; E_Exc_dir=NULL;
vv=NULL; vv=NULL;
vi=NULL; vi=NULL;
Curr_Exc_delay = NULL;
Curr_Exc_amp=NULL;
Curr_Exc_dir=NULL;
iv=NULL; iv=NULL;
ii=NULL; ii=NULL;
for (int n=0;n<3;++n) for (int n=0;n<3;++n)
{ {
discLines[n]=NULL; discLines[n]=NULL;
E_Exc_index[n]=NULL; E_Exc_index[n]=NULL;
Curr_Exc_index[n]=NULL;
} }
MainOp=NULL; MainOp=NULL;
@ -69,14 +74,21 @@ void Operator::Init()
EC_L[n]=NULL; EC_L[n]=NULL;
EC_R[n]=NULL; EC_R[n]=NULL;
} }
for (int n=0;n<6;++n)
m_BC[n]=0;
} }
void Operator::Reset() void Operator::Reset()
{ {
delete[] ExciteSignal; delete[] ExciteSignal_volt;
delete[] ExciteSignal_curr;
delete[] E_Exc_delay; delete[] E_Exc_delay;
delete[] E_Exc_dir; delete[] E_Exc_dir;
delete[] E_Exc_amp; delete[] E_Exc_amp;
delete[] Curr_Exc_delay;
delete[] Curr_Exc_dir;
delete[] Curr_Exc_amp;
Delete_N_3DArray(vv,numLines); Delete_N_3DArray(vv,numLines);
Delete_N_3DArray(vi,numLines); Delete_N_3DArray(vi,numLines);
Delete_N_3DArray(iv,numLines); Delete_N_3DArray(iv,numLines);
@ -85,6 +97,7 @@ void Operator::Reset()
{ {
delete[] discLines[n]; delete[] discLines[n];
delete[] E_Exc_index[n]; delete[] E_Exc_index[n];
delete[] Curr_Exc_index[n];
} }
delete MainOp; delete MainOp;
delete DualOp; delete DualOp;
@ -310,9 +323,9 @@ unsigned int Operator::GetMaxExcitationTimestep() const
unsigned int maxStep=0; unsigned int maxStep=0;
for (unsigned int n=1;n<ExciteLength+1;++n) for (unsigned int n=1;n<ExciteLength+1;++n)
{ {
if (fabs(ExciteSignal[n])>maxAmp) if (fabs(ExciteSignal_volt[n])>maxAmp)
{ {
maxAmp = fabs(ExciteSignal[n]); maxAmp = fabs(ExciteSignal_volt[n]);
maxStep = n; maxStep = n;
} }
} }
@ -325,13 +338,18 @@ unsigned int Operator::CalcGaussianPulsExcitation(double f0, double fc)
ExciteLength = (unsigned int)(2.0 * 9.0/(2.0*PI*fc) / dT); ExciteLength = (unsigned int)(2.0 * 9.0/(2.0*PI*fc) / dT);
// cerr << "Operator::CalcGaussianPulsExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl; // cerr << "Operator::CalcGaussianPulsExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl;
delete[] ExciteSignal; delete[] ExciteSignal_volt;
ExciteSignal = new FDTD_FLOAT[ExciteLength+1]; delete[] ExciteSignal_curr;
ExciteSignal[0]=0.0; ExciteSignal_volt = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal_curr = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal_volt[0]=0.0;
ExciteSignal_curr[0]=0.0;
for (unsigned int n=1;n<ExciteLength+1;++n) for (unsigned int n=1;n<ExciteLength+1;++n)
{ {
ExciteSignal[n] = cos(2.0*PI*f0*(n*dT-9.0/(2.0*PI*fc)))*exp(-1*pow(2.0*PI*fc*n*dT/3.0-3,2)); double t = (n-1)*dT;
// cerr << ExciteSignal[n] << endl; ExciteSignal_volt[n] = cos(2.0*PI*f0*(t-9.0/(2.0*PI*fc)))*exp(-1*pow(2.0*PI*fc*t/3.0-3,2));
t += 0.5*dT;
ExciteSignal_curr[n] = cos(2.0*PI*f0*(t-9.0/(2.0*PI*fc)))*exp(-1*pow(2.0*PI*fc*t/3.0-3,2));
} }
return CalcNyquistNum(f0+fc); return CalcNyquistNum(f0+fc);
} }
@ -342,10 +360,14 @@ unsigned int Operator::CalcDiracPulsExcitation()
ExciteLength = 1; ExciteLength = 1;
// cerr << "Operator::CalcDiracPulsExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl; // cerr << "Operator::CalcDiracPulsExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl;
delete[] ExciteSignal; delete[] ExciteSignal_volt;
ExciteSignal = new FDTD_FLOAT[ExciteLength+1]; delete[] ExciteSignal_curr;
ExciteSignal[0]=0.0; ExciteSignal_volt = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal[1]=1.0; ExciteSignal_curr = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal_volt[0]=0.0;
ExciteSignal_volt[1]=1.0;
ExciteSignal_curr[0]=0.0;
ExciteSignal_curr[1]=1.0;
return 1; return 1;
} }
@ -355,10 +377,14 @@ unsigned int Operator::CalcStepExcitation()
if (dT==0) return 0; if (dT==0) return 0;
ExciteLength = 1; ExciteLength = 1;
delete[] ExciteSignal; delete[] ExciteSignal_volt;
ExciteSignal = new FDTD_FLOAT[ExciteLength+1]; delete[] ExciteSignal_curr;
ExciteSignal[0]=1.0; ExciteSignal_volt = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal[1]=1.0; ExciteSignal_curr = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal_volt[0]=1.0;
ExciteSignal_volt[1]=1.0;
ExciteSignal_curr[0]=0.0;
ExciteSignal_curr[1]=1.0;
return 1; return 1;
} }
@ -370,9 +396,12 @@ unsigned int Operator::CalcCustomExcitation(double f0, int nTS, string signal)
ExciteLength = (unsigned int)(nTS); ExciteLength = (unsigned int)(nTS);
// cerr << "Operator::CalcSinusExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl; // cerr << "Operator::CalcSinusExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl;
delete[] ExciteSignal; delete[] ExciteSignal_volt;
ExciteSignal = new FDTD_FLOAT[ExciteLength+1]; delete[] ExciteSignal_curr;
ExciteSignal[0]=0.0; ExciteSignal_volt = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal_curr = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal_volt[0]=0.0;
ExciteSignal_curr[0]=0.0;
FunctionParser fParse; FunctionParser fParse;
fParse.AddConstant("pi", 3.14159265358979323846); fParse.AddConstant("pi", 3.14159265358979323846);
fParse.AddConstant("e", 2.71828182845904523536); fParse.AddConstant("e", 2.71828182845904523536);
@ -386,8 +415,9 @@ unsigned int Operator::CalcCustomExcitation(double f0, int nTS, string signal)
for (unsigned int n=1;n<ExciteLength+1;++n) for (unsigned int n=1;n<ExciteLength+1;++n)
{ {
vars[0] = (n-1)*GetTimestep(); vars[0] = (n-1)*GetTimestep();
ExciteSignal[n] = fParse.Eval(vars); ExciteSignal_volt[n] = fParse.Eval(vars);
// cerr << ExciteSignal[n] << endl; vars[0] += 0.5*GetTimestep();
ExciteSignal_curr[n] = fParse.Eval(vars);
} }
return CalcNyquistNum(f0); return CalcNyquistNum(f0);
} }
@ -399,13 +429,18 @@ unsigned int Operator::CalcSinusExcitation(double f0, int nTS)
ExciteLength = (unsigned int)(nTS); ExciteLength = (unsigned int)(nTS);
// cerr << "Operator::CalcSinusExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl; // cerr << "Operator::CalcSinusExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl;
delete[] ExciteSignal; delete[] ExciteSignal_volt;
ExciteSignal = new FDTD_FLOAT[ExciteLength+1]; delete[] ExciteSignal_curr;
ExciteSignal[0]=0.0; ExciteSignal_volt = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal_curr = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal_volt[0]=0.0;
ExciteSignal_curr[0]=0.0;
for (unsigned int n=1;n<ExciteLength+1;++n) for (unsigned int n=1;n<ExciteLength+1;++n)
{ {
ExciteSignal[n] = sin(2.0*PI*f0*n*dT); double t = (n-1)*dT;
// cerr << ExciteSignal[n] << endl; ExciteSignal_volt[n] = sin(2.0*PI*f0*t);
t += 0.5*dT;
ExciteSignal_curr[n] = sin(2.0*PI*f0*t);
} }
return CalcNyquistNum(f0); return CalcNyquistNum(f0);
} }
@ -574,9 +609,14 @@ int Operator::CalcECOperator()
bool PEC[6]={1,1,1,1,1,1}; bool PEC[6]={1,1,1,1,1,1};
ApplyElectricBC(PEC); ApplyElectricBC(PEC);
if (CalcEFieldExcitation()==false) return -1; if (CalcFieldExcitation()==false) return -1;
CalcPEC(); CalcPEC();
bool PMC[6];
for (int n=0;n<6;++n)
PMC[n] = m_BC[n]==1;
ApplyMagneticBC(PMC);
return 0; return 0;
} }
@ -883,18 +923,25 @@ double Operator::CalcTimestep()
return 0; return 0;
} }
bool Operator::CalcEFieldExcitation() bool Operator::CalcFieldExcitation()
{ {
if (dT==0) return false; if (dT==0) return false;
vector<unsigned int> vIndex[3];
vector<FDTD_FLOAT> vExcit;
vector<unsigned int> vDelay;
vector<unsigned int> vDir;
unsigned int pos[3]; unsigned int pos[3];
double coord[3];
double delta[3]; double delta[3];
double amp=0; double amp=0;
vector<unsigned int> volt_vIndex[3];
vector<FDTD_FLOAT> volt_vExcit;
vector<unsigned int> volt_vDelay;
vector<unsigned int> volt_vDir;
double volt_coord[3];
vector<unsigned int> curr_vIndex[3];
vector<FDTD_FLOAT> curr_vExcit;
vector<unsigned int> curr_vDelay;
vector<unsigned int> curr_vDir;
double curr_coord[3];
for (pos[2]=0;pos[2]<numLines[2];++pos[2]) for (pos[2]=0;pos[2]<numLines[2];++pos[2])
{ {
delta[2]=fabs(MainOp->GetIndexDelta(2,pos[2])); delta[2]=fabs(MainOp->GetIndexDelta(2,pos[2]));
@ -904,29 +951,31 @@ bool Operator::CalcEFieldExcitation()
for (pos[0]=0;pos[0]<numLines[0];++pos[0]) for (pos[0]=0;pos[0]<numLines[0];++pos[0])
{ {
delta[0]=fabs(MainOp->GetIndexDelta(0,pos[0])); delta[0]=fabs(MainOp->GetIndexDelta(0,pos[0]));
//electric field excite
for (int n=0;n<3;++n) for (int n=0;n<3;++n)
{ {
coord[0] = discLines[0][pos[0]]; volt_coord[0] = discLines[0][pos[0]];
coord[1] = discLines[1][pos[1]]; volt_coord[1] = discLines[1][pos[1]];
coord[2] = discLines[2][pos[2]]; volt_coord[2] = discLines[2][pos[2]];
coord[n]+=delta[n]*0.5; volt_coord[n]+=delta[n]*0.5;
CSProperties* prop = CSX->GetPropertyByCoordPriority(coord,(CSProperties::PropertyType)(CSProperties::ELECTRODE)); CSProperties* prop = CSX->GetPropertyByCoordPriority(volt_coord,(CSProperties::PropertyType)(CSProperties::ELECTRODE));
if (prop) if (prop)
{ {
CSPropElectrode* elec = prop->ToElectrode(); CSPropElectrode* elec = prop->ToElectrode();
if (elec!=NULL) if (elec!=NULL)
{ {
if ((elec->GetActiveDir(n)) )//&& (pos[n]<numLines[n]-1)) if ((elec->GetActiveDir(n)) && ( (elec->GetExcitType()==0) || (elec->GetExcitType()==1) ))//&& (pos[n]<numLines[n]-1))
{ {
amp = elec->GetWeightedExcitation(n,coord)*GetMeshDelta(n,pos);// delta[n]*gridDelta; amp = elec->GetWeightedExcitation(n,volt_coord)*GetMeshDelta(n,pos);// delta[n]*gridDelta;
if (amp!=0) if (amp!=0)
{ {
vExcit.push_back(amp); volt_vExcit.push_back(amp);
vDelay.push_back((unsigned int)(elec->GetDelay()/dT)); volt_vDelay.push_back((unsigned int)(elec->GetDelay()/dT));
vDir.push_back(n); volt_vDir.push_back(n);
vIndex[0].push_back(pos[0]); volt_vIndex[0].push_back(pos[0]);
vIndex[1].push_back(pos[1]); volt_vIndex[1].push_back(pos[1]);
vIndex[2].push_back(pos[2]); volt_vIndex[2].push_back(pos[2]);
} }
if (elec->GetExcitType()==1) //hard excite if (elec->GetExcitType()==1) //hard excite
{ {
@ -937,6 +986,45 @@ bool Operator::CalcEFieldExcitation()
} }
} }
} }
//magnetic field excite
for (int n=0;n<3;++n)
{
int nP = (n+1)%3;
int nPP = (n+2)%3;
curr_coord[0] = discLines[0][pos[0]];
curr_coord[1] = discLines[1][pos[1]];
curr_coord[2] = discLines[2][pos[2]];
curr_coord[nP] +=delta[nP]*0.5;
curr_coord[nPP] +=delta[nPP]*0.5;
CSProperties* prop = CSX->GetPropertyByCoordPriority(curr_coord,(CSProperties::PropertyType)(CSProperties::ELECTRODE));
if (prop)
{
CSPropElectrode* elec = prop->ToElectrode();
if (elec!=NULL)
{
if ((elec->GetActiveDir(n)) && ( (elec->GetExcitType()==2) || (elec->GetExcitType()==3) ))//&& (pos[n]<numLines[n]-1))
{
amp = elec->GetWeightedExcitation(n,curr_coord)*GetMeshDelta(n,pos,true);// delta[n]*gridDelta;
if (amp!=0)
{
curr_vExcit.push_back(amp);
curr_vDelay.push_back((unsigned int)(elec->GetDelay()/dT));
curr_vDir.push_back(n);
curr_vIndex[0].push_back(pos[0]);
curr_vIndex[1].push_back(pos[1]);
curr_vIndex[2].push_back(pos[2]);
}
if (elec->GetExcitType()==3) //hard excite
{
ii[n][pos[0]][pos[1]][pos[2]] = 0;
iv[n][pos[0]][pos[1]][pos[2]] = 0;
}
}
}
}
}
} }
} }
} }
@ -972,24 +1060,24 @@ bool Operator::CalcEFieldExcitation()
pos[2] = path.posPath[2].at(t); pos[2] = path.posPath[2].at(t);
MainOp->SetPos(pos[0],pos[1],pos[2]); MainOp->SetPos(pos[0],pos[1],pos[2]);
deltaN=fabs(MainOp->GetIndexDelta(n,pos[n])); deltaN=fabs(MainOp->GetIndexDelta(n,pos[n]));
coord[0] = discLines[0][pos[0]]; volt_coord[0] = discLines[0][pos[0]];
coord[1] = discLines[1][pos[1]]; volt_coord[1] = discLines[1][pos[1]];
coord[2] = discLines[2][pos[2]]; volt_coord[2] = discLines[2][pos[2]];
coord[n] += 0.5*deltaN; volt_coord[n] += 0.5*deltaN;
// cerr << n << " " << coord[0] << " " << coord[1] << " " << coord[2] << endl; // cerr << n << " " << coord[0] << " " << coord[1] << " " << coord[2] << endl;
if (elec!=NULL) if (elec!=NULL)
{ {
if ((elec->GetActiveDir(n)) && (pos[n]<numLines[n]-1)) if ((elec->GetActiveDir(n)) && (pos[n]<numLines[n]-1) && ( (elec->GetExcitType()==0) || (elec->GetExcitType()==1) ))
{ {
amp = elec->GetWeightedExcitation(n,coord)*deltaN*gridDelta; amp = elec->GetWeightedExcitation(n,volt_coord)*deltaN*gridDelta;
if (amp!=0) if (amp!=0)
{ {
vExcit.push_back(amp); volt_vExcit.push_back(amp);
vDelay.push_back((unsigned int)(elec->GetDelay()/dT)); volt_vDelay.push_back((unsigned int)(elec->GetDelay()/dT));
vDir.push_back(n); volt_vDir.push_back(n);
vIndex[0].push_back(pos[0]); volt_vIndex[0].push_back(pos[0]);
vIndex[1].push_back(pos[1]); volt_vIndex[1].push_back(pos[1]);
vIndex[2].push_back(pos[2]); volt_vIndex[2].push_back(pos[2]);
} }
if (elec->GetExcitType()==1) //hard excite if (elec->GetExcitType()==1) //hard excite
{ {
@ -1004,16 +1092,18 @@ bool Operator::CalcEFieldExcitation()
} }
} }
E_Exc_Count = vExcit.size();
cerr << "Operator::CalcEFieldExcitation: Found number of excitations points: " << E_Exc_Count << endl; // set voltage excitations
E_Exc_Count = volt_vExcit.size();
cerr << "Operator::CalcFieldExcitation: Number of voltage excitation points: " << E_Exc_Count << endl;
if (E_Exc_Count==0) if (E_Exc_Count==0)
cerr << "No E-Field excitation found!" << endl; cerr << "No E-Field/voltage excitation found!" << endl;
for (int n=0;n<3;++n) for (int n=0;n<3;++n)
{ {
delete[] E_Exc_index[n]; delete[] E_Exc_index[n];
E_Exc_index[n] = new unsigned int[E_Exc_Count]; E_Exc_index[n] = new unsigned int[E_Exc_Count];
for (unsigned int i=0;i<E_Exc_Count;++i) for (unsigned int i=0;i<E_Exc_Count;++i)
E_Exc_index[n][i]=vIndex[n].at(i); E_Exc_index[n][i]=volt_vIndex[n].at(i);
} }
delete[] E_Exc_delay; delete[] E_Exc_delay;
E_Exc_delay = new unsigned int[E_Exc_Count]; E_Exc_delay = new unsigned int[E_Exc_Count];
@ -1023,10 +1113,37 @@ bool Operator::CalcEFieldExcitation()
E_Exc_dir = new unsigned short[E_Exc_Count]; E_Exc_dir = new unsigned short[E_Exc_Count];
for (unsigned int i=0;i<E_Exc_Count;++i) for (unsigned int i=0;i<E_Exc_Count;++i)
{ {
E_Exc_delay[i]=vDelay.at(i); E_Exc_delay[i]=volt_vDelay.at(i);
E_Exc_amp[i]=vExcit.at(i); E_Exc_amp[i]=volt_vExcit.at(i);
E_Exc_dir[i]=vDir.at(i); E_Exc_dir[i]=volt_vDir.at(i);
} }
// set current excitations
Curr_Exc_Count = curr_vExcit.size();
cerr << "Operator::CalcFieldExcitation: Number of current excitation points: " << Curr_Exc_Count << endl;
if (Curr_Exc_Count==0)
cerr << "No H-Field/current excitation found!" << endl;
for (int n=0;n<3;++n)
{
delete[] Curr_Exc_index[n];
Curr_Exc_index[n] = new unsigned int[Curr_Exc_Count];
for (unsigned int i=0;i<Curr_Exc_Count;++i)
Curr_Exc_index[n][i]=curr_vIndex[n].at(i);
}
delete[] Curr_Exc_delay;
Curr_Exc_delay = new unsigned int[Curr_Exc_Count];
delete[] Curr_Exc_amp;
Curr_Exc_amp = new FDTD_FLOAT[Curr_Exc_Count];
delete[] Curr_Exc_dir;
Curr_Exc_dir = new unsigned short[Curr_Exc_Count];
for (unsigned int i=0;i<Curr_Exc_Count;++i)
{
Curr_Exc_delay[i]=curr_vDelay.at(i);
Curr_Exc_amp[i]=curr_vExcit.at(i);
Curr_Exc_dir[i]=curr_vDir.at(i);
}
return true; return true;
} }

21
FDTD/operator.h
View File

@ -35,6 +35,7 @@ public:
virtual ~Operator(); virtual ~Operator();
virtual bool SetGeometryCSX(ContinuousStructure* geo); virtual bool SetGeometryCSX(ContinuousStructure* geo);
virtual ContinuousStructure* GetGeometryCSX() {return CSX;}
virtual int CalcECOperator(); virtual int CalcECOperator();
@ -52,6 +53,7 @@ public:
//! Get the excitation timestep with the (first) max amplitude //! Get the excitation timestep with the (first) max amplitude
virtual unsigned int GetMaxExcitationTimestep() const; virtual unsigned int GetMaxExcitationTimestep() const;
virtual void SetBoundaryCondition(int* BCs) {for (int n=0;n<6;++n) m_BC[n]=BCs[n];}
virtual void ApplyElectricBC(bool* dirs); //applied by default to all boundaries virtual void ApplyElectricBC(bool* dirs); //applied by default to all boundaries
virtual void ApplyMagneticBC(bool* dirs); virtual void ApplyMagneticBC(bool* dirs);
@ -103,9 +105,10 @@ protected:
ContinuousStructure* CSX; ContinuousStructure* CSX;
//E-Field Excitation int m_BC[6];
//! Calc the electric field excitation.
virtual bool CalcEFieldExcitation(); //! Calculate the field excitations.
virtual bool CalcFieldExcitation();
virtual bool CalcPEC(); virtual bool CalcPEC();
@ -144,14 +147,22 @@ public:
//Excitation time-signal //Excitation time-signal
unsigned int ExciteLength; unsigned int ExciteLength;
FDTD_FLOAT* ExciteSignal; FDTD_FLOAT* ExciteSignal_volt;
FDTD_FLOAT* ExciteSignal_curr;
//E-Field Excitation //E-Field/voltage Excitation
unsigned int E_Exc_Count; unsigned int E_Exc_Count;
unsigned int* E_Exc_index[3]; unsigned int* E_Exc_index[3];
unsigned short* E_Exc_dir; unsigned short* E_Exc_dir;
FDTD_FLOAT* E_Exc_amp; //represented as edge-voltages!! FDTD_FLOAT* E_Exc_amp; //represented as edge-voltages!!
unsigned int* E_Exc_delay; unsigned int* E_Exc_delay;
//H-Field/current Excitation
unsigned int Curr_Exc_Count;
unsigned int* Curr_Exc_index[3];
unsigned short* Curr_Exc_dir;
FDTD_FLOAT* Curr_Exc_amp; //represented as edge-currents!!
unsigned int* Curr_Exc_delay;
}; };
#endif // OPERATOR_H #endif // OPERATOR_H

69
FDTD/operator_ext_mur_abc.cpp
View File

@ -22,23 +22,35 @@
Operator_Ext_Mur_ABC::Operator_Ext_Mur_ABC(Operator* op) : Operator_Extension(op) Operator_Ext_Mur_ABC::Operator_Ext_Mur_ABC(Operator* op) : Operator_Extension(op)
{ {
m_Mur_Coeff = NULL;
m_ny = -1; m_ny = -1;
m_nyP = -1; m_nyP = -1;
m_nyPP = -1; m_nyPP = -1;
m_LineNr = 0; m_LineNr = 0;
m_LineNr_Shift = 0; m_LineNr_Shift = 0;
m_Mur_Coeff = 0;
m_Mur_Coeff_nyP = NULL;
m_Mur_Coeff_nyPP = NULL;
m_numLines[0]=0;
m_numLines[1]=0;
} }
Operator_Ext_Mur_ABC::~Operator_Ext_Mur_ABC() Operator_Ext_Mur_ABC::~Operator_Ext_Mur_ABC()
{ {
Delete2DArray(m_Mur_Coeff_nyP,m_numLines);
m_Mur_Coeff_nyP = NULL;
Delete2DArray(m_Mur_Coeff_nyPP,m_numLines);
m_Mur_Coeff_nyPP = NULL;
} }
void Operator_Ext_Mur_ABC::SetDirection(int ny, bool top_ny) void Operator_Ext_Mur_ABC::SetDirection(int ny, bool top_ny)
{ {
if ((ny<0) || (ny>2)) if ((ny<0) || (ny>2))
return; return;
Delete2DArray(m_Mur_Coeff_nyP,m_numLines);
Delete2DArray(m_Mur_Coeff_nyPP,m_numLines);
m_ny = ny; m_ny = ny;
m_nyP = (ny+1)%3; m_nyP = (ny+1)%3;
m_nyPP = (ny+2)%3; m_nyPP = (ny+2)%3;
@ -55,6 +67,10 @@ void Operator_Ext_Mur_ABC::SetDirection(int ny, bool top_ny)
m_numLines[0] = m_Op->GetNumberOfLines(m_nyP); m_numLines[0] = m_Op->GetNumberOfLines(m_nyP);
m_numLines[1] = m_Op->GetNumberOfLines(m_nyPP); m_numLines[1] = m_Op->GetNumberOfLines(m_nyPP);
m_Mur_Coeff_nyP = Create2DArray(m_numLines);
m_Mur_Coeff_nyPP = Create2DArray(m_numLines);
} }
bool Operator_Ext_Mur_ABC::BuildExtension() bool Operator_Ext_Mur_ABC::BuildExtension()
@ -68,8 +84,53 @@ bool Operator_Ext_Mur_ABC::BuildExtension()
unsigned int pos[] = {0,0,0}; unsigned int pos[] = {0,0,0};
pos[m_ny] = m_LineNr; pos[m_ny] = m_LineNr;
double delta = fabs(m_Op->GetMeshDelta(m_ny,pos)); double delta = fabs(m_Op->GetMeshDelta(m_ny,pos));
m_Mur_Coeff = (__C0__ * dT - delta) / (__C0__ * dT + delta); double coord[] = {0,0,0};
// cerr << "Operator_Ext_Mur_ABC::BuildExtension(): " << m_Mur_Coeff << endl; coord[0] = m_Op->GetDiscLine(0,pos[0]);
coord[1] = m_Op->GetDiscLine(1,pos[1]);
coord[2] = m_Op->GetDiscLine(2,pos[2]);
double eps,mue;
double c0t;
if (m_LineNr==0)
coord[m_ny] = m_Op->GetDiscLine(m_ny,pos[m_ny]) + delta/2 / m_Op->GetGridDelta();
else
coord[m_ny] = m_Op->GetDiscLine(m_ny,pos[m_ny]) - delta/2 / m_Op->GetGridDelta();
for (pos[m_nyP]=0;pos[m_nyP]<m_numLines[0];++pos[m_nyP])
{
coord[m_nyP] = m_Op->GetDiscLine(m_nyP,pos[m_nyP]);
for (pos[m_nyPP]=0;pos[m_nyPP]<m_numLines[1];++pos[m_nyPP])
{
coord[m_nyPP] = m_Op->GetDiscLine(m_nyPP,pos[m_nyPP]);
CSProperties* prop = m_Op->GetGeometryCSX()->GetPropertyByCoordPriority(coord,CSProperties::MATERIAL);
if (prop)
{
CSPropMaterial* mat = prop->ToMaterial();
//nP
eps = mat->GetEpsilonWeighted(m_nyP,coord);
mue = mat->GetMueWeighted(m_nyP,coord);
c0t = __C0__ * dT / sqrt(eps*mue);
m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]] = (c0t - delta) / (c0t + delta);
//nPP
eps = mat->GetEpsilonWeighted(m_nyPP,coord);
mue = mat->GetMueWeighted(m_nyPP,coord);
c0t = __C0__ * dT / sqrt(eps*mue);
m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] = (c0t - delta) / (c0t + delta);
}
else
{
c0t = __C0__ * dT;
m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]] = (c0t - delta) / (c0t + delta);
m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] = m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]];
}
// cerr << m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]] << " : " << m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]] << endl;
}
}
// cerr << "Operator_Ext_Mur_ABC::BuildExtension(): " << m_ny << " @ " << m_LineNr << endl;
return true; return true;
} }

3
FDTD/operator_ext_mur_abc.h
View File

@ -43,7 +43,8 @@ protected:
unsigned int m_numLines[2]; unsigned int m_numLines[2];
FDTD_FLOAT m_Mur_Coeff; FDTD_FLOAT** m_Mur_Coeff_nyP;
FDTD_FLOAT** m_Mur_Coeff_nyPP;
}; };
#endif // OPERATOR_EXT_MUR_ABC_H #endif // OPERATOR_EXT_MUR_ABC_H

1
FDTD/operator_extension.h
View File

@ -24,6 +24,7 @@ class Engine_Extension;
//! Abstract base-class for all operator extensions //! Abstract base-class for all operator extensions
class Operator_Extension class Operator_Extension
{ {
friend class Engine_Extension;
public: public:
virtual bool BuildExtension() {return true;} virtual bool BuildExtension() {return true;}

7
FDTD/processfields.cpp
View File

@ -307,7 +307,7 @@ bool ProcessFields::DumpMultiScalarArray2VTK(ofstream &file, string names[], FDT
return true; return true;
} }
bool ProcessFields::DumpVectorArray2HDF5(string filename, string name, FDTD_FLOAT const* const* const* const* array, unsigned int const* numLines) bool ProcessFields::DumpVectorArray2HDF5(string filename, string name, FDTD_FLOAT const* const* const* const* array, unsigned int const* numLines, float time)
{ {
const H5std_string FILE_NAME(filename); const H5std_string FILE_NAME(filename);
const H5std_string DATASET_NAME( name ); const H5std_string DATASET_NAME( name );
@ -329,6 +329,11 @@ bool ProcessFields::DumpVectorArray2HDF5(string filename, string name, FDTD_FLOA
// datatype.setOrder( H5T_ORDER_LE ); // datatype.setOrder( H5T_ORDER_LE );
H5::DataSet dataset = group.createDataSet( DATASET_NAME, datatype, dataspace ); H5::DataSet dataset = group.createDataSet( DATASET_NAME, datatype, dataspace );
hsize_t t_dimsf[] = {1};
H5::DataSpace t_dataspace( 1, t_dimsf );
H5::Attribute attr = dataset.createAttribute("time",H5::PredType::NATIVE_FLOAT,t_dataspace);
attr.write( H5::PredType::NATIVE_FLOAT , &time);
// I have not the slightest idea why this array-copy action is necessary... but it's the only way hdf5 does what it is supposed to do anyway!! // I have not the slightest idea why this array-copy action is necessary... but it's the only way hdf5 does what it is supposed to do anyway!!
// at least it is save in case FDTD_FLOAT was defined as double... // at least it is save in case FDTD_FLOAT was defined as double...
// why does hdf5 write the dimensions backwards??? or matlab??? // why does hdf5 write the dimensions backwards??? or matlab???

2
FDTD/processfields.h
View File

@ -57,7 +57,7 @@ public:
static bool DumpScalarArray2VTK(ofstream &file, string name, FDTD_FLOAT const* const* const* array, double const* const* discLines, unsigned int const* numLines, unsigned int precision=12); static bool DumpScalarArray2VTK(ofstream &file, string name, FDTD_FLOAT const* const* const* array, double const* const* discLines, unsigned int const* numLines, unsigned int precision=12);
static bool DumpMultiScalarArray2VTK(ofstream &file, string names[], FDTD_FLOAT const* const* const* const* array, unsigned int numFields, double const* const* discLines, unsigned int const* numLines, unsigned int precision=12); static bool DumpMultiScalarArray2VTK(ofstream &file, string names[], FDTD_FLOAT const* const* const* const* array, unsigned int numFields, double const* const* discLines, unsigned int const* numLines, unsigned int precision=12);
static bool DumpVectorArray2HDF5(string filename, string name, FDTD_FLOAT const* const* const* const* array, unsigned int const* numLines); static bool DumpVectorArray2HDF5(string filename, string name, FDTD_FLOAT const* const* const* const* array, unsigned int const* numLines, float time=0);
double CalcTotalEnergy() const; double CalcTotalEnergy() const;

12
FDTD/processfields_td.cpp
View File

@ -44,10 +44,10 @@ void ProcessFieldsTD::DumpCellInterpol(string filename)
OpPos[0]=start[0]+pos[0]*subSample[0]; OpPos[0]=start[0]+pos[0]*subSample[0];
for (pos[1]=0;pos[1]<numDLines[1];++pos[1]) for (pos[1]=0;pos[1]<numDLines[1];++pos[1])
{ {
OpPos[1]=start[1]+pos[1]*subSample[0]; OpPos[1]=start[1]+pos[1]*subSample[1];
for (pos[2]=0;pos[2]<numDLines[2];++pos[2]) for (pos[2]=0;pos[2]<numDLines[2];++pos[2])
{ {
OpPos[2]=start[2]+pos[2]*subSample[0]; OpPos[2]=start[2]+pos[2]*subSample[2];
//in x //in x
delta = Op->GetMeshDelta(0,OpPos); //Op->discLines[0][OpPos[0]+1] - Op->discLines[0][OpPos[0]]; delta = Op->GetMeshDelta(0,OpPos); //Op->discLines[0][OpPos[0]+1] - Op->discLines[0][OpPos[0]];
if (delta) if (delta)
@ -84,7 +84,7 @@ void ProcessFieldsTD::DumpCellInterpol(string filename)
{ {
stringstream ss; stringstream ss;
ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps();
DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),E_T,numDLines); DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),E_T,numDLines,Eng->GetNumberOfTimesteps()*Op->GetTimestep());
} }
else else
cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl; cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl;
@ -148,7 +148,7 @@ void ProcessFieldsTD::DumpCellInterpol(string filename)
{ {
stringstream ss; stringstream ss;
ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps();
DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),H_T,numDLines); DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),H_T,numDLines,(0.5+Eng->GetNumberOfTimesteps())*Op->GetTimestep());
} }
else else
cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl; cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl;
@ -198,7 +198,7 @@ void ProcessFieldsTD::DumpNoInterpol(string filename)
{ {
stringstream ss; stringstream ss;
ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps();
DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),E_T,numLines); DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),E_T,numLines,Eng->GetNumberOfTimesteps()*Op->GetTimestep());
} }
else else
cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl; cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl;
@ -244,7 +244,7 @@ void ProcessFieldsTD::DumpNoInterpol(string filename)
{ {
stringstream ss; stringstream ss;
ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps();
DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),H_T,numLines); DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),H_T,numLines,(0.5+Eng->GetNumberOfTimesteps())*Op->GetTimestep());
} }
else else
cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl; cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl;

17
TESTSUITE/helperscripts/invoke_openEMS.m
View File

@ -1,17 +0,0 @@
function invoke_openEMS( opts )
if nargin < 1
opts = '';
end
% opts = [opts ' --disable-dumps'];
% opts = [opts ' --debug-material'];
% opts = [opts ' --engine=multithreaded'];
% opts = [opts ' --engine=sse'];
filename = mfilename('fullpath');
dir = fileparts( filename );
openEMS_Path = [dir '/../../'];
command = [openEMS_Path 'openEMS.sh ' opts];
disp(command);
system(command);

3
matlab/ReadHDF5FieldData.m
View File

@ -15,11 +15,12 @@ for n=1:numel(info.GroupHierarchy.Groups)
if strcmp(info.GroupHierarchy.Groups(n).Name,'/FieldData') if strcmp(info.GroupHierarchy.Groups(n).Name,'/FieldData')
for m=1:numel(info.GroupHierarchy.Groups(n).Datasets) for m=1:numel(info.GroupHierarchy.Groups(n).Datasets)
names{m} = info.GroupHierarchy.Groups(n).Datasets(m).Name; names{m} = info.GroupHierarchy.Groups(n).Datasets(m).Name;
hdf_fielddata.time(m) = double(info.GroupHierarchy.Groups(n).Datasets(m).Attributes.Value);
end end
end end
end end
hdf_fielddata.names = names; hdf_fielddata.names = names;
for n=1:numel(names) for n=1:numel(names)
hdf_fielddata.values{n} = hdf5read(file,names{n}); hdf_fielddata.values{n} = double(hdf5read(file,names{n}));
end end

20
matlab/examples/Circ_Waveguide.m
View File

@ -1,5 +1,5 @@
close all; close all
clear all; clear
clc clc
%% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -18,7 +18,7 @@ f0 = 400e6;
p11 = 1.841; p11 = 1.841;
kc = p11 / rad /unit; kc = p11 / rad /unit;
k = 2*pi*f0/C0; k = 2*pi*f0/C0;
fc = C0*kc/2/pi fc = C0*kc/2/pi;
beta = sqrt(k^2 - kc^2); beta = sqrt(k^2 - kc^2);
kc = kc*unit; kc = kc*unit;
@ -28,8 +28,7 @@ func_Ea = [ num2str(1/kc) '*sin(a)*0.5*(j0(' num2str(kc) '*rho)-jn(2,' num2str
func_Ex = [func_Er '*cos(a) - ' func_Ea '*sin(a)']; func_Ex = [func_Er '*cos(a) - ' func_Ea '*sin(a)'];
func_Ey = [func_Er '*sin(a) + ' func_Ea '*cos(a)']; func_Ey = [func_Er '*sin(a) + ' func_Ea '*cos(a)'];
%% define file pathes and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% define and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
openEMS_Path = [pwd() '/../../']
openEMS_opts = ''; openEMS_opts = '';
% openEMS_opts = [openEMS_opts ' --disable-dumps']; % openEMS_opts = [openEMS_opts ' --disable-dumps'];
% openEMS_opts = [openEMS_opts ' --debug-material']; % openEMS_opts = [openEMS_opts ' --debug-material'];
@ -110,21 +109,14 @@ CSX = AddProbe(CSX,'ut1',0);
start = [ -rad 0 0/2 ];stop = [ rad 0 0/2 ]; start = [ -rad 0 0/2 ];stop = [ rad 0 0/2 ];
CSX = AddBox(CSX,'ut1', 0 ,start,stop); CSX = AddBox(CSX,'ut1', 0 ,start,stop);
%current calc
CSX = AddProbe(CSX,'it1',1);
mid = 0.5*(coax_rad_i+coax_rad_ai);
start = [ -mid -mid length/2 ];stop = [ mid mid length/2 ];
CSX = AddBox(CSX,'it1', 0 ,start,stop);
%% Write openEMS compatoble xml-file %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% Write openEMS compatoble xml-file %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX); WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
%% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
savePath = pwd(); savePath = pwd();
cd(Sim_Path); %cd to working dir cd(Sim_Path); %cd to working dir
command = [openEMS_Path 'openEMS.sh ' Sim_CSX ' ' openEMS_opts]; args = [Sim_CSX ' ' openEMS_opts];
disp(command); invoke_openEMS(args);
system(command)
cd(savePath); cd(savePath);
%% do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

14
matlab/examples/Circ_Waveguide_CylinderCoords.m
View File

@ -1,5 +1,5 @@
close all; close all
clear all; clear
clc clc
%% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -22,15 +22,14 @@ f0 = 400e6;
p11 = 1.841; p11 = 1.841;
kc = p11 / rad /unit; kc = p11 / rad /unit;
k = 2*pi*f0/C0; k = 2*pi*f0/C0;
fc = C0*kc/2/pi fc = C0*kc/2/pi;
beta = sqrt(k^2 - kc^2); beta = sqrt(k^2 - kc^2);
kc = kc*unit; kc = kc*unit;
func_Er = [ num2str(-1/kc^2,15) '/rho*cos(a)*j1(' num2str(kc,15) '*rho)']; func_Er = [ num2str(-1/kc^2,15) '/rho*cos(a)*j1(' num2str(kc,15) '*rho)'];
func_Ea = [ num2str(1/kc,15) '*sin(a)*0.5*(j0(' num2str(kc,15) '*rho)-jn(2,' num2str(kc,15) '*rho))']; func_Ea = [ num2str(1/kc,15) '*sin(a)*0.5*(j0(' num2str(kc,15) '*rho)-jn(2,' num2str(kc,15) '*rho))'];
%% define file pathes and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% define openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
openEMS_Path = [pwd() '/../../']
openEMS_opts = ''; openEMS_opts = '';
% openEMS_opts = [openEMS_opts ' --disable-dumps']; % openEMS_opts = [openEMS_opts ' --disable-dumps'];
% openEMS_opts = [openEMS_opts ' --debug-material']; % openEMS_opts = [openEMS_opts ' --debug-material'];
@ -117,9 +116,8 @@ WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
%% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
savePath = pwd(); savePath = pwd();
cd(Sim_Path); %cd to working dir cd(Sim_Path); %cd to working dir
command = [openEMS_Path 'openEMS.sh ' Sim_CSX ' ' openEMS_opts]; args = [Sim_CSX ' ' openEMS_opts];
disp(command); invoke_openEMS(args);
system(command)
cd(savePath); cd(savePath);
%% do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

12
matlab/examples/Coax.m
View File

@ -1,5 +1,5 @@
close all; close all
clear all; clear
clc clc
%% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -18,8 +18,7 @@ Z0 = sqrt(MUE0/EPS0);
f0 = 0.5e9; f0 = 0.5e9;
epsR = 1; epsR = 1;
%% define file pathes and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% define openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
openEMS_Path = [pwd() '/../../']
openEMS_opts = ''; openEMS_opts = '';
openEMS_opts = [openEMS_opts ' --disable-dumps']; openEMS_opts = [openEMS_opts ' --disable-dumps'];
% openEMS_opts = [openEMS_opts ' --debug-material']; % openEMS_opts = [openEMS_opts ' --debug-material'];
@ -98,9 +97,8 @@ WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
%% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
savePath = pwd(); savePath = pwd();
cd(Sim_Path); %cd to working dir cd(Sim_Path); %cd to working dir
command = [openEMS_Path 'openEMS.sh ' Sim_CSX ' ' openEMS_opts]; args = [Sim_CSX ' ' openEMS_opts];
disp(command); invoke_openEMS(args)
system(command)
cd(savePath); cd(savePath);
%% postproc & do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% postproc & do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

14
matlab/examples/Coax_CylinderCoords.m
View File

@ -1,5 +1,5 @@
close all; close all
clear all; clear
clc clc
%% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -21,8 +21,7 @@ max_alpha = max_mesh;
N_alpha = ceil(rad_a * 2*pi / max_alpha); N_alpha = ceil(rad_a * 2*pi / max_alpha);
mesh_res = [max_mesh 2*pi/N_alpha max_mesh]; mesh_res = [max_mesh 2*pi/N_alpha max_mesh];
%% define file pathes and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% define and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
openEMS_Path = [pwd() '/../../'];
openEMS_opts = ''; openEMS_opts = '';
openEMS_opts = [openEMS_opts ' --disable-dumps']; openEMS_opts = [openEMS_opts ' --disable-dumps'];
% openEMS_opts = [openEMS_opts ' --debug-material']; % openEMS_opts = [openEMS_opts ' --debug-material'];
@ -107,9 +106,8 @@ WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
%% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
savePath = pwd(); savePath = pwd();
cd(Sim_Path); %cd to working dir cd(Sim_Path); %cd to working dir
command = [openEMS_Path 'openEMS.sh ' Sim_CSX ' ' openEMS_opts]; args = [Sim_CSX ' ' openEMS_opts];
disp(command); invoke_openEMS(args);
system(command)
cd(savePath); cd(savePath);
%% postproc & do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% postproc & do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -129,5 +127,3 @@ plot(UI.FD{1}.f,real(Z),'Linewidth',2);
plot(UI.FD{1}.f,imag(Z),'r','Linewidth',2); plot(UI.FD{1}.f,imag(Z),'r','Linewidth',2);
xlim([0 2*f0]); xlim([0 2*f0]);
legend('Z_L','\Re\{Z\}','\Im\{Z\}','Location','Best'); legend('Z_L','\Re\{Z\}','\Im\{Z\}','Location','Best');

10
matlab/examples/Coax_partial_CylinderCoords.m
View File

@ -1,5 +1,5 @@
close all; close all
clear all; clear
clc clc
EPS0 = 8.85418781762e-12; EPS0 = 8.85418781762e-12;
@ -21,7 +21,6 @@ max_alpha = max_mesh;
N_alpha = ceil(rad_a * 2*pi * partial / max_alpha); N_alpha = ceil(rad_a * 2*pi * partial / max_alpha);
mesh_res = [max_mesh 2*pi*partial/N_alpha max_mesh]; mesh_res = [max_mesh 2*pi*partial/N_alpha max_mesh];
openEMS_Path = [pwd() '/../../']
openEMS_opts = ''; openEMS_opts = '';
openEMS_opts = [openEMS_opts ' --disable-dumps']; openEMS_opts = [openEMS_opts ' --disable-dumps'];
% openEMS_opts = [openEMS_opts ' --debug-material']; % openEMS_opts = [openEMS_opts ' --debug-material'];
@ -105,9 +104,8 @@ WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
%cd to working dir and run openEMS %cd to working dir and run openEMS
savePath = pwd(); savePath = pwd();
cd(Sim_Path); %cd to working dir cd(Sim_Path); %cd to working dir
command = [openEMS_Path 'openEMS.sh ' Sim_CSX ' ' openEMS_opts]; args = [Sim_CSX ' ' openEMS_opts];
disp(command); invoke_openEMS(args)
system(command)
cd(savePath); cd(savePath);
UI = ReadUI({'ut1_1','ut1_2','it1'},'tmp/'); UI = ReadUI({'ut1_1','ut1_2','it1'},'tmp/');

18
matlab/examples/Helix.m
View File

@ -1,5 +1,5 @@
close all; close all
clear; clear
clc clc
%% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -16,8 +16,7 @@ port_resist = 1000;
f_max = 100e6; f_max = 100e6;
f_excite = 300e6; f_excite = 300e6;
%% define file pathes and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% define openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
openEMS_Path = [pwd() '/../../']
openEMS_opts = ''; openEMS_opts = '';
% openEMS_opts = [openEMS_opts ' --disable-dumps']; % openEMS_opts = [openEMS_opts ' --disable-dumps'];
openEMS_opts = [openEMS_opts ' --debug-material']; openEMS_opts = [openEMS_opts ' --debug-material'];
@ -27,8 +26,8 @@ openEMS_opts = [openEMS_opts ' --debug-boxes'];
Sim_Path = 'tmp'; Sim_Path = 'tmp';
Sim_CSX = 'helix.xml'; Sim_CSX = 'helix.xml';
rmdir(Sim_Path,'s'); [status, message, messageid] = rmdir(Sim_Path,'s');
mkdir(Sim_Path); [status,message,messageid] = mkdir(Sim_Path);
%% setup FDTD parameter & excitation function %%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% setup FDTD parameter & excitation function %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FDTD = InitFDTD(30000,1e-6); FDTD = InitFDTD(30000,1e-6);
@ -120,11 +119,10 @@ CSX = AddBox(CSX,'Ht_',0 , start,stop);
WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX); WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
%% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
savePath = pwd(); savePath = pwd;
cd(Sim_Path); %cd to working dir cd(Sim_Path); %cd to working dir
command = [openEMS_Path 'openEMS.sh ' Sim_CSX ' ' openEMS_opts]; args = [Sim_CSX ' ' openEMS_opts];
disp(command); invoke_openEMS(args);
system(command)
cd(savePath); cd(savePath);
%% postproc & do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% postproc & do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

12
matlab/examples/MSL.m
View File

@ -1,5 +1,5 @@
close all; close all
clear all; clear
clc clc
%% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -14,8 +14,7 @@ mesh_res = [5 5 10];
EPS0 = 8.85418781762e-12; EPS0 = 8.85418781762e-12;
MUE0 = 1.256637062e-6; MUE0 = 1.256637062e-6;
%% define file pathes and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% define openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
openEMS_Path = [pwd() '/../../']
openEMS_opts = ''; openEMS_opts = '';
% openEMS_opts = [openEMS_opts ' --disable-dumps']; % openEMS_opts = [openEMS_opts ' --disable-dumps'];
% openEMS_opts = [openEMS_opts ' --debug-material']; % openEMS_opts = [openEMS_opts ' --debug-material'];
@ -86,8 +85,7 @@ WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
%% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
savePath = pwd(); savePath = pwd();
cd(Sim_Path); %cd to working dir cd(Sim_Path); %cd to working dir
command = [openEMS_Path 'openEMS.sh ' Sim_CSX ' ' openEMS_opts]; args = [Sim_CSX ' ' openEMS_opts];
disp(command); invoke_openEMS(args);
system(command)
cd(savePath); cd(savePath);

15
matlab/examples/PlaneWave.m
View File

@ -1,5 +1,5 @@
close all; close all
clear all; clear
clc clc
%% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -12,8 +12,7 @@ mesh_res = 25;
EPS0 = 8.85418781762e-12; EPS0 = 8.85418781762e-12;
MUE0 = 1.256637062e-6; MUE0 = 1.256637062e-6;
%% define file pathes and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% define openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
openEMS_Path = [pwd() '/../../']
openEMS_opts = ''; openEMS_opts = '';
% openEMS_opts = [openEMS_opts ' --disable-dumps']; % openEMS_opts = [openEMS_opts ' --disable-dumps'];
% openEMS_opts = [openEMS_opts ' --debug-material']; % openEMS_opts = [openEMS_opts ' --debug-material'];
@ -25,7 +24,7 @@ Sim_CSX = 'plane_wave.xml';
mkdir(Sim_Path); mkdir(Sim_Path);
%% setup FDTD parameter & excitation function %%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% setup FDTD parameter & excitation function %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FDTD = InitFDTD(5e5,1e-6,'OverSampling',10); FDTD = InitFDTD(5e5,1e-5,'OverSampling',10);
FDTD = SetGaussExcite(FDTD,0.5e9,0.5e9); FDTD = SetGaussExcite(FDTD,0.5e9,0.5e9);
BC = [1 1 0 0 0 0]; BC = [1 1 0 0 0 0];
FDTD = SetBoundaryCond(FDTD,BC); FDTD = SetBoundaryCond(FDTD,BC);
@ -71,9 +70,8 @@ WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
%% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
savePath = pwd(); savePath = pwd();
cd(Sim_Path); %cd to working dir cd(Sim_Path); %cd to working dir
command = [openEMS_Path 'openEMS.sh ' Sim_CSX ' ' openEMS_opts]; args = [Sim_CSX ' ' openEMS_opts];
disp(command); invoke_openEMS(args);
system(command)
cd(savePath); cd(savePath);
%% do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -83,4 +81,3 @@ PlotArgs.component=2;
PlotArgs.Limit = 'auto'; PlotArgs.Limit = 'auto';
PlotHDF5FieldData('tmp/Et.h5',PlotArgs) PlotHDF5FieldData('tmp/Et.h5',PlotArgs)

15
matlab/examples/Rect_Waveguide.m
View File

@ -1,5 +1,5 @@
close all; close all
clear all; clear
clc clc
%% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% setup the simulation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -24,14 +24,13 @@ f0 = 1e9;
k = 2*pi*f0/C0; k = 2*pi*f0/C0;
kc = sqrt((m*pi/a/unit)^2 + (n*pi/b/unit)^2); kc = sqrt((m*pi/a/unit)^2 + (n*pi/b/unit)^2);
fc = C0*kc/2/pi fc = C0*kc/2/pi;
beta = sqrt(k^2 - kc^2); beta = sqrt(k^2 - kc^2);
func_Ex = [num2str(n/b/unit) '*cos(' num2str(m*pi/a) '*x)*sin(' num2str(n*pi/b) '*y)']; func_Ex = [num2str(n/b/unit) '*cos(' num2str(m*pi/a) '*x)*sin(' num2str(n*pi/b) '*y)'];
func_Ey = [num2str(m/a/unit) '*sin(' num2str(m*pi/a) '*x)*cos(' num2str(n*pi/b) '*y)']; func_Ey = [num2str(m/a/unit) '*sin(' num2str(m*pi/a) '*x)*cos(' num2str(n*pi/b) '*y)'];
%% define file pathes and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% define and openEMS options %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
openEMS_Path = [pwd() '/../../']
openEMS_opts = ''; openEMS_opts = '';
% openEMS_opts = [openEMS_opts ' --disable-dumps']; % openEMS_opts = [openEMS_opts ' --disable-dumps'];
% openEMS_opts = [openEMS_opts ' --debug-material']; % openEMS_opts = [openEMS_opts ' --debug-material'];
@ -95,9 +94,8 @@ WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
%% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% cd to working dir and run openEMS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
savePath = pwd(); savePath = pwd();
cd(Sim_Path); %cd to working dir cd(Sim_Path); %cd to working dir
command = [openEMS_Path 'openEMS.sh ' Sim_CSX ' ' openEMS_opts]; args = [Sim_CSX ' ' openEMS_opts];
disp(command); invoke_openEMS(args);
system(command)
cd(savePath); cd(savePath);
%% do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% do the plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -107,4 +105,3 @@ PlotArgs.component=2;
PlotArgs.Limit = 'auto'; PlotArgs.Limit = 'auto';
PlotHDF5FieldData('tmp/Et.h5',PlotArgs) PlotHDF5FieldData('tmp/Et.h5',PlotArgs)

10
openems.cpp
View File

@ -244,10 +244,6 @@ int openEMS::SetupFDTD(const char* file)
return(-2); return(-2);
} }
bool PMC[6];
for (int n=0;n<6;++n)
PMC[n]=(bounds[n]==1);
//*************** setup operator ************// //*************** setup operator ************//
cout << "Create Operator..." << endl; cout << "Create Operator..." << endl;
if (CylinderCoords) if (CylinderCoords)
@ -266,8 +262,10 @@ int openEMS::SetupFDTD(const char* file)
if (FDTD_Op->SetGeometryCSX(&CSX)==false) return(2); if (FDTD_Op->SetGeometryCSX(&CSX)==false) return(2);
FDTD_Op->SetBoundaryCondition(bounds); //operator only knows about PEC and PMC, everything else is defined by extensions (see below)
/**************************** create all operator/engine extensions here !!!! **********************************/ /**************************** create all operator/engine extensions here !!!! **********************************/
//Mur-ABC //Mur-ABC, defined as extension to the operator
for (int n=0;n<6;++n) for (int n=0;n<6;++n)
{ {
if (bounds[n]==2) if (bounds[n]==2)
@ -295,8 +293,6 @@ int openEMS::SetupFDTD(const char* file)
FDTD_Op->ShowStat(); FDTD_Op->ShowStat();
FDTD_Op->ApplyMagneticBC(PMC);
cout << "Creation time for operator: " << difftime(OpDoneTime,startTime) << " s" << endl; cout << "Creation time for operator: " << difftime(OpDoneTime,startTime) << " s" << endl;
//create FDTD engine //create FDTD engine