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

operator: new steady state detection operator extension 

Signed-off-by: Thorsten Liebig <Thorsten.Liebig@gmx.de>
pull/13/head
Thorsten Liebig 2015-05-04 20:47:19 +02:00
parent 846e7c960e
commit cd1db5d21b
14 changed files with 443 additions and 36 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

49
FDTD/excitation.cpp
View File

@ -28,18 +28,16 @@ Excitation::Excitation()
Signal_volt = 0;
Signal_curr = 0;
m_Excit_Type = -1;
this->Reset(0);
dT = 0;
m_nyquistTS = 0;
m_f_max = 0;
m_foi = 0;
m_Excite_Elem = NULL;
m_Excit_Type = -1;
m_SignalPeriod = 0;
}
Excitation::~Excitation()
{
delete[] Signal_volt;
delete[] Signal_curr;
this->Reset(0);
}
void Excitation::Reset( double timestep )
@ -49,15 +47,13 @@ void Excitation::Reset( double timestep )
delete[] Signal_curr;
Signal_curr = 0;
m_Excit_Type = -1;
dT = timestep;
m_nyquistTS = 0;
m_f_max = 0;
m_foi = 0;
}
bool Excitation::setupExcitation( TiXmlElement* Excite, unsigned int maxTS )
bool Excitation::setupExcitation( TiXmlElement* Excite)
{
if (!Excite)
{
@ -65,6 +61,23 @@ bool Excitation::setupExcitation( TiXmlElement* Excite, unsigned int maxTS )
return false;
}
m_Excite_Elem = Excite;
double f0=0;
m_Excite_Elem->QueryIntAttribute("Type",&m_Excit_Type);
m_SignalPeriod = 0;
switch (m_Excit_Type)
{
case 1: // sinusoidal excite
m_Excite_Elem->QueryDoubleAttribute("f0",&f0);
m_SignalPeriod = 1/f0;
break;
}
return true;
}
bool Excitation::buildExcitationSignal(unsigned int maxTS)
{
if (dT<=0)
{
cerr << "Excitation::setupExcitation: Error, invalid timestep... " << endl;
@ -73,17 +86,15 @@ bool Excitation::setupExcitation( TiXmlElement* Excite, unsigned int maxTS )
double f0=0;
double fc=0;
Excite->QueryIntAttribute("Type",&m_Excit_Type);
switch (m_Excit_Type)
{
case 0:
Excite->QueryDoubleAttribute("f0",&f0);
Excite->QueryDoubleAttribute("fc",&fc);
m_Excite_Elem->QueryDoubleAttribute("f0",&f0);
m_Excite_Elem->QueryDoubleAttribute("fc",&fc);
CalcGaussianPulsExcitation(f0,fc,maxTS);
break;
case 1:
Excite->QueryDoubleAttribute("f0",&f0);
m_Excite_Elem->QueryDoubleAttribute("f0",&f0);
CalcSinusExcitation(f0,maxTS);
break;
case 2:
@ -93,18 +104,18 @@ bool Excitation::setupExcitation( TiXmlElement* Excite, unsigned int maxTS )
CalcStepExcitation();
break;
case 10:
Excite->QueryDoubleAttribute("f0",&f0);
CalcCustomExcitation(f0,maxTS,Excite->Attribute("Function"));
m_Excite_Elem->QueryDoubleAttribute("f0",&f0);
CalcCustomExcitation(f0,maxTS,m_Excite_Elem->Attribute("Function"));
break;
default:
cerr << "Excitation::setupExcitation: Unknown excitation type: \"" << m_Excit_Type<< "\" !!" << endl;
cerr << "Excitation::buildExcitationSignal: Unknown excitation type: \"" << m_Excit_Type<< "\" !!" << endl;
m_Excit_Type = -1;
return false;
}
if (GetNyquistNum() == 0)
{
cerr << "Excitation::setupExcitation: Unknown error... excitation setup failed!!" << endl;
cerr << "Excitation::buildExcitationSignal: Unknown error... excitation setup failed!!" << endl;
return false;
}

9
FDTD/excitation.h
View File

@ -32,7 +32,8 @@ public:
virtual void Reset( double timestep );
bool setupExcitation( TiXmlElement* Excite, unsigned int maxTS );
bool setupExcitation(TiXmlElement* Excite);
bool buildExcitationSignal(unsigned int maxTS);
//! Get the excitation timestep with the (first) max amplitude
virtual unsigned int GetMaxExcitationTimestep() const;
@ -61,14 +62,20 @@ public:
//! Get the frequency of interest
double GetFrequencyOfInterest() const {return m_foi;}
//! Get the signal period, 0 if not a periodical signal
double GetSignalPeriod() const {return m_SignalPeriod;}
FDTD_FLOAT* GetVoltageSignal() const {return Signal_volt;}
FDTD_FLOAT* GetCurrentSignal() const {return Signal_curr;}
protected:
double dT;
unsigned int m_nyquistTS;
double m_SignalPeriod;
int m_Excit_Type;
TiXmlElement* m_Excite_Elem;
//Excitation time-signal
unsigned int Length;
FDTD_FLOAT* Signal_volt;

2
FDTD/extensions/CMakeLists.txt
View File

@ -21,6 +21,8 @@ set(SOURCES
engine_ext_excitation.cpp
operator_ext_tfsf.cpp
engine_ext_tfsf.cpp
operator_ext_steadystate.cpp
engine_ext_steadystate.cpp
)
# FDTD/extensions lib

104
FDTD/extensions/engine_ext_steadystate.cpp Normal file
View File

@ -0,0 +1,104 @@
/*
* Copyright (C) 2015 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 <http://www.gnu.org/licenses/>.
*/
#include "engine_ext_steadystate.h"
#include "operator_ext_steadystate.h"
#include "FDTD/engine_sse.h"
#include "FDTD/engine_interface_fdtd.h"
Engine_Ext_SteadyState::Engine_Ext_SteadyState(Operator_Ext_SteadyState* op_ext): Engine_Extension(op_ext)
{
m_Op_SS = op_ext;
m_Priority = ENG_EXT_PRIO_STEADYSTATE;
for (size_t n=0;n<m_Op_SS->m_E_probe_dir.size();++n)
{
double* rec = new double[m_Op_SS->m_TS_period*2];
m_E_records.push_back(rec);
}
m_last_max_diff = 1;
last_total_energy = 0;
m_Eng_Interface = NULL;
}
Engine_Ext_SteadyState::~Engine_Ext_SteadyState()
{
for (size_t n=0;n<m_E_records.size();++n)
{
delete[] m_E_records.at(n);
m_E_records.at(n) = NULL;
}
m_E_records.clear();
delete m_Eng_Interface;
m_Eng_Interface = NULL;
}
void Engine_Ext_SteadyState::Apply2Voltages()
{
unsigned int p = m_Op_SS->m_TS_period;
unsigned int TS = m_Eng->GetNumberOfTimesteps();
unsigned int rel_pos = m_Eng->GetNumberOfTimesteps()%(2*p);
for (size_t n=0;n<m_E_records.size();++n)
m_E_records.at(n)[rel_pos] = m_Eng->GetVolt(m_Op_SS->m_E_probe_dir.at(n), m_Op_SS->m_E_probe_pos[0].at(n), m_Op_SS->m_E_probe_pos[1].at(n), m_Op_SS->m_E_probe_pos[2].at(n));
if ((TS%(m_Op_SS->m_TS_period)==0) && (TS>=2*p))
{
bool no_valid = true;
m_last_max_diff = 0;
double curr_total_energy = m_Eng_Interface->CalcFastEnergy();
if (last_total_energy>0)
{
m_last_max_diff = abs(curr_total_energy-last_total_energy)/last_total_energy;
no_valid = false;
}
// cerr << curr_total_energy << "/" << last_total_energy << "=" << abs(curr_total_energy-last_total_energy)/last_total_energy << endl;
last_total_energy = curr_total_energy;
unsigned int old_pos = 0;
unsigned int new_pos = p;
if (rel_pos<=p)
{
new_pos = 0;
old_pos = p;
}
//cerr << TS << "/" << rel_pos << ": one period complete, new_pos" << new_pos << " old pos: " << old_pos << endl;
for (size_t n=0;n<m_E_records.size();++n)
{
double *buf = m_E_records.at(n);
double curr_pow = 0;
double diff_pow = 0;
for (unsigned int nt=0;nt<p;++nt)
{
curr_pow += buf[nt+new_pos]*buf[nt+new_pos];
diff_pow += (buf[nt+old_pos]-buf[nt+new_pos])*(buf[nt+old_pos]-buf[nt+new_pos]);
}
if (curr_pow>0)
{
m_last_max_diff = max(m_last_max_diff, diff_pow/curr_pow);
no_valid = false;
}
// cerr << "curr_pow: " << curr_pow << " diff_pow: " << diff_pow << " diff: " << diff_pow/curr_pow << endl;
// cerr << m_last_max_diff << endl;
}
if ((no_valid) || (m_last_max_diff>1))
m_last_max_diff = 1;
// cerr << m_last_max_diff << endl;
}
}
void Engine_Ext_SteadyState::Apply2Current()
{
}

51
FDTD/extensions/engine_ext_steadystate.h Normal file
View File

@ -0,0 +1,51 @@
/*
* Copyright (C) 2015 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 <http://www.gnu.org/licenses/>.
*/
#ifndef ENGINE_EXT_STEADYSTATE_H
#define ENGINE_EXT_STEADYSTATE_H
#include "engine_extension.h"
#include "FDTD/engine.h"
#include "FDTD/operator.h"
class Operator_Ext_SteadyState;
class Engine_Interface_FDTD;
class Engine_Ext_SteadyState : public Engine_Extension
{
public:
Engine_Ext_SteadyState(Operator_Ext_SteadyState* op_ext);
virtual ~Engine_Ext_SteadyState();
virtual void Apply2Voltages();
virtual void Apply2Current();
void SetEngineInterface(Engine_Interface_FDTD* eng_if) {m_Eng_Interface=eng_if;}
double GetLastDiff() {return m_last_max_diff;}
protected:
Operator_Ext_SteadyState* m_Op_SS;
double m_last_max_diff;
vector<double*> m_E_records;
vector<double*> m_H_records;
double last_total_energy;
Engine_Interface_FDTD* m_Eng_Interface;
};
#endif // ENGINE_EXT_STEADYSTATE_H

1
FDTD/extensions/engine_extension.h
View File

@ -21,6 +21,7 @@
#define ENG_EXT_PRIO_DEFAULT 0 //default engine extension priority
// priority definitions for some important extensions
#define ENG_EXT_PRIO_STEADYSTATE +2e6 //steady state extension priority
#define ENG_EXT_PRIO_UPML +1e6 //unaxial pml extension priority
#define ENG_EXT_PRIO_CYLINDER +1e5 //cylindrial extension priority
#define ENG_EXT_PRIO_TFSF +5e4 //total-field/scattered-field extension priority

104
FDTD/extensions/operator_ext_steadystate.cpp Normal file
View File

@ -0,0 +1,104 @@
/*
* Copyright (C) 2015 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 <http://www.gnu.org/licenses/>.
*/
#include "operator_ext_steadystate.h"
#include "engine_ext_steadystate.h"
Operator_Ext_SteadyState::Operator_Ext_SteadyState(Operator* op, double period): Operator_Extension(op)
{
this->Reset();
m_T_period = period;
}
Operator_Ext_SteadyState::Operator_Ext_SteadyState(Operator* op, Operator_Ext_SteadyState* op_ext): Operator_Extension(op, op_ext)
{
this->Reset();
m_T_period = op_ext->m_T_period;
}
Operator_Ext_SteadyState::~Operator_Ext_SteadyState()
{
}
Operator_Extension* Operator_Ext_SteadyState::Clone(Operator* op)
{
//disable cloning, only the main operator needs to have a steady state detection
UNUSED(op);
return NULL;
}
bool Operator_Ext_SteadyState::BuildExtension()
{
double dT = m_Op->GetTimestep();
m_TS_period = round(m_T_period/dT);
m_T_period = m_TS_period*dT;
return true;
}
void Operator_Ext_SteadyState::Reset()
{
for (int n=0;n<3;++n)
{
m_E_probe_pos[n].clear();
m_H_probe_pos[n].clear();
}
m_E_probe_dir.clear();
m_H_probe_dir.clear();
m_T_period = 0;
m_TS_period = 0;
}
bool Operator_Ext_SteadyState::Add_E_Probe(unsigned int pos[3], int dir)
{
if ((dir<0) || (dir>2))
return false;
for (int n=0;n<3;++n)
if (pos[n]>=m_Op->GetNumberOfLines(n))
return false;
for (int n=0;n<3;++n)
m_E_probe_pos[n].push_back(pos[n]);
m_E_probe_dir.push_back(dir);
return true;
}
bool Operator_Ext_SteadyState::Add_H_Probe(unsigned int pos[3], int dir)
{
if ((dir<0) || (dir>2))
return false;
for (int n=0;n<3;++n)
if (pos[n]>=m_Op->GetNumberOfLines(n))
return false;
for (int n=0;n<3;++n)
m_H_probe_pos[n].push_back(pos[n]);
m_H_probe_dir.push_back(dir);
return true;
}
Engine_Extension* Operator_Ext_SteadyState::CreateEngineExtention()
{
m_Eng_Ext = new Engine_Ext_SteadyState(this);
return m_Eng_Ext;
}
void Operator_Ext_SteadyState::ShowStat(ostream &ostr) const
{
Operator_Extension::ShowStat(ostr);
cout << "Period time (s): " << m_T_period << "\t Period TS: " << m_TS_period << endl;
cout << "Number of E probes\t: " << m_E_probe_dir.size() << endl;
cout << "Number of H probes\t: " << m_H_probe_dir.size() << endl;
}

61
FDTD/extensions/operator_ext_steadystate.h Normal file
View File

@ -0,0 +1,61 @@
/*
* Copyright (C) 2015 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 <http://www.gnu.org/licenses/>.
*/
#ifndef OPERATOR_EXT_STEADYSTATE_H
#define OPERATOR_EXT_STEADYSTATE_H
#include "operator_extension.h"
#include "FDTD/operator.h"
class Engine_Ext_SteadyState;
class Operator_Ext_SteadyState : public Operator_Extension
{
friend class Engine_Ext_SteadyState;
public:
Operator_Ext_SteadyState(Operator* op, double period);
virtual ~Operator_Ext_SteadyState();
virtual Operator_Extension* Clone(Operator* op);
virtual bool BuildExtension();
virtual Engine_Extension* CreateEngineExtention();
virtual bool IsCylinderCoordsSave(bool closedAlpha, bool R0_included) const {UNUSED(closedAlpha); UNUSED(R0_included); return true;}
virtual bool IsCylindricalMultiGridSave(bool child) const {UNUSED(child); return true;}
virtual bool IsMPISave() const {return true;}
virtual string GetExtensionName() const {return string("Steady-State Detection Extension");}
virtual void ShowStat(ostream &ostr) const;
virtual void Reset();
bool Add_E_Probe(unsigned int pos[3], int dir);
bool Add_H_Probe(unsigned int pos[3], int dir);
protected:
Operator_Ext_SteadyState(Operator* op, Operator_Ext_SteadyState* op_ext);
double m_T_period;
unsigned int m_TS_period;
vector<unsigned int> m_E_probe_pos[3];
vector<unsigned int> m_E_probe_dir;
vector<unsigned int> m_H_probe_pos[3];
vector<unsigned int> m_H_probe_dir;
};
#endif // OPERATOR_EXT_STEADYSTATE_H

2
FDTD/extensions/operator_extension.cpp
View File

@ -28,6 +28,7 @@ Operator_Extension::Operator_Extension(Operator* op)
m_Op_Cyl = dynamic_cast<Operator_Cylinder*>(op);
if (m_Op_Cyl)
m_CC_R0_included=m_Op_Cyl->GetR0Included();
m_Eng_Ext = NULL;
}
Operator_Extension::~Operator_Extension()
@ -42,6 +43,7 @@ Operator_Extension::Operator_Extension(Operator* op, Operator_Extension* op_ext)
m_Active = op_ext->m_Active;
if (m_Op_Cyl)
m_CC_R0_included=m_Op_Cyl->GetR0Included();
m_Eng_Ext = NULL;
}
void Operator_Extension::ShowStat(ostream &ostr) const

