continued (internal) reorganization of excitation

FDTD/excitation.cpp

@@ -23,14 +23,14 @@
 #include "tinyxml.h"
 #include "excitation.h"
 
-Excitation::Excitation( double timestep )
+Excitation::Excitation()
 {
 	Signal_volt = 0;
 	Signal_curr = 0;
 
 	m_Excit_Type = -1;
 
-	dT = timestep;
+	dT = 0;
 	m_nyquistTS = 0;
 }
 
@@ -57,7 +57,13 @@ bool Excitation::setupExcitation( TiXmlElement* Excite, unsigned int maxTS )
 {
 	if (!Excite)
 	{
-		cerr << "Can't read openEMS excitation settings... " << endl;
+		cerr << "Excitation::setupExcitation: Error, can't read openEMS excitation settings... " << endl;
+		return false;
+	}
+
+	if (dT<=0)
+	{
+		cerr << "Excitation::setupExcitation: Error, invalid timestep... " << endl;
 		return false;
 	}
 

FDTD/excitation.h

@@ -27,7 +27,7 @@ class TiXmlElement;
 class Excitation
 {
 public:
-	Excitation( double timestep );
+	Excitation();
 	virtual ~Excitation();
 
 	virtual void Reset( double timestep );

FDTD/extensions/engine_ext_mur_abc.cpp

@@ -42,7 +42,6 @@ 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;
-	Excitation* Exc = m_Op_mur->m_Op->GetExcitationSignal();
 	Operator_Ext_Excitation* Exc_ext = m_Op_mur->m_Op->GetExcitationExtension();
 	for (unsigned int n=0; n<Exc_ext->GetVoltCount(); ++n)
 	{

FDTD/extensions/operator_ext_excitation.cpp

@@ -22,6 +22,23 @@
 
 Operator_Ext_Excitation::Operator_Ext_Excitation(Operator* op) : Operator_Extension(op)
 {
+	Init();
+}
+
+Operator_Ext_Excitation::~Operator_Ext_Excitation()
+{
+	Reset();
+}
+
+Operator_Extension* Operator_Ext_Excitation::Clone(Operator* op)
+{
+	Operator_Ext_Excitation* clone = new Operator_Ext_Excitation(op, this);
+	return clone;
+}
+
+void Operator_Ext_Excitation::Init()
+{
+	Operator_Extension::Init();
 	Volt_delay = 0;
 	Volt_amp = 0;
 	Volt_dir = 0;
@@ -38,17 +55,12 @@ Operator_Ext_Excitation::Operator_Ext_Excitation(Operator* op) : Operator_Extens
 		Volt_Count_Dir[n] = 0;
 		Curr_Count_Dir[n] = 0;
 	}
-
+	m_Exc = 0;
-	m_Exc = m_Op->m_Exc;
 }
 
-Operator_Ext_Excitation::~Operator_Ext_Excitation()
+void Operator_Ext_Excitation::Reset()
-{
-	Reset();
-}
-
-void Operator_Ext_Excitation::Reset( )
 {
+	Operator_Extension::Reset();
 	delete[] Volt_delay;
 	Volt_delay = 0;
 	delete[] Volt_dir;
@@ -80,17 +92,19 @@ void Operator_Ext_Excitation::Reset( )
 
 Operator_Ext_Excitation::Operator_Ext_Excitation(Operator* op, Operator_Ext_Excitation* op_ext) : Operator_Extension(op, op_ext)
 {
-	m_Exc = NULL;
+	Init();
 }
 
 bool Operator_Ext_Excitation::BuildExtension()
 {
+	m_Exc = m_Op->GetExcitationSignal();
 	double dT = m_Op->GetTimestep();
 	if (dT==0)
 		return false;
 	if (m_Exc==0)
 		return false;
 
+	Reset();
 	ContinuousStructure* CSX = m_Op->GetGeometryCSX();
 
 	unsigned int pos[3];

FDTD/extensions/operator_ext_excitation.h

@@ -32,7 +32,7 @@ public:
 	Operator_Ext_Excitation(Operator* op);
 	~Operator_Ext_Excitation();
 
-	virtual Operator_Extension* Clone(Operator* op) {UNUSED(op);return NULL;}
+	virtual Operator_Extension* Clone(Operator* op);
 
 	virtual bool BuildExtension();
 
@@ -45,6 +45,7 @@ public:
 
 	virtual void ShowStat(ostream &ostr) const;
 
+	virtual void Init();
 	virtual void Reset();
 
 	unsigned int GetVoltCount() const {return Volt_Count;}

FDTD/extensions/operator_extension.cpp

@@ -22,6 +22,7 @@
 Operator_Extension::Operator_Extension(Operator* op)
 {
 	m_Op = op;
+	m_Active = true;
 
 	m_CC_R0_included = false;
 	m_Op_Cyl = dynamic_cast<Operator_Cylinder*>(op);

FDTD/extensions/operator_extension.h

@@ -60,6 +60,10 @@ public:
 
 	virtual void ShowStat(ostream &ostr) const;
 
+	virtual bool IsActive() const {return m_Active;}
+	virtual void SetActive(bool active=true) {m_Active=active;}
+
+	virtual void Init() {}
 	virtual void Reset() {}
 
 protected:
@@ -67,6 +71,8 @@ protected:
 	//! Copy constructor
 	Operator_Extension(Operator* op, Operator_Extension* op_ext);
 
+	bool m_Active;
+
 	//FDTD Operator
 	Operator* m_Op;
 

FDTD/operator.cpp

@@ -104,8 +104,6 @@ void Operator::Delete()
 		delete[] EC_R[n];EC_R[n]=0;
 	}
 
-	delete m_Exc;m_Exc=0;
-
 	Delete_N_3DArray(m_epsR,numLines);
 	m_epsR=0;
 	Delete_N_3DArray(m_kappa,numLines);
@@ -489,23 +487,26 @@ void Operator::DumpOperator2File(string filename)
 
 	vtk_Writer->SetNativeDump(true);
 
-	if (m_Op_Ext_Exc)
+	//find excitation extension
+	Operator_Ext_Excitation* Op_Ext_Exc=GetExcitationExtension();
+
+	if (Op_Ext_Exc)
 	{
 		FDTD_FLOAT**** exc = NULL;
-		if (m_Op_Ext_Exc->Volt_Count>0)
+		if (Op_Ext_Exc->Volt_Count>0)
 		{
 			exc = Create_N_3DArray<FDTD_FLOAT>(numLines);
-			for (unsigned int n=0; n<m_Op_Ext_Exc->Volt_Count; ++n)
+			for (unsigned int n=0; n<  Op_Ext_Exc->Volt_Count; ++n)
-				exc[m_Op_Ext_Exc->Volt_dir[n]][m_Op_Ext_Exc->Volt_index[0][n]][m_Op_Ext_Exc->Volt_index[1][n]][m_Op_Ext_Exc->Volt_index[2][n]] = m_Op_Ext_Exc->Volt_amp[n];
+				exc[  Op_Ext_Exc->Volt_dir[n]][  Op_Ext_Exc->Volt_index[0][n]][  Op_Ext_Exc->Volt_index[1][n]][  Op_Ext_Exc->Volt_index[2][n]] =   Op_Ext_Exc->Volt_amp[n];
 			vtk_Writer->AddVectorField("exc_volt",exc);
 			Delete_N_3DArray(exc,numLines);
 		}
 
-		if (m_Op_Ext_Exc->Curr_Count>0)
+		if (  Op_Ext_Exc->Curr_Count>0)
 		{
 			exc = Create_N_3DArray<FDTD_FLOAT>(numLines);
-			for (unsigned int n=0; n<m_Op_Ext_Exc->Curr_Count; ++n)
+			for (unsigned int n=0; n<  Op_Ext_Exc->Curr_Count; ++n)
-				exc[m_Op_Ext_Exc->Curr_dir[n]][m_Op_Ext_Exc->Curr_index[0][n]][m_Op_Ext_Exc->Curr_index[1][n]][m_Op_Ext_Exc->Curr_index[2][n]] = m_Op_Ext_Exc->Curr_amp[n];
+				exc[  Op_Ext_Exc->Curr_dir[n]][  Op_Ext_Exc->Curr_index[0][n]][  Op_Ext_Exc->Curr_index[1][n]][  Op_Ext_Exc->Curr_index[2][n]] =   Op_Ext_Exc->Curr_amp[n];
 			vtk_Writer->AddVectorField("exc_curr",exc);
 			Delete_N_3DArray(exc,numLines);
 		}
@@ -829,16 +830,9 @@ double Operator::GetDiscMaterial(int type, int n, const unsigned int pos[3]) con
 	return 0;
 }
 
-void Operator::InitExcitation()
+void Operator::SetExcitationSignal(Excitation* exc)
 {
-	if (m_Exc!=NULL)
+	m_Exc=exc;
-		m_Exc->Reset(dT);
-	else
-	{
-		m_Exc = new Excitation( dT );
-		m_Op_Ext_Exc = new Operator_Ext_Excitation(this);
-		this->AddExtension(m_Op_Ext_Exc);
-	}
 }
 
 void Operator::Calc_ECOperatorPos(int n, unsigned int* pos)
@@ -897,6 +891,8 @@ int Operator::CalcECOperator( DebugFlags debugFlags )
 	else
 		CalcTimestep();
 
+	m_Exc->Reset(dT);
+
 	InitOperator();
 
 	unsigned int pos[3];
@@ -932,12 +928,23 @@ int Operator::CalcECOperator( DebugFlags debugFlags )
 		PMC[n] = m_BC[n]==1;
 	ApplyMagneticBC(PMC);
 
-	InitExcitation();
-
 	//all information available for extension... create now...
 	for (size_t n=0; n<m_Op_exts.size(); ++n)
 		m_Op_exts.at(n)->BuildExtension();
 
+	//remove inactive extensions
+	vector<Operator_Extension*>::iterator it = m_Op_exts.begin();
+	while (it!=m_Op_exts.end())
+	{
+		if ( (*it)->IsActive() == false)
+		{
+			m_Op_exts.erase(it);
+			it = m_Op_exts.begin(); //restart search for inactive extension
+		}
+		else
+			++it;
+	}
+
 	if (debugFlags & debugMaterial)
 		DumpMaterial2File( "material_dump" );
 	if (debugFlags & debugOperator)
@@ -1422,12 +1429,6 @@ bool Operator::Calc_LumpedElements()
 	return true;
 }
 
-void Operator::DumpExciationSignals()
-{
-	m_Exc->DumpVoltageExcite("et");
-	m_Exc->DumpCurrentExcite("ht");
-}
-
 void Operator::Init_EC()
 {
 	for (int n=0; n<3; ++n)
@@ -1710,6 +1711,19 @@ void Operator::CalcPEC_Curves()
 	}
 }
 
+Operator_Ext_Excitation* Operator::GetExcitationExtension() const
+{
+	//search for excitation extension
+	Operator_Ext_Excitation* Op_Ext_Exc=0;
+	for (size_t n=0; n<m_Op_exts.size(); ++n)
+	{
+		Op_Ext_Exc = dynamic_cast<Operator_Ext_Excitation*>(m_Op_exts.at(n));
+		if (Op_Ext_Exc)
+			break;
+	}
+	return Op_Ext_Exc;
+}
+
 void Operator::AddExtension(Operator_Extension* op_ext)
 {
 	m_Op_exts.push_back(op_ext);

FDTD/operator.h

@@ -52,10 +52,6 @@ public:
 
 	virtual int CalcECOperator( DebugFlags debugFlags = None );
 
-	virtual bool SetupExcitation(TiXmlElement* Excite, unsigned int maxTS) {return m_Exc->setupExcitation(Excite,maxTS);}
-
-	virtual void DumpExciationSignals();
-
 	// the next four functions need to be reimplemented in a derived class
 	inline virtual FDTD_FLOAT GetVV( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return vv[n][x][y][z]; }
 	inline virtual FDTD_FLOAT GetVI( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return vi[n][x][y][z]; }
@@ -146,9 +142,10 @@ public:
 
 	virtual double GetDiscMaterial(int type, int ny, const unsigned int pos[3]) const;
 
-	Excitation* GetExcitationSignal() const {return m_Exc;}
+	virtual void SetExcitationSignal(Excitation* exc);
+	virtual Excitation* GetExcitationSignal() const {return m_Exc;}
 
-	Operator_Ext_Excitation* GetExcitationExtension() const {return m_Op_Ext_Exc;}
+	Operator_Ext_Excitation* GetExcitationExtension() const;
 
 protected:
 	//! use New() for creating a new Operator
@@ -159,7 +156,6 @@ protected:
 	virtual void Reset();
 	virtual void InitOperator();
 	virtual void InitDataStorage();
-	virtual void InitExcitation();
 
 	virtual bool SetupCSXGrid(CSRectGrid* grid);
 
@@ -236,7 +232,7 @@ protected:
 
 	// excitation classes
 	Excitation* m_Exc; // excitation time signal class
-	Operator_Ext_Excitation* m_Op_Ext_Exc; // excitation extension
+//	Operator_Ext_Excitation* m_Op_Ext_Exc; // excitation extension
 
 	// engine/post-proc needs access
 public:

FDTD/operator_cylindermultigrid.cpp

@@ -445,11 +445,10 @@ void Operator_CylinderMultiGrid::SetupInterpolation()
 	}
 }
 
-bool Operator_CylinderMultiGrid::SetupExcitation(TiXmlElement* Excite, unsigned int maxTS)
+void Operator_CylinderMultiGrid::SetExcitationSignal(Excitation* exc)
 {
-	if (!m_InnerOp->SetupExcitation(Excite,maxTS))
+	m_InnerOp->SetExcitationSignal(exc);
-		return false;
+	Operator_Cylinder::SetExcitationSignal(exc);
-	return m_Exc->setupExcitation(Excite,maxTS);
 }
 
 void Operator_CylinderMultiGrid::Delete()

FDTD/operator_cylindermultigrid.h

@@ -46,14 +46,14 @@ public:
 
 	virtual unsigned int GetSplitPos() const {return m_Split_Pos;}
 
-	virtual bool SetupExcitation(TiXmlElement* Excite, unsigned int maxTS);
-
 	virtual void SetBoundaryCondition(int* BCs);
 
 	virtual void AddExtension(Operator_Extension* op_ext);
 
 	Operator_Cylinder* GetInnerOperator() const {return m_InnerOp;}
 
+	virtual void SetExcitationSignal(Excitation* exc);
+
 	virtual void ShowStat() const;
 
 #ifdef MPI_SUPPORT

FDTD/operator_mpi.cpp

@@ -98,12 +98,6 @@ void Operator_MPI::SetNeighborDown(int ny, int id)
 	m_NeighborDown[ny]=id;
 }
 
