diff --git a/FDTD/engine.cpp b/FDTD/engine.cpp index 21c6c3a..cfbbe06 100644 --- a/FDTD/engine.cpp +++ b/FDTD/engine.cpp @@ -113,19 +113,19 @@ void Engine::ApplyVoltageExcite() { int exc_pos; //soft voltage excitation here (E-field excite) - for (unsigned int n=0;nE_Exc_Count;++n) + for (unsigned int n=0;nExc->E_Count;++n) { - exc_pos = (int)numTS - (int)Op->E_Exc_delay[n]; - exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->ExciteLength); + exc_pos = (int)numTS - (int)Op->Exc->E_delay[n]; + exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->Exc->Length); // if (n==0) cerr << numTS << " => " << Op->ExciteSignal[exc_pos] << endl; - GetVolt(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]; + GetVolt(Op->Exc->E_dir[n],Op->Exc->E_index[0][n],Op->Exc->E_index[1][n],Op->Exc->E_index[2][n]) += Op->Exc->E_amp[n]*Op->Exc->Signal_volt[exc_pos]; } // write the first excitation into the file "et1" - if (Op->E_Exc_Count >= 1) { - exc_pos = (int)numTS - (int)Op->E_Exc_delay[0]; - exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->ExciteLength); - 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 + if (Op->Exc->E_Count >= 1) { + exc_pos = (int)numTS - (int)Op->Exc->E_delay[0]; + exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->Exc->Length); + file_et1 << numTS * Op->GetTimestep() << "\t" << Op->Exc->E_amp[0]*Op->Exc->Signal_volt[exc_pos] << "\n"; // do not use std::endl here, because it will do an implicit flush } } @@ -159,12 +159,12 @@ void Engine::ApplyCurrentExcite() { int exc_pos; //soft current excitation here (H-field excite) - for (unsigned int n=0;nCurr_Exc_Count;++n) + for (unsigned int n=0;nExc->Curr_Count;++n) { - exc_pos = (int)numTS - (int)Op->Curr_Exc_delay[n]; - exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->ExciteLength); + exc_pos = (int)numTS - (int)Op->Exc->Curr_delay[n]; + exc_pos *= (exc_pos>0 && exc_pos<=(int)Op->Exc->Length); // if (n==0) cerr << numTS << " => " << Op->ExciteSignal[exc_pos] << endl; - GetCurr(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]; + GetCurr(Op->Exc->Curr_dir[n],Op->Exc->Curr_index[0][n],Op->Exc->Curr_index[1][n],Op->Exc->Curr_index[2][n]) += Op->Exc->Curr_amp[n]*Op->Exc->Signal_curr[exc_pos]; } } diff --git a/FDTD/engine_ext_mur_abc.cpp b/FDTD/engine_ext_mur_abc.cpp index a5c3c19..ad77d78 100644 --- a/FDTD/engine_ext_mur_abc.cpp +++ b/FDTD/engine_ext_mur_abc.cpp @@ -39,18 +39,18 @@ Engine_Ext_Mur_ABC::Engine_Ext_Mur_ABC(Operator_Ext_Mur_ABC* op_ext) : Engine_Ex //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;nm_Op->E_Exc_Count;++n) + for (unsigned int n=0;nm_Op->Exc->E_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 ( ((m_Op_mur->m_Op->Exc->E_dir[n]==m_nyP) || (m_Op_mur->m_Op->Exc->E_dir[n]==m_nyPP)) && (m_Op_mur->m_Op->Exc->E_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]; + if ((int)m_Op_mur->m_Op->Exc->E_delay[n]>maxDelay) + maxDelay = (int)m_Op_mur->m_Op->Exc->E_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 + m_start_TS = maxDelay + m_Op_mur->m_Op->Exc->Length + 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; } } diff --git a/FDTD/excitation.cpp b/FDTD/excitation.cpp new file mode 100644 index 0000000..6be7717 --- /dev/null +++ b/FDTD/excitation.cpp @@ -0,0 +1,296 @@ +/* +* Copyright (C) 2010 Thorsten Liebig (Thorsten.Liebig@gmx.de) +* +* This program is free software: you can redistribute it and/or modify +* it under the terms of the GNU General Public License as published by +* the Free Software Foundation, either version 3 of the License, or +* (at your option) any later version. +* +* This program is distributed in the hope that it will be useful, +* but WITHOUT ANY WARRANTY; without even the implied warranty of +* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +* GNU General Public License for more details. +* +* You should have received a copy of the GNU General Public License +* along with this program. If not, see . +*/ + +#include "tools/array_ops.h" +#include "fparser.hh" +#include "tinyxml.h" +#include "excitation.h" + +Excitation::Excitation( double timestep ) +{ + Signal_volt = 0; + Signal_curr = 0; + E_delay = 0; + E_amp = 0; + E_dir = 0; + Curr_delay = 0; + Curr_amp = 0; + Curr_dir = 0; + + for (int n=0;n<3;++n) { + E_index[n] = 0; + Curr_index[n] = 0; + } + + dT = timestep; + m_nyquistTS = nan(""); +} + +Excitation::~Excitation() +{ + delete[] Signal_volt; + delete[] Signal_curr; + delete[] E_delay; + delete[] E_dir; + delete[] E_amp; + delete[] Curr_delay; + delete[] Curr_dir; + delete[] Curr_amp; + for (int n=0;n<3;++n) { + delete[] E_index[n]; + delete[] Curr_index[n]; + } + +} + +bool Excitation::setupExcitation( TiXmlElement* Excite, unsigned int maxTS ) +{ + if (!Excite) { + cerr << "Can't read openEMS excitation settings... " << endl; + return false; + } + + int Excit_Type=0; + double f0=0; + double fc=0; + Excite->QueryIntAttribute("Type",&Excit_Type); + + switch (Excit_Type) + { + case 0: + Excite->QueryDoubleAttribute("f0",&f0); + Excite->QueryDoubleAttribute("fc",&fc); + CalcGaussianPulsExcitation(f0,fc); + break; + case 1: + Excite->QueryDoubleAttribute("f0",&f0); + CalcSinusExcitation(f0,maxTS); + break; + case 2: + CalcDiracPulsExcitation(); + break; + case 3: + CalcStepExcitation(); + break; + case 10: + Excite->QueryDoubleAttribute("f0",&f0); + CalcCustomExcitation(f0,maxTS,Excite->Attribute("Function")); + break; + } + + if (GetNyquistNum() == 0) { + cerr << "openEMS: excitation setup failed!!" << endl; + return false; + } + + return true; +} + +unsigned int Excitation::CalcNyquistNum(double fmax) const +{ + if (dT==0) return 1; + double T0 = 1/fmax; + return floor(T0/2/dT); +} + +unsigned int Excitation::GetMaxExcitationTimestep() const +{ + FDTD_FLOAT maxAmp=0; + unsigned int maxStep=0; + for (unsigned int n=1;nmaxAmp) + { + maxAmp = fabs(Signal_volt[n]); + maxStep = n; + } + } + return maxStep; +} + +void Excitation::CalcGaussianPulsExcitation(double f0, double fc) +{ + if (dT==0) return; + + Length = (unsigned int)(2.0 * 9.0/(2.0*PI*fc) / dT); +// cerr << "Operator::CalcGaussianPulsExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl; + delete[] Signal_volt; + delete[] Signal_curr; + Signal_volt = new FDTD_FLOAT[Length+1]; + Signal_curr = new FDTD_FLOAT[Length+1]; + Signal_volt[0]=0.0; + Signal_curr[0]=0.0; + for (unsigned int n=1;n const volt_vIndex[3], vector const& volt_vExcit, + vector const& volt_vDelay, vector const& volt_vDir ) +{ + E_Count = volt_vIndex[0].size(); + for (int n=0; n<3; n++) { + delete[] E_index[n]; + E_index[n] = new unsigned int[E_Count]; + } + delete[] E_delay; + delete[] E_amp; + delete[] E_dir; + E_delay = new unsigned int[E_Count]; + E_amp = new FDTD_FLOAT[E_Count]; + E_dir = new unsigned short[E_Count]; + + cerr << "Excitation::setupVoltageExcitation(): Number of voltage excitation points: " << E_Count << endl; +// if (E_Count==0) +// cerr << "No E-Field/voltage excitation found!" << endl; + for (int n=0; n<3; n++) + for (unsigned int i=0; i const curr_vIndex[3], vector const& curr_vExcit, + vector const& curr_vDelay, vector const& curr_vDir ) +{ + Curr_Count = curr_vIndex[0].size(); + for (int n=0; n<3; n++) { + delete[] Curr_index[n]; + Curr_index[n] = new unsigned int[Curr_Count]; + } + delete[] Curr_delay; + delete[] Curr_amp; + delete[] Curr_dir; + Curr_delay = new unsigned int[Curr_Count]; + Curr_amp = new FDTD_FLOAT[Curr_Count]; + Curr_dir = new unsigned short[Curr_Count]; + + cerr << "Excitation::setupCurrentExcitation(): Number of current excitation points: " << Curr_Count << endl; +// if (Curr_Count==0) +// cerr << "No H-Field/current excitation found!" << endl; + for (int n=0;n<3;++n) + for (unsigned int i=0; i. +*/ + +#ifndef EXCITATION_H +#define EXCITATION_H + +#include +#include +#include "tools/constants.h" + +class TiXmlElement; + +class Excitation +{ +public: + Excitation( double timestep ); + virtual ~Excitation(); + + bool setupExcitation( TiXmlElement* Excite, unsigned int maxTS ); + + //! Get the excitation timestep with the (first) max amplitude + virtual unsigned int GetMaxExcitationTimestep() const; + unsigned int CalcNyquistNum(double fmax) const; + + void setupVoltageExcitation( vector const volt_vIndex[3], vector const& volt_vExcit, + vector const& volt_vDelay, vector const& volt_vDir ); + void setupCurrentExcitation( vector const curr_vIndex[3], vector const& curr_vExcit, + vector const& curr_vDelay, vector const& curr_vDir ); + + void SetNyquistNum(unsigned int nyquist) {m_nyquistTS=nyquist;} + unsigned int GetNyquistNum() const {return m_nyquistTS;} + + //Excitation time-signal + unsigned int Length; + FDTD_FLOAT* Signal_volt; + FDTD_FLOAT* Signal_curr; + + //E-Field/voltage Excitation + unsigned int E_Count; + unsigned int* E_index[3]; + unsigned short* E_dir; + FDTD_FLOAT* E_amp; //represented as edge-voltages!! + unsigned int* E_delay; + + //H-Field/current Excitation + unsigned int Curr_Count; + unsigned int* Curr_index[3]; + unsigned short* Curr_dir; + FDTD_FLOAT* Curr_amp; //represented as edge-currents!! + unsigned int* Curr_delay; + +protected: + double dT; + unsigned int m_nyquistTS; + + //! Calculate a custom signal + virtual void CalcCustomExcitation(double f0, int nTS, string signal); + //! Calculate an excitation with center of f0 and the half bandwidth fc + virtual void CalcGaussianPulsExcitation(double f0, double fc); + //! Calculate a sinusoidal excitation with frequency f0 and a duration of nTS number of timesteps + virtual void CalcSinusExcitation(double f0, int nTS); + //! Calculate a dirac impuls excitation + virtual void CalcDiracPulsExcitation(); + //! Calculate a step excitation + virtual void CalcStepExcitation(); +}; + +#endif // EXCITATION_H diff --git a/FDTD/operator.cpp b/FDTD/operator.cpp index 032276a..22b69c9 100644 --- a/FDTD/operator.cpp +++ b/FDTD/operator.cpp @@ -22,6 +22,7 @@ #include "tools/array_ops.h" #include "fparser.hh" + Operator* Operator::New() { Operator* op = new Operator(); @@ -31,6 +32,7 @@ Operator* Operator::New() Operator::Operator() { + Exc = 0; } Operator::~Operator() @@ -45,24 +47,12 @@ void Operator::Init() { CSX = NULL; - ExciteSignal_volt = NULL; - ExciteSignal_curr = NULL; - E_Exc_delay = NULL; - E_Exc_amp=NULL; - E_Exc_dir=NULL; vv=NULL; vi=NULL; - Curr_Exc_delay = NULL; - Curr_Exc_amp=NULL; - Curr_Exc_dir=NULL; iv=NULL; ii=NULL; for (int n=0;n<3;++n) - { discLines[n]=NULL; - E_Exc_index[n]=NULL; - Curr_Exc_index[n]=NULL; - } MainOp=NULL; DualOp=NULL; @@ -77,28 +67,18 @@ void Operator::Init() for (int n=0;n<6;++n) m_BC[n]=0; + + Exc = 0; } void Operator::Reset() { - delete[] ExciteSignal_volt; - delete[] ExciteSignal_curr; - delete[] E_Exc_delay; - delete[] E_Exc_dir; - delete[] E_Exc_amp; - delete[] Curr_Exc_delay; - delete[] Curr_Exc_dir; - delete[] Curr_Exc_amp; Delete_N_3DArray(vv,numLines); Delete_N_3DArray(vi,numLines); Delete_N_3DArray(iv,numLines); Delete_N_3DArray(ii,numLines); for (int n=0;n<3;++n) - { delete[] discLines[n]; - delete[] E_Exc_index[n]; - delete[] Curr_Exc_index[n]; - } delete MainOp; delete DualOp; for (int n=0;n<3;++n) @@ -109,14 +89,9 @@ void Operator::Reset() delete[] EC_R[n]; } - Init(); -} + delete Exc; -unsigned int Operator::CalcNyquistNum(double fmax) -{ - if (dT==0) return 1; - double T0 = 1/fmax; - return floor(T0/2/dT); + Init(); } string Operator::GetDirName(int ny) const @@ -310,140 +285,13 @@ void Operator::ShowStat() const cout << "Size of Field-Data : " << FieldSize << " Byte (" << (double)FieldSize/MBdiff << " MB) " << endl; cout << "-----------------------------------" << endl; cout << "Timestep (s) : " << dT << endl; - cout << "Nyquist criteria (TS) : " << m_nyquistTS << endl; - cout << "Nyquist criteria (s) : " << m_nyquistTS*dT << endl; - cout << "Excitation Length (TS) : " << ExciteLength << endl; - cout << "Excitation Length (s) : " << ExciteLength*dT << endl; + cout << "Nyquist criteria (TS) : " << Exc->GetNyquistNum() << endl; + cout << "Nyquist criteria (s) : " << Exc->GetNyquistNum()*dT << endl; + cout << "Excitation Length (TS) : " << Exc->Length << endl; + cout << "Excitation Length (s) : " << Exc->Length*dT << endl; cout << "-----------------------------------" << endl; } -unsigned int Operator::GetMaxExcitationTimestep() const -{ - FDTD_FLOAT maxAmp=0; - unsigned int maxStep=0; - for (unsigned int n=1;nmaxAmp) - { - maxAmp = fabs(ExciteSignal_volt[n]); - maxStep = n; - } - } - return maxStep; -} - -unsigned int Operator::CalcGaussianPulsExcitation(double f0, double fc) -{ - if (dT==0) return 0; - - ExciteLength = (unsigned int)(2.0 * 9.0/(2.0*PI*fc) / dT); -// cerr << "Operator::CalcGaussianPulsExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl; - delete[] ExciteSignal_volt; - delete[] ExciteSignal_curr; - 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;nE_Count;++n) + exc[Exc->E_dir[n]][Exc->E_index[0][n]][Exc->E_index[1][n]][Exc->E_index[2][n]] = Exc->E_amp[n]; } string names[] = {"vv", "vi", "iv" , "ii", "exc"}; @@ -557,6 +405,12 @@ void Operator::InitOperator() ii = Create_N_3DArray(numLines); } +void Operator::InitExcitation() +{ + delete Exc; + Exc = new Excitation( dT ); +} + void Operator::Calc_ECOperatorPos(int n, unsigned int* pos) { unsigned int i = MainOp->SetPos(pos[0],pos[1],pos[2]); @@ -576,6 +430,8 @@ int Operator::CalcECOperator() InitOperator(); + InitExcitation(); + unsigned int pos[3]; for (int n=0;n<3;++n) @@ -925,7 +781,11 @@ double Operator::CalcTimestep() bool Operator::CalcFieldExcitation() { - if (dT==0) return false; + if (dT==0) + return false; + if (Exc==0) + return false; + unsigned int pos[3]; double delta[3]; double amp=0; @@ -1094,59 +954,16 @@ bool Operator::CalcFieldExcitation() // 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) - cerr << "No E-Field/voltage excitation found!" << endl; - for (int n=0;n<3;++n) - { - delete[] E_Exc_index[n]; - E_Exc_index[n] = new unsigned int[E_Exc_Count]; - for (unsigned int i=0;isetupVoltageExcitation( volt_vIndex, volt_vExcit, volt_vDelay, volt_vDir ); // 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;isetupCurrentExcitation( curr_vIndex, curr_vExcit, curr_vDelay, curr_vDir ); return true; } + + bool Operator::CalcPEC() { unsigned int pos[3]; diff --git a/FDTD/operator.h b/FDTD/operator.h index 22fe808..5865b17 100644 --- a/FDTD/operator.h +++ b/FDTD/operator.h @@ -21,8 +21,7 @@ #include "ContinuousStructure.h" #include "tools/AdrOp.h" #include "tools/constants.h" - -#define FDTD_FLOAT float +#include "excitation.h" class Operator_Extension; @@ -45,21 +44,6 @@ public: inline virtual FDTD_FLOAT& GetII( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return ii[n][x][y][z]; } inline virtual FDTD_FLOAT& GetIV( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return iv[n][x][y][z]; } - - //! Calculate a custom signal \return number of Nyquist timesteps defined by f0 - virtual unsigned int CalcCustomExcitation(double f0, int nTS, string signal); - //! Calculate an excitation with center of f0 and the half bandwidth fc \return number of Nyquist timesteps - virtual unsigned int CalcGaussianPulsExcitation(double f0, double fc); - //! Calculate a sinusoidal excitation with frequency f0 and a duration of nTS number of timesteps \return number of Nyquist timesteps - virtual unsigned int CalcSinusExcitation(double f0, int nTS); - //! Calculate a dirac impuls excitation \return number of Nyquist timesteps - virtual unsigned int CalcDiracPulsExcitation(); - //! Calculate a step excitation \return number of Nyquist timesteps - virtual unsigned int CalcStepExcitation(); - - //! Get the excitation timestep with the (first) max amplitude - 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 ApplyMagneticBC(bool* dirs); @@ -69,10 +53,6 @@ public: virtual unsigned int GetNumberOfLines(int ny) const {return numLines[ny];} - void SetNyquistNum(unsigned int nyquist) {m_nyquistTS=nyquist;} - unsigned int GetNyquistNum() const {return m_nyquistTS;} - unsigned int CalcNyquistNum(double fmax); - void ShowStat() const; void DumpOperator2File(string filename); @@ -102,6 +82,7 @@ protected: virtual void Init(); virtual void Reset(); virtual void InitOperator(); + virtual void InitExcitation(); struct Grid_Path { @@ -122,7 +103,6 @@ protected: //Calc timestep only internal use virtual double CalcTimestep(); double dT; //FDTD timestep! - unsigned int m_nyquistTS; //! Calc operator at certain pos virtual void Calc_ECOperatorPos(int n, unsigned int* pos); @@ -152,24 +132,7 @@ public: FDTD_FLOAT**** ii; //calc new current from old current FDTD_FLOAT**** iv; //calc new current from old voltage - //Excitation time-signal - unsigned int ExciteLength; - FDTD_FLOAT* ExciteSignal_volt; - FDTD_FLOAT* ExciteSignal_curr; - - //E-Field/voltage Excitation - unsigned int E_Exc_Count; - unsigned int* E_Exc_index[3]; - unsigned short* E_Exc_dir; - FDTD_FLOAT* E_Exc_amp; //represented as edge-voltages!! - 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; + Excitation* Exc; }; #endif // OPERATOR_H diff --git a/FDTD/operator_sse.cpp b/FDTD/operator_sse.cpp index 1d097a1..8b09d33 100644 --- a/FDTD/operator_sse.cpp +++ b/FDTD/operator_sse.cpp @@ -49,8 +49,11 @@ void Operator_sse::Reset() Delete_N_3DArray_v4sf(f4_vi,numLines); Delete_N_3DArray_v4sf(f4_iv,numLines); Delete_N_3DArray_v4sf(f4_ii,numLines); + f4_vv = 0; + f4_vi = 0; + f4_iv = 0; + f4_ii = 0; Operator::Reset(); -// Init(); // FIXME this calls Operator::Init() twice... } void Operator_sse::InitOperator() diff --git a/openEMS.pro b/openEMS.pro index 07c8b28..9c38386 100644 --- a/openEMS.pro +++ b/openEMS.pro @@ -46,12 +46,13 @@ SOURCES += main.cpp \ FDTD/engine_multithread.cpp \ FDTD/operator_cylinder.cpp \ FDTD/engine_cylinder.cpp \ - FDTD/engine_sse.cpp \ - FDTD/operator_sse.cpp \ + FDTD/engine_sse.cpp \ + FDTD/operator_sse.cpp \ FDTD/operator_extension.cpp \ FDTD/engine_extension.cpp \ FDTD/engine_ext_mur_abc.cpp \ - FDTD/operator_ext_mur_abc.cpp + FDTD/operator_ext_mur_abc.cpp \ + FDTD/excitation.cpp HEADERS += tools/ErrorMsg.h \ tools/AdrOp.h \ tools/constants.h \ @@ -67,71 +68,73 @@ HEADERS += tools/ErrorMsg.h \ openems.h \ FDTD/engine_multithread.h \ FDTD/operator_cylinder.h \ - FDTD/engine_sse.h \ - FDTD/operator_sse.h \ + FDTD/engine_sse.h \ + FDTD/operator_sse.h \ FDTD/engine_cylinder.h \ FDTD/operator_extension.h \ FDTD/engine_extension.h \ FDTD/engine_ext_mur_abc.h \ - FDTD/operator_ext_mur_abc.h + FDTD/operator_ext_mur_abc.h \ + FDTD/excitation.h QMAKE_CXXFLAGS_RELEASE = -O2 \ -g \ - -march=native + -march=native -msse -msse2 QMAKE_CXXFLAGS_DEBUG = -O0 \ -g \ - -march=native + -march=native -msse -msse2 - -# # to use ABI2 target: # qmake CONFIG+="ABI2 bits64" -o Makefile.ABI2-64 openEMS.pro # make -fMakefile.ABI2-64 -# - -ABI2 { - CONFIG-=debug debug_and_release - CONFIG+=release - QMAKE_CFLAGS_RELEASE=-O2 -fabi-version=2 - QMAKE_CXXFLAGS_RELEASE=-O2 -fabi-version=2 - QMAKE_CC = apgcc - QMAKE_CXX = apg++ - QMAKE_LINK = apg++ - QMAKE_LINK_SHLIB = apg++ - QMAKE_LFLAGS_RPATH = - QMAKE_LFLAGS = \'-Wl,-rpath,\$$ORIGIN/lib\' +ABI2 { + CONFIG -= debug \ + debug_and_release + CONFIG += release + QMAKE_CFLAGS_RELEASE = -O2 \ + -fabi-version=2 + QMAKE_CXXFLAGS_RELEASE = -O2 \ + -fabi-version=2 + QMAKE_CC = apgcc + QMAKE_CXX = apg++ + QMAKE_LINK = apg++ + QMAKE_LINK_SHLIB = apg++ + QMAKE_LFLAGS_RPATH = + QMAKE_LFLAGS = \'-Wl,-rpath,\$$ORIGIN/lib\' } - -bits64 { - QMAKE_CXXFLAGS_RELEASE+=-m64 -march=athlon64 - QMAKE_LFLAGS_RELEASE+=-m64 -march=athlon64 - OBJECTS_DIR = ABI2-64 - LIBS = ../CSXCAD/ABI2-64/libCSXCAD.so - LIBS += ../fparser/ABI2-64/libfparser.so - LIBS += ../tinyxml/ABI2-64/libtinyxml.so - LIBS += ../boost-64/lib/libboost_thread.so - LIBS += ../hdf5-64/lib/libhdf5.so - LIBS += ../hdf5-64/lib/libhdf5_cpp.so -lpthread - INCLUDEPATH += ../hdf5-64/include - INCLUDEPATH += ../boost-64/include +bits64 { + QMAKE_CXXFLAGS_RELEASE += -m64 \ + -march=athlon64 + QMAKE_LFLAGS_RELEASE += -m64 \ + -march=athlon64 + OBJECTS_DIR = ABI2-64 + LIBS = ../CSXCAD/ABI2-64/libCSXCAD.so + LIBS += ../fparser/ABI2-64/libfparser.so + LIBS += ../tinyxml/ABI2-64/libtinyxml.so + LIBS += ../boost-64/lib/libboost_thread.so + LIBS += ../hdf5-64/lib/libhdf5.so + LIBS += ../hdf5-64/lib/libhdf5_cpp.so \ + -lpthread + INCLUDEPATH += ../hdf5-64/include + INCLUDEPATH += ../boost-64/include } - -bits32 { - QMAKE_CXXFLAGS_RELEASE+=-m32 -march=i686 - QMAKE_LFLAGS_RELEASE+=-m32 -march=i686 - OBJECTS_DIR = ABI2-32 - LIBS = ../CSXCAD/ABI2-32/libCSXCAD.so - LIBS += ../fparser/ABI2-32/libfparser.so - LIBS += ../tinyxml/ABI2-32/libtinyxml.so - LIBS += ../boost-32/lib/libboost_thread.so - LIBS += ../hdf5-32/lib/libhdf5.so - LIBS += ../hdf5-32/lib/libhdf5_cpp.so - INCLUDEPATH += ../hdf5-32/include - INCLUDEPATH += ../boost-32/include +bits32 { + QMAKE_CXXFLAGS_RELEASE += -m32 \ + -march=i686 + QMAKE_LFLAGS_RELEASE += -m32 \ + -march=i686 + OBJECTS_DIR = ABI2-32 + LIBS = ../CSXCAD/ABI2-32/libCSXCAD.so + LIBS += ../fparser/ABI2-32/libfparser.so + LIBS += ../tinyxml/ABI2-32/libtinyxml.so + LIBS += ../boost-32/lib/libboost_thread.so + LIBS += ../hdf5-32/lib/libhdf5.so + LIBS += ../hdf5-32/lib/libhdf5_cpp.so + INCLUDEPATH += ../hdf5-32/include + INCLUDEPATH += ../boost-32/include } - -ABI2 { - DESTDIR = $$OBJECTS_DIR - MOC_DIR = $$OBJECTS_DIR - UI_DIR = $$OBJECTS_DIR - RCC_DIR = $$OBJECTS_DIR +ABI2 { + DESTDIR = $$OBJECTS_DIR + MOC_DIR = $$OBJECTS_DIR + UI_DIR = $$OBJECTS_DIR + RCC_DIR = $$OBJECTS_DIR } diff --git a/openems.cpp b/openems.cpp index cab2722..21c45db 100644 --- a/openems.cpp +++ b/openems.cpp @@ -124,51 +124,6 @@ bool openEMS::parseCommandLineArgument( const char *argv ) return false; } -void openEMS::SetupExcitation(TiXmlElement* Excite) -{ - if (Excite==NULL) - { - cerr << "Can't read openEMS Excitation Settings... " << endl; - exit(-2); - } - - int Excit_Type=0; - double f0=0; - double fc=0; - Excite->QueryIntAttribute("Type",&Excit_Type); - - unsigned int Nyquist = 0; - switch (Excit_Type) - { - case 0: - Excite->QueryDoubleAttribute("f0",&f0); - Excite->QueryDoubleAttribute("fc",&fc); - Nyquist = FDTD_Op->CalcGaussianPulsExcitation(f0,fc); - break; - case 1: - Excite->QueryDoubleAttribute("f0",&f0); - Nyquist = FDTD_Op->CalcSinusExcitation(f0,NrTS); - break; - case 2: - Nyquist = FDTD_Op->CalcDiracPulsExcitation(); - break; - case 3: - Nyquist = FDTD_Op->CalcStepExcitation(); - break; - case 10: - Excite->QueryDoubleAttribute("f0",&f0); - Nyquist = FDTD_Op->CalcCustomExcitation(f0,NrTS,Excite->Attribute("Function")); - break; - } - - if (!Nyquist) - { - cerr << "openEMS: excitation setup failed!!" << endl; - exit(2); - } - FDTD_Op->SetNyquistNum(Nyquist); -} - int openEMS::SetupFDTD(const char* file) { if (file==NULL) return -1; @@ -278,7 +233,8 @@ int openEMS::SetupFDTD(const char* file) FDTD_Op->CalcECOperator(); - SetupExcitation(FDTD_Opts->FirstChildElement("Excitation")); + if (!FDTD_Op->Exc->setupExcitation( FDTD_Opts->FirstChildElement("Excitation"), NrTS )) + exit(2); if (DebugMat) { @@ -318,7 +274,7 @@ int openEMS::SetupFDTD(const char* file) //*************** setup processing ************// cout << "Setting up processing..." << endl; - unsigned int Nyquist = FDTD_Op->GetNyquistNum(); + unsigned int Nyquist = FDTD_Op->Exc->GetNyquistNum(); PA = new ProcessingArray(Nyquist); double start[3]; @@ -414,9 +370,9 @@ void openEMS::RunFDTD() double maxE=0,currE=0; //add all timesteps to end-crit field processing with max excite amplitude - unsigned int maxExcite = FDTD_Op->GetMaxExcitationTimestep(); - for (unsigned int n=0;nE_Exc_Count;++n) - ProcField->AddStep(FDTD_Op->E_Exc_delay[n]+maxExcite); + unsigned int maxExcite = FDTD_Op->Exc->GetMaxExcitationTimestep(); + for (unsigned int n=0;nExc->E_Count;++n) + ProcField->AddStep(FDTD_Op->Exc->E_delay[n]+maxExcite); double change=1; int prevTS=0,currTS=0; diff --git a/openems.h b/openems.h index 2ccd7f3..c710dda 100644 --- a/openems.h +++ b/openems.h @@ -45,8 +45,6 @@ public: void DebugBox() {m_debugBox=true;} protected: - void SetupExcitation(TiXmlElement* Excite); - bool CylinderCoords; //! Number of Timesteps diff --git a/tools/constants.h b/tools/constants.h index ab65d4f..9707fdf 100644 --- a/tools/constants.h +++ b/tools/constants.h @@ -18,6 +18,8 @@ #ifndef CONSTANTS_H #define CONSTANTS_H +#define FDTD_FLOAT float + #define __EPS0__ 8.85418781762e-12 #define __MUE0__ 1.256637062e-6 #define __C0__ 299792458