2
FDTD/extensions/operator_extension.h
View File

@ -49,6 +49,7 @@ public:
virtual bool BuildExtension() {return true;}
virtual Engine_Extension* CreateEngineExtention() {return 0;}
virtual Engine_Extension* GetEngineExtention() {return m_Eng_Ext;}
//! The cylindrical operator will check whether the extension is save to use. Default is false. Derive this method to override.
virtual bool IsCylinderCoordsSave(bool closedAlpha, bool R0_included) const {UNUSED(closedAlpha); UNUSED(R0_included); return false;}
@ -78,6 +79,7 @@ protected:
//FDTD Operator
Operator* m_Op;
Engine_Extension* m_Eng_Ext;
//Cylindrical FDTD Operator (not NULL if a cylindrical FDTD is used)
Operator_Cylinder* m_Op_Cyl;

8
FDTD/operator.cpp
View File

@ -999,6 +999,14 @@ int Operator::CalcECOperator( DebugFlags debugFlags )
dT*=m_TimeStepFactor;
if (m_Exc->GetSignalPeriod()>0)
{
unsigned int TS = floor(m_Exc->GetSignalPeriod()/dT);
double new_dT = m_Exc->GetSignalPeriod()/TS;
cerr << "Operartor::CalcECOperator: Increasing timestep by " << round((new_dT-dT)/dT*1000)/10.0 << "% to " << new_dT << " (" << dT << ") to match periodic signal" << endl;
dT = new_dT;
}
m_Exc->Reset(dT);
InitOperator();

