diff --git a/FDTD/engine.cpp b/FDTD/engine.cpp index 436377a..dc8b10a 100644 --- a/FDTD/engine.cpp +++ b/FDTD/engine.cpp @@ -31,6 +31,8 @@ Engine* Engine::New(const Operator* op) Engine::Engine(const Operator* op) { + m_type = BASIC; + numTS = 0; Op = op; for (int n=0;n<3;++n) numLines[n] = Op->GetOriginalNumLines(n); diff --git a/FDTD/engine.h b/FDTD/engine.h index d28858a..b0f52ac 100644 --- a/FDTD/engine.h +++ b/FDTD/engine.h @@ -26,6 +26,11 @@ class Engine_Extension; class Engine { public: + enum EngineType + { + BASIC, SSE, UNKNOWN + }; + static Engine* New(const Operator* op); virtual ~Engine(); @@ -37,11 +42,11 @@ public: virtual unsigned int GetNumberOfTimesteps() {return numTS;}; - inline virtual FDTD_FLOAT& GetVolt( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return volt[n][x][y][z]; } - inline virtual FDTD_FLOAT& GetCurr( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return curr[n][x][y][z]; } - - inline virtual FDTD_FLOAT& GetVolt( unsigned int n, unsigned int pos[] ) { return GetVolt(n,pos[0],pos[1],pos[2]); } - inline virtual FDTD_FLOAT& GetCurr( unsigned int n, unsigned int pos[] ) { return GetCurr(n,pos[0],pos[1],pos[2]); } + //this access functions muss be overloaded by any new engine using a different storage model + inline virtual FDTD_FLOAT& GetVolt( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) { return volt[n][x][y][z]; } + inline virtual FDTD_FLOAT& GetVolt( unsigned int n, unsigned int pos[3] ) { return volt[n][pos[0]][pos[1]][pos[2]]; } + inline virtual FDTD_FLOAT& GetCurr( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) { return curr[n][x][y][z]; } + inline virtual FDTD_FLOAT& GetCurr( unsigned int n, unsigned int pos[3] ) { return curr[n][pos[0]][pos[1]][pos[2]]; } virtual void UpdateVoltages(unsigned int startX, unsigned int numX); virtual void ApplyVoltageExcite(); @@ -51,7 +56,11 @@ public: inline size_t GetExtensionCount() {return m_Eng_exts.size();} inline Engine_Extension* GetExtension(size_t nr) {return m_Eng_exts.at(nr);} + EngineType GetType() const {return m_type;} + protected: + EngineType m_type; + Engine(const Operator* op); const Operator* Op; diff --git a/FDTD/engine_ext_cylinder.cpp b/FDTD/engine_ext_cylinder.cpp index 72fc021..1ac704e 100644 --- a/FDTD/engine_ext_cylinder.cpp +++ b/FDTD/engine_ext_cylinder.cpp @@ -21,11 +21,13 @@ Engine_Ext_Cylinder::Engine_Ext_Cylinder(Operator_Ext_Cylinder* op_ext) : Engine_Extension(op_ext) { + cyl_Op = op_ext; + CC_closedAlpha = op_ext->CC_closedAlpha; CC_R0_included = op_ext->CC_R0_included; for (int n=0;n<3;++n) - numLines[n] = op_ext->m_Op->GetNumberOfLines(n); + numLines[n] = op_ext->m_Op->GetOriginalNumLines(n); } void Engine_Ext_Cylinder::Apply2Voltages() diff --git a/FDTD/engine_ext_mur_abc.cpp b/FDTD/engine_ext_mur_abc.cpp index ad77d78..23c97c8 100644 --- a/FDTD/engine_ext_mur_abc.cpp +++ b/FDTD/engine_ext_mur_abc.cpp @@ -18,6 +18,7 @@ #include "engine_ext_mur_abc.h" #include "operator_ext_mur_abc.h" #include "engine.h" +#include "engine_sse.h" #include "tools/array_ops.h" Engine_Ext_Mur_ABC::Engine_Ext_Mur_ABC(Operator_Ext_Mur_ABC* op_ext) : Engine_Extension(op_ext) @@ -72,15 +73,50 @@ void Engine_Ext_Mur_ABC::DoPreVoltageUpdates() pos[m_ny] = m_LineNr; pos_shift[m_ny] = m_LineNr_Shift; - for (pos[m_nyP]=0;pos[m_nyP]GetType()) { - pos_shift[m_nyP] = pos[m_nyP]; - for (pos[m_nyPP]=0;pos[m_nyPP]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_Op_mur->m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->GetVolt(m_nyPP,pos); + for (pos[m_nyP]=0;pos[m_nyP]Engine::GetVolt(m_nyP,pos_shift) - m_Op_mur->m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->Engine::GetVolt(m_nyP,pos); + m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]] = m_Eng->Engine::GetVolt(m_nyPP,pos_shift) - m_Op_mur->m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->Engine::GetVolt(m_nyPP,pos); + } + } + break; } + case Engine::SSE: + { + Engine_sse* eng_sse = (Engine_sse*) m_Eng; + for (pos[m_nyP]=0;pos[m_nyP]Engine_sse::GetVolt(m_nyP,pos_shift) - m_Op_mur->m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]] * eng_sse->Engine_sse::GetVolt(m_nyP,pos); + m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]] = eng_sse->Engine_sse::GetVolt(m_nyPP,pos_shift) - m_Op_mur->m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] * eng_sse->Engine_sse::GetVolt(m_nyPP,pos); + } + } + break; + } + default: + for (pos[m_nyP]=0;pos[m_nyP]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_Op_mur->m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->GetVolt(m_nyPP,pos); + } + } + break; } } @@ -93,15 +129,52 @@ void Engine_Ext_Mur_ABC::DoPostVoltageUpdates() pos[m_ny] = m_LineNr; pos_shift[m_ny] = m_LineNr_Shift; - for (pos[m_nyP]=0;pos[m_nyP]GetType()) { - pos_shift[m_nyP] = pos[m_nyP]; - for (pos[m_nyPP]=0;pos[m_nyPP]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_Op_mur->m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->GetVolt(m_nyPP,pos_shift); + for (pos[m_nyP]=0;pos[m_nyP]m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->Engine::GetVolt(m_nyP,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->Engine::GetVolt(m_nyPP,pos_shift); + } + } + break; } + + case Engine::SSE: + { + Engine_sse* eng_sse = (Engine_sse*) m_Eng; + for (pos[m_nyP]=0;pos[m_nyP]m_Mur_Coeff_nyP[pos[m_nyP]][pos[m_nyPP]] * eng_sse->Engine_sse::GetVolt(m_nyP,pos_shift); + m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]] += m_Op_mur->m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] * eng_sse->Engine_sse::GetVolt(m_nyPP,pos_shift); + } + } + break; + } + + default: + for (pos[m_nyP]=0;pos[m_nyP]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_Op_mur->m_Mur_Coeff_nyPP[pos[m_nyP]][pos[m_nyPP]] * m_Eng->GetVolt(m_nyPP,pos_shift); + } + } + break; } } @@ -112,12 +185,46 @@ void Engine_Ext_Mur_ABC::Apply2Voltages() unsigned int pos[] = {0,0,0}; pos[m_ny] = m_LineNr; - for (pos[m_nyP]=0;pos[m_nyP]GetType()) { - for (pos[m_nyPP]=0;pos[m_nyPP]GetVolt(m_nyP,pos) = m_volt_nyP[pos[m_nyP]][pos[m_nyPP]]; - m_Eng->GetVolt(m_nyPP,pos) = m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]]; + for (pos[m_nyP]=0;pos[m_nyP]Engine::GetVolt(m_nyP,pos) = m_volt_nyP[pos[m_nyP]][pos[m_nyPP]]; + m_Eng->Engine::GetVolt(m_nyPP,pos) = m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]]; + } + } + break; } + + case Engine::SSE: + { + Engine_sse* eng_sse = (Engine_sse*) m_Eng; + for (pos[m_nyP]=0;pos[m_nyP]Engine_sse::GetVolt(m_nyP,pos) = m_volt_nyP[pos[m_nyP]][pos[m_nyPP]]; + eng_sse->Engine_sse::GetVolt(m_nyPP,pos) = m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]]; + } + } + break; + } + + default: + for (pos[m_nyP]=0;pos[m_nyP]GetVolt(m_nyP,pos) = m_volt_nyP[pos[m_nyP]][pos[m_nyPP]]; + m_Eng->GetVolt(m_nyPP,pos) = m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]]; + } + } + break; } + } diff --git a/FDTD/engine_multithread.cpp b/FDTD/engine_multithread.cpp index e614789..041bcc4 100644 --- a/FDTD/engine_multithread.cpp +++ b/FDTD/engine_multithread.cpp @@ -45,6 +45,7 @@ Engine_Multithread* Engine_Multithread::New(const Operator* op, unsigned int num Engine_Multithread::Engine_Multithread(const Operator* op) : Engine(op) { + m_type = UNKNOWN; m_RunEngine = NULL; } diff --git a/FDTD/engine_multithread.h b/FDTD/engine_multithread.h index af50eda..bd5afcb 100644 --- a/FDTD/engine_multithread.h +++ b/FDTD/engine_multithread.h @@ -82,8 +82,11 @@ public: static Engine_Multithread* New(const Operator* op, unsigned int numThreads = 0); virtual ~Engine_Multithread(); - inline virtual FDTD_FLOAT& GetVolt( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return m_RunEngine->GetVolt(n,x,y,z); } - inline virtual FDTD_FLOAT& GetCurr( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return m_RunEngine->GetCurr(n,x,y,z);} + //this access functions muss be overloaded by any new engine using a different storage model + inline virtual FDTD_FLOAT& GetVolt( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) { return m_RunEngine->GetVolt(n,x,y,z); } + inline virtual FDTD_FLOAT& GetVolt( unsigned int n, unsigned int pos[3] ) { return m_RunEngine->GetVolt(n,pos); } + inline virtual FDTD_FLOAT& GetCurr( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) { return m_RunEngine->GetCurr(n,x,y,z); } + inline virtual FDTD_FLOAT& GetCurr( unsigned int n, unsigned int pos[3] ) { return m_RunEngine->GetCurr(n,pos);} virtual void setNumThreads( unsigned int numThreads ); virtual void Init(); diff --git a/FDTD/engine_sse.cpp b/FDTD/engine_sse.cpp index d7fb842..d5a8bc7 100644 --- a/FDTD/engine_sse.cpp +++ b/FDTD/engine_sse.cpp @@ -28,6 +28,7 @@ Engine_sse* Engine_sse::New(const Operator_sse* op) Engine_sse::Engine_sse(const Operator_sse* op) : Engine(op) { + m_type = SSE; Op = op; f4_volt = 0; f4_curr = 0; diff --git a/FDTD/engine_sse.h b/FDTD/engine_sse.h index 8cb517d..00cd324 100644 --- a/FDTD/engine_sse.h +++ b/FDTD/engine_sse.h @@ -32,8 +32,11 @@ public: virtual unsigned int GetNumberOfTimesteps() {return numTS;}; - inline virtual FDTD_FLOAT& GetVolt( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return f4_volt[n][x][y][z%numVectors].f[z/numVectors]; } - inline virtual FDTD_FLOAT& GetCurr( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return f4_curr[n][x][y][z%numVectors].f[z/numVectors]; } + //this access functions muss be overloaded by any new engine using a different storage model + inline virtual FDTD_FLOAT& GetVolt( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) { return f4_volt[n][x][y][z%numVectors].f[z/numVectors]; } + inline virtual FDTD_FLOAT& GetVolt( unsigned int n, unsigned int pos[3] ) { return f4_volt[n][pos[0]][pos[1]][pos[2]%numVectors].f[pos[2]/numVectors]; } + inline virtual FDTD_FLOAT& GetCurr( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) { return f4_curr[n][x][y][z%numVectors].f[z/numVectors]; } + inline virtual FDTD_FLOAT& GetCurr( unsigned int n, unsigned int pos[3] ) { return f4_curr[n][pos[0]][pos[1]][pos[2]%numVectors].f[pos[2]/numVectors]; } protected: Engine_sse(const Operator_sse* op); diff --git a/FDTD/operator_ext_cylinder.cpp b/FDTD/operator_ext_cylinder.cpp index df5860b..81e6b6b 100644 --- a/FDTD/operator_ext_cylinder.cpp +++ b/FDTD/operator_ext_cylinder.cpp @@ -29,18 +29,18 @@ bool Operator_Ext_Cylinder::BuildExtension() if ((CC_R0_included==false) || (CC_closedAlpha==false)) return true; - vv_R0 = new FDTD_FLOAT[m_Op->GetNumberOfLines(2)]; - vi_R0 = new FDTD_FLOAT[m_Op->GetNumberOfLines(2)]; + vv_R0 = new FDTD_FLOAT[m_Op->GetOriginalNumLines(2)]; + vi_R0 = new FDTD_FLOAT[m_Op->GetOriginalNumLines(2)]; unsigned int pos[3]; double inEC[4]; double dT = m_Op->GetTimestep(); pos[0]=0; - for (pos[2]=0;pos[2]GetNumberOfLines(2);++pos[2]) + for (pos[2]=0;pos[2]GetOriginalNumLines(2);++pos[2]) { double C=0; double G=0; - for (pos[1]=0;pos[1]GetNumberOfLines(1)-1;++pos[1]) + for (pos[1]=0;pos[1]GetOriginalNumLines(1)-1;++pos[1]) { m_Op_Cyl->Calc_ECPos(2,pos,inEC); C+=inEC[0]*0.5; diff --git a/FDTD/processfields.cpp b/FDTD/processfields.cpp index a4d62b2..cd01c6a 100644 --- a/FDTD/processfields.cpp +++ b/FDTD/processfields.cpp @@ -58,16 +58,21 @@ void ProcessFields::InitProcess() unsigned int* NrLines; double** Lines; - if (m_DumpMode==CELL_INTERPOLATE) - { - NrLines = numDLines; - Lines = discDLines; - } - else if (m_DumpMode==NO_INTERPOLATION) + if (m_DumpMode==NO_INTERPOLATION) { NrLines = numLines; Lines = discLines; } + else if (m_DumpMode==NODE_INTERPOLATE) + { + NrLines = numLines; + Lines = discLines; + } + else if (m_DumpMode==CELL_INTERPOLATE) + { + NrLines = numDLines; + Lines = discDLines; + } else return; diff --git a/FDTD/processfields_td.cpp b/FDTD/processfields_td.cpp index 0d0f2a7..22f2099 100644 --- a/FDTD/processfields_td.cpp +++ b/FDTD/processfields_td.cpp @@ -388,6 +388,8 @@ int ProcessFieldsTD::Process() { if (m_DumpMode==NO_INTERPOLATION) DumpNoInterpol(m_fileName); + if (m_DumpMode==NODE_INTERPOLATE) + DumpNodeInterpol(m_fileName); if (m_DumpMode==CELL_INTERPOLATE) DumpCellInterpol(m_fileName); }