-void Operator_MPI::DumpExciationSignals()
-{
-	if (m_MyID==0)
-		Operator_SSE_Compressed::DumpExciationSignals();
-}
-
 void Operator_MPI::Init()
 {
 	Operator_SSE_Compressed::Init();

FDTD/operator_mpi.h

@@ -46,8 +46,6 @@ public:
 	virtual void SetNeighborUp(int ny, int id);
 	virtual void SetNeighborDown(int ny, int id);
 
-	virtual void DumpExciationSignals();
-
 	//! Set the lower original mesh position
 	virtual void SetSplitPos(int ny, unsigned int pos) {m_SplitPos[ny]=pos;}
 	virtual void SetOriginalMesh(CSRectGrid* orig_Mesh);

openems.cpp

@@ -23,6 +23,7 @@
 #include "FDTD/operator_cylindermultigrid.h"
 #include "FDTD/engine_multithread.h"
 #include "FDTD/operator_multithread.h"
+#include "FDTD/extensions/operator_ext_excitation.h"
 #include "FDTD/extensions/operator_ext_mur_abc.h"
 #include "FDTD/extensions/operator_ext_pml_sf.h"
 #include "FDTD/extensions/operator_ext_upml.h"
@@ -71,6 +72,7 @@ openEMS::openEMS()
 	m_engine_numThreads = 0;
 
 	m_Abort = false;
+	m_Exc = 0;
 }
 
 openEMS::~openEMS()