8
openEMS.pro
View File

@ -237,7 +237,9 @@ SOURCES += FDTD/extensions/engine_extension.cpp \
FDTD/extensions/operator_ext_excitation.cpp \
FDTD/extensions/engine_ext_excitation.cpp \
FDTD/extensions/operator_ext_tfsf.cpp \
FDTD/extensions/engine_ext_tfsf.cpp
FDTD/extensions/engine_ext_tfsf.cpp \
FDTD/extensions/operator_ext_steadystate.cpp \
FDTD/extensions/engine_ext_steadystate.cpp
# Common source files
SOURCES += Common/operator_base.cpp \
@ -303,7 +305,9 @@ HEADERS += FDTD/extensions/operator_extension.h \
FDTD/extensions/operator_ext_excitation.h \
FDTD/extensions/engine_ext_excitation.h \
FDTD/extensions/operator_ext_tfsf.h \
FDTD/extensions/engine_ext_tfsf.h
FDTD/extensions/engine_ext_tfsf.h \
FDTD/extensions/operator_ext_steadystate.h \
FDTD/extensions/engine_ext_steadystate.h
# Common header files
HEADERS += Common/operator_base.h \

76
openems.cpp
View File

@ -31,6 +31,8 @@
#include "FDTD/extensions/operator_ext_upml.h"
#include "FDTD/extensions/operator_ext_lorentzmaterial.h"
#include "FDTD/extensions/operator_ext_conductingsheet.h"
#include "FDTD/extensions/operator_ext_steadystate.h"
#include "FDTD/extensions/engine_ext_steadystate.h"
#include "FDTD/engine_interface_fdtd.h"
#include "FDTD/engine_interface_cylindrical_fdtd.h"
#include "Common/processvoltage.h"
@ -61,6 +63,7 @@ openEMS::openEMS()
{
FDTD_Op=NULL;
FDTD_Eng=NULL;
Eng_Ext_SSD=NULL;
m_CSX=NULL;
PA=NULL;
CylinderCoords = false;
@ -313,7 +316,7 @@ Engine_Interface_FDTD* openEMS::NewEngineInterface(int multigridlevel)
if (mgl==multigridlevel)
{
if (g_settings.GetVerboseLevel()>0)
cerr << __func__ << ": Operator with requested multi-grid level found." << endl;
cout << __func__ << ": Operator with requested multi-grid level found." << endl;
return new Engine_Interface_Cylindrical_FDTD(op_cyl_mg);
}
Operator_Cylinder* op_cyl_inner = op_cyl_mg->GetInnerOperator();
@ -321,7 +324,7 @@ Engine_Interface_FDTD* openEMS::NewEngineInterface(int multigridlevel)
if (op_cyl_mg==NULL) //inner most operator reached
{
if (g_settings.GetVerboseLevel()>0)
cerr << __func__ << ": Operator with highest multi-grid level chosen." << endl;
cout << __func__ << ": Operator with highest multi-grid level chosen." << endl;
return new Engine_Interface_Cylindrical_FDTD(op_cyl_inner);
}
// try next level
@ -641,9 +644,6 @@ int openEMS::SetupFDTD(const char* file)
if (m_OverSampling<2)
m_OverSampling=2;
double maxTime=0;
FDTD_Opts->QueryDoubleAttribute("MaxTime",&maxTime);
TiXmlElement* BC = FDTD_Opts->FirstChildElement("BoundaryCond");
if (BC==NULL)
{
@ -688,10 +688,12 @@ int openEMS::SetupFDTD(const char* file)
{
FDTD_Op->SetCellConstantMaterial();
if (g_settings.GetVerboseLevel()>0)
cerr << "Enabling constant cell material assumption." << endl;
cout << "Enabling constant cell material assumption." << endl;
}
m_Exc = new Excitation();
if (!m_Exc->setupExcitation(FDTD_Opts->FirstChildElement("Excitation")))
exit(2);
FDTD_Op->SetExcitationSignal(m_Exc);
FDTD_Op->AddExtension(new Operator_Ext_Excitation(FDTD_Op));
if (!CylinderCoords)
@ -713,6 +715,36 @@ int openEMS::SetupFDTD(const char* file)
if (FDTD_Opts->QueryDoubleAttribute("TimeStepFactor",&timestepfactor)==TIXML_SUCCESS)
FDTD_Op->SetTimestepFactor(timestepfactor);
// Is a steady state detection requested
Operator_Ext_SteadyState* Op_Ext_SSD = NULL;
if (m_Exc->GetSignalPeriod()>0)
{
cout << "Create a steady state detection using a period of " << m_Exc->GetSignalPeriod() << " s" << endl;
Op_Ext_SSD = new Operator_Ext_SteadyState(FDTD_Op, m_Exc->GetSignalPeriod());
unsigned int pos[3];
for (int p=0;p<3;++p)
pos[p] = FDTD_Op->GetNumberOfLines(p)/2;
Op_Ext_SSD->Add_E_Probe(pos, 0);
Op_Ext_SSD->Add_E_Probe(pos, 1);
Op_Ext_SSD->Add_E_Probe(pos, 2);
for (int n=0;n<3;++n)
{
for (int p=0;p<3;++p)
pos[p] = FDTD_Op->GetNumberOfLines(p)/2;
pos[n] = 10;
Op_Ext_SSD->Add_E_Probe(pos, 0);
Op_Ext_SSD->Add_E_Probe(pos, 1);
Op_Ext_SSD->Add_E_Probe(pos, 2);
pos[n] = FDTD_Op->GetNumberOfLines(n)-10;
Op_Ext_SSD->Add_E_Probe(pos, 0);
Op_Ext_SSD->Add_E_Probe(pos, 1);
Op_Ext_SSD->Add_E_Probe(pos, 2);
}
FDTD_Op->AddExtension(Op_Ext_SSD);
}
if ((m_CSX->GetQtyPropertyType(CSProperties::LORENTZMATERIAL)>0) || (m_CSX->GetQtyPropertyType(CSProperties::DEBYEMATERIAL)>0))
FDTD_Op->AddExtension(new Operator_Ext_LorentzMaterial(FDTD_Op));
if (m_CSX->GetQtyPropertyType(CSProperties::CONDUCTINGSHEET)>0)
@ -746,11 +778,13 @@ int openEMS::SetupFDTD(const char* file)
FDTD_Op->SetMaterialStoreFlags(2,false);
FDTD_Op->SetMaterialStoreFlags(3,false);
double maxTime=0;
FDTD_Opts->QueryDoubleAttribute("MaxTime",&maxTime);
unsigned int maxTime_TS = (unsigned int)(maxTime/FDTD_Op->GetTimestep());
if ((maxTime_TS>0) && (maxTime_TS<NrTS))
NrTS = maxTime_TS;
if (!m_Exc->setupExcitation(FDTD_Opts->FirstChildElement("Excitation"),NrTS))
if (!m_Exc->buildExcitationSignal(NrTS))
exit(2);
m_Exc->DumpVoltageExcite("et");
m_Exc->DumpCurrentExcite("ht");
@ -784,6 +818,12 @@ int openEMS::SetupFDTD(const char* file)
//create FDTD engine
FDTD_Eng = FDTD_Op->CreateEngine();
if (Op_Ext_SSD)
{
Eng_Ext_SSD = dynamic_cast<Engine_Ext_SteadyState*>(Op_Ext_SSD->GetEngineExtention());
Eng_Ext_SSD->SetEngineInterface(this->NewEngineInterface());
}
//setup all processing classes
if (SetupProcessing()==false)
return 2;
@ -879,7 +919,7 @@ void openEMS::RunFDTD()
FDTD_Eng->IterateTS(step);
step=PA->Process();
if (ProcField->CheckTimestep())
if ((Eng_Ext_SSD==NULL) && ProcField->CheckTimestep())
{
currE = ProcField->CalcTotalEnergyEstimate();
if (currE>maxE)
@ -896,16 +936,24 @@ void openEMS::RunFDTD()
if (t_diff>4)
{
currE = ProcField->CalcTotalEnergyEstimate();
if (currE>maxE)
maxE=currE;
t_run = CalcDiffTime(currTime,startTime);
speed = numCells*(currTS-prevTS)/t_diff;
cout << "[@" << FormatTime(t_run) << "] Timestep: " << setw(12) << currTS ;
cout << " || Speed: " << setw(6) << setprecision(1) << std::fixed << speed*1e-6 << " MC/s (" << setw(4) << setprecision(3) << std::scientific << t_diff/(currTS-prevTS) << " s/TS)" ;
if (maxE)
change = currE/maxE;
cout << " || Energy: ~" << setw(6) << setprecision(2) << std::scientific << currE << " (-" << setw(5) << setprecision(2) << std::fixed << fabs(10.0*log10(change)) << "dB)" << endl;
if (Eng_Ext_SSD==NULL)
{
currE = ProcField->CalcTotalEnergyEstimate();
if (currE>maxE)
maxE=currE;
if (maxE)
change = currE/maxE;
cout << " || Energy: ~" << setw(6) << setprecision(2) << std::scientific << currE << " (-" << setw(5) << setprecision(2) << std::fixed << fabs(10.0*log10(change)) << "dB)" << endl;
}
else
{
change = Eng_Ext_SSD->GetLastDiff();
cout << " || SteadyState: " << setw(6) << setprecision(2) << std::fixed << 10.0*log10(change) << " dB" << endl;
}
prevTime=currTime;
prevTS=currTS;

2
openems.h
View File

@ -35,6 +35,7 @@ class TiXmlElement;
class ContinuousStructure;
class Engine_Interface_FDTD;
class Excitation;
class Engine_Ext_SteadyState;
double CalcDiffTime(timeval t1, timeval t2);
string FormatTime(int sec);
@ -91,6 +92,7 @@ protected:
bool m_CellConstantMaterial;
Operator* FDTD_Op;
Engine* FDTD_Eng;
Engine_Ext_SteadyState* Eng_Ext_SSD;
ProcessingArray* PA;
Excitation* m_Exc;