@@ -89,6 +91,8 @@ void openEMS::Reset()
 	FDTD_Op=0;
 	delete m_CSX;
 	m_CSX=0;
+	delete m_Exc;
+	m_Exc=0;
 }
 
 //! \brief processes a command line argument
@@ -624,6 +628,10 @@ int openEMS::SetupFDTD(const char* file)
 	if (SetupOperator(FDTD_Opts)==false)
 		return 2;
 
+	m_Exc = new Excitation();
+	FDTD_Op->SetExcitationSignal(m_Exc);
+	FDTD_Op->AddExtension(new Operator_Ext_Excitation(FDTD_Op));
+
 	if (FDTD_Op->SetGeometryCSX(m_CSX)==false) return(2);
 
 	SetupBoundaryConditions(BC);
@@ -667,10 +675,10 @@ int openEMS::SetupFDTD(const char* file)
 	if ((maxTime_TS>0) && (maxTime_TS<NrTS))
 		NrTS = maxTime_TS;
 
-	if (!FDTD_Op->SetupExcitation( FDTD_Opts->FirstChildElement("Excitation"), NrTS ))
+	if (!m_Exc->setupExcitation(FDTD_Opts->FirstChildElement("Excitation"),NrTS))
 		exit(2);
-
+	m_Exc->DumpVoltageExcite("et");
-	FDTD_Op->DumpExciationSignals();
+	m_Exc->DumpCurrentExcite("ht");
 
 	timeval OpDoneTime;
 	gettimeofday(&OpDoneTime,NULL);
@@ -681,15 +689,14 @@ int openEMS::SetupFDTD(const char* file)
 		FDTD_Op->ShowExtStat();
 		cout << "Creation time for operator: " << CalcDiffTime(OpDoneTime,startTime) << " s" << endl;
 	}
-	Excitation* Exc=FDTD_Op->GetExcitationSignal();
 	cout << "FDTD simulation size: " << FDTD_Op->GetNumberOfLines(0) << "x" << FDTD_Op->GetNumberOfLines(1) << "x" << FDTD_Op->GetNumberOfLines(2) << " --> "  << FDTD_Op->GetNumberCells() << " FDTD cells " << endl;
-	cout << "FDTD timestep is: " <<FDTD_Op->GetTimestep()  << " s; Nyquist rate: " << Exc->GetNyquistNum() << " timesteps @" << CalcNyquistFrequency(Exc->GetNyquistNum(),FDTD_Op->GetTimestep()) << " Hz" << endl;
+	cout << "FDTD timestep is: " <<FDTD_Op->GetTimestep()  << " s; Nyquist rate: " <<  m_Exc->GetNyquistNum() << " timesteps @" << CalcNyquistFrequency(m_Exc->GetNyquistNum(),FDTD_Op->GetTimestep()) << " Hz" << endl;
-	if (Exc->GetNyquistNum()>1000)
+	if (m_Exc->GetNyquistNum()>1000)
 		cerr << "openEMS::SetupFDTD: Warning, the timestep seems to be very small --> long simulation. Check your mesh!?" << endl;
 
-	cout << "Excitation signal length is: " << Exc->GetLength() << " timesteps (" << Exc->GetLength()*FDTD_Op->GetTimestep() << "s)" << endl;
+	cout << "Excitation signal length is: " <<  m_Exc->GetLength() << " timesteps (" <<  m_Exc->GetLength()*FDTD_Op->GetTimestep() << "s)" << endl;
-	cout << "Max. number of timesteps: " << NrTS << " ( --> " << (double)NrTS/(double)(Exc->GetLength()) << " * Excitation signal length)" << endl;
+	cout << "Max. number of timesteps: " << NrTS << " ( --> " << (double)NrTS/(double)(m_Exc->GetLength()) << " * Excitation signal length)" << endl;
-	if ( ((double)NrTS/(double)Exc->GetLength() < 3) && (Exc->GetExciteType()==0))
+	if ( ((double)NrTS/(double)m_Exc->GetLength() < 3) && (m_Exc->GetExciteType()==0))
 		cerr << "openEMS::SetupFDTD: Warning, max. number of timesteps is smaller than three times the excitation. " << endl << \
 				"\tYou may want to choose a higher number of max. timesteps... " << endl;
 

openems.h

@@ -31,6 +31,7 @@ class ProcessingArray;
 class TiXmlElement;
 class ContinuousStructure;
 class Engine_Interface_FDTD;
+class Excitation;
 
 double CalcDiffTime(timeval t1, timeval t2);
 string FormatTime(int sec);
@@ -84,6 +85,8 @@ protected:
 	Engine* FDTD_Eng;
 	ProcessingArray* PA;
 
+	Excitation* m_Exc;
+
 	bool m_Abort;
 
 #ifdef MPI_SUPPORT