diff --git a/FDTD/engine.cpp b/FDTD/engine.cpp index cfbbe06..4b8d569 100644 --- a/FDTD/engine.cpp +++ b/FDTD/engine.cpp @@ -34,7 +34,8 @@ Engine::Engine(const Operator* op) Op = op; for (int n=0;n<3;++n) { - numLines[n] = Op->GetNumberOfLines(n); +// numLines[n] = Op->GetNumberOfLines(n); + numLines[n] = Op->numLines[n]; } for (size_t n=0;nGetNumberOfExtentions();++n) diff --git a/FDTD/engine_cylinder.cpp b/FDTD/engine_cylinder.cpp deleted file mode 100644 index 99d349a..0000000 --- a/FDTD/engine_cylinder.cpp +++ /dev/null @@ -1,136 +0,0 @@ -/* -* 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 "engine_cylinder.h" -#include "engine_extension.h" -#include "operator_extension.h" - -Engine_Cylinder* Engine_Cylinder::New(const Operator_Cylinder* op) -{ - Engine_Cylinder* e = new Engine_Cylinder(op); - e->Init(); - return e; -} - -Engine_Cylinder::Engine_Cylinder(const Operator_Cylinder* op) : Engine(op) -{ - cyl_Op = op; - if (cyl_Op->GetClosedAlpha()) - { - ++numLines[1]; //necessary for dobled voltage and current line in alpha-dir, operator will return one smaller for correct post-processing - } -} - -inline void Engine_Cylinder::CloseAlphaVoltages() -{ - unsigned int pos[3]; - // copy voltages from last alpha-plane to first - unsigned int last_A_Line = numLines[1]-1; - for (pos[0]=0;pos[0]vv_R0[pos[2]]; - for (pos[1]=0;pos[1]GetClosedAlpha();++pos[1]) - { - volt[2][0][0][pos[2]] += cyl_Op->vi_R0[pos[2]] * curr[1][0][pos[1]][pos[2]]; - } - } - for (pos[1]=0;pos[1]GetClosedAlpha()==false) - return Engine::IterateTS(iterTS); - - for (unsigned int iter=0;iterDoPreVoltageUpdates(); - - UpdateVoltages(); - - for (size_t n=0;nDoPostVoltageUpdates(); - for (size_t n=0;nApply2Voltages(); - - if (cyl_Op->GetR0Included()) - R0IncludeVoltages(); - - ApplyVoltageExcite(); - - CloseAlphaVoltages(); - - //current updates with extensions - for (size_t n=0;nDoPreCurrentUpdates(); - - UpdateCurrents(); - - for (size_t n=0;nDoPostCurrentUpdates(); - for (size_t n=0;nApply2Current(); - - ApplyCurrentExcite(); - - CloseAlphaCurrents(); - - ++numTS; - } - - return true; -} diff --git a/FDTD/engine_ext_cylinder.cpp b/FDTD/engine_ext_cylinder.cpp new file mode 100644 index 0000000..72fc021 --- /dev/null +++ b/FDTD/engine_ext_cylinder.cpp @@ -0,0 +1,89 @@ +/* +* 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 "engine.h" +#include "engine_ext_cylinder.h" +#include "operator_ext_cylinder.h" + +Engine_Ext_Cylinder::Engine_Ext_Cylinder(Operator_Ext_Cylinder* op_ext) : Engine_Extension(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); +} + +void Engine_Ext_Cylinder::Apply2Voltages() +{ + if (CC_closedAlpha==false) return; + + if (CC_R0_included) + { + unsigned int pos[3]; + pos[0] = 0; + for (pos[2]=0;pos[2]GetVolt(2,0,0,pos[2]) *= cyl_Op->vv_R0[pos[2]]; + for (pos[1]=0;pos[1]GetVolt(2,0,0,pos[2]) += cyl_Op->vi_R0[pos[2]] * m_Eng->GetCurr(1,0,pos[1],pos[2]); + } + } + for (pos[1]=0;pos[1]GetVolt(1,0,pos[1],pos[2]) = 0; //no voltage in alpha-direction at r=0 + m_Eng->GetVolt(2,0,pos[1],pos[2]) = m_Eng->GetVolt(2,0,0,pos[2]); + } + } + } + + //close alpha + unsigned int pos[3]; + // copy voltages from last alpha-plane to first + unsigned int last_A_Line = numLines[1]-1; + for (pos[0]=0;pos[0]GetVolt(0,pos[0],0,pos[2]) = m_Eng->GetVolt(0,pos[0],last_A_Line,pos[2]); + m_Eng->GetVolt(1,pos[0],0,pos[2]) = m_Eng->GetVolt(1,pos[0],last_A_Line,pos[2]); + m_Eng->GetVolt(2,pos[0],0,pos[2]) = m_Eng->GetVolt(2,pos[0],last_A_Line,pos[2]); + } + } +} + +void Engine_Ext_Cylinder::Apply2Current() +{ + if (CC_closedAlpha==false) return; + + //close alpha + unsigned int pos[3]; + // copy currents from first alpha-plane to last + for (pos[0]=0;pos[0]GetCurr(0,pos[0],last_A_Line,pos[2]) = m_Eng->GetCurr(0,pos[0],0,pos[2]); + m_Eng->GetCurr(1,pos[0],last_A_Line,pos[2]) = m_Eng->GetCurr(1,pos[0],0,pos[2]); + m_Eng->GetCurr(2,pos[0],last_A_Line,pos[2]) = m_Eng->GetCurr(2,pos[0],0,pos[2]); + } + } +} diff --git a/FDTD/engine_cylinder.h b/FDTD/engine_ext_cylinder.h similarity index 69% rename from FDTD/engine_cylinder.h rename to FDTD/engine_ext_cylinder.h index df93b30..a3ec90f 100644 --- a/FDTD/engine_cylinder.h +++ b/FDTD/engine_ext_cylinder.h @@ -19,25 +19,27 @@ #define ENGINE_CYLINDER_H #include "engine.h" +#include "engine_extension.h" #include "operator_cylinder.h" -class Engine_Cylinder : public Engine +class Operator_Ext_Cylinder; + +class Engine_Ext_Cylinder : public Engine_Extension { public: - static Engine_Cylinder* New(const Operator_Cylinder* op); + Engine_Ext_Cylinder(Operator_Ext_Cylinder* op_ext); - //!Iterate a number of timesteps - virtual bool IterateTS(unsigned int iterTS); + virtual void Apply2Voltages(); + + virtual void Apply2Current(); protected: - Engine_Cylinder(const Operator_Cylinder* op); + Operator_Ext_Cylinder* cyl_Op; - virtual inline void CloseAlphaVoltages(); - virtual inline void CloseAlphaCurrents(); + unsigned int numLines[3]; - virtual inline void R0IncludeVoltages(); - - const Operator_Cylinder* cyl_Op; + bool CC_closedAlpha; + bool CC_R0_included; }; #endif // ENGINE_CYLINDER_H diff --git a/FDTD/operator.h b/FDTD/operator.h index ce00963..b924923 100644 --- a/FDTD/operator.h +++ b/FDTD/operator.h @@ -28,6 +28,7 @@ class Operator_Extension; //! Abstract base-class for the FDTD-operator class Operator { + friend class Engine; public: //! Create a new operator static Operator* New(); diff --git a/FDTD/operator_cylinder.cpp b/FDTD/operator_cylinder.cpp index 579be89..e7ba875 100644 --- a/FDTD/operator_cylinder.cpp +++ b/FDTD/operator_cylinder.cpp @@ -17,6 +17,7 @@ #include "operator_cylinder.h" #include "operator_extension.h" +#include "operator_ext_cylinder.h" Operator_Cylinder* Operator_Cylinder::New() { @@ -38,26 +39,20 @@ void Operator_Cylinder::Init() { CC_closedAlpha = false; CC_R0_included = false; - vv_R0 = NULL; - vi_R0 = NULL; Operator::Init(); } void Operator_Cylinder::Reset() { Operator::Reset(); - delete[] vv_R0;vv_R0=NULL; - delete[] vi_R0;vi_R0=NULL; } void Operator_Cylinder::InitOperator() { - if (CC_R0_included) - { - vv_R0 = new FDTD_FLOAT[numLines[2]]; - vi_R0 = new FDTD_FLOAT[numLines[2]]; - } Operator::InitOperator(); + + if (CC_closedAlpha || CC_R0_included) + this->AddExtension(new Operator_Ext_Cylinder(this)); } inline unsigned int Operator_Cylinder::GetNumberOfLines(int ny) const @@ -123,43 +118,12 @@ bool Operator_Cylinder::SetGeometryCSX(ContinuousStructure* geo) else if (discLines[0][0]==0.0) { cout << "Operator_Cylinder::SetGeometryCSX: r=0 included..." << endl; - CC_R0_included=true; + CC_R0_included= true; //also needed for correct ec-calculation } return true; } -int Operator_Cylinder::CalcECOperator() -{ - int val = Operator::CalcECOperator(); - if (val) - return val; - - //if r=0 is not included -> obviously no special treatment for r=0 - //if alpha direction is not closed, PEC-BC at r=0 necessary and already set... - if ((CC_R0_included==false) || (CC_closedAlpha==false)) - return val; - - unsigned int pos[3]; - double inEC[4]; - pos[0]=0; - for (pos[2]=0;pos[2]GetR0Included(); + CC_closedAlpha=m_Op_Cyl->GetClosedAlpha(); + + vv_R0 = NULL; + vi_R0 = NULL; +} + +Operator_Ext_Cylinder::~Operator_Ext_Cylinder() +{ + delete[] vv_R0;vv_R0=NULL; + delete[] vi_R0;vi_R0=NULL; +} + +bool Operator_Ext_Cylinder::BuildExtension() +{ + delete[] vv_R0;vv_R0=NULL; + delete[] vi_R0;vi_R0=NULL; + + //if r=0 is not included -> obviously no special treatment for r=0 + //if alpha direction is not closed, PEC-BC at r=0 necessary and already set... + 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)]; + + 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]) + { + double C=0; + double G=0; + for (pos[1]=0;pos[1]GetNumberOfLines(1)-1;++pos[1]) + { + m_Op_Cyl->Calc_ECPos(2,pos,inEC); + C+=inEC[0]*0.5; + G+=inEC[1]*0.5; + } + m_Op->GetVV(2,0,0,pos[2]) = 1; + vv_R0[pos[2]] = (1-dT*G/2/C)/(1+dT*G/2/C); + vi_R0[pos[2]] = (dT/C)/(1+dT*G/2/C); + } + return true; +} + +Engine_Extension* Operator_Ext_Cylinder::CreateEngineExtention() +{ + Engine_Ext_Cylinder* eng_ext = new Engine_Ext_Cylinder(this); + return eng_ext; +} diff --git a/FDTD/operator_ext_cylinder.h b/FDTD/operator_ext_cylinder.h new file mode 100644 index 0000000..e306913 --- /dev/null +++ b/FDTD/operator_ext_cylinder.h @@ -0,0 +1,36 @@ +#ifndef OPERATOR_EXT_CYLINDER_H +#define OPERATOR_EXT_CYLINDER_H + +#include "operator_extension.h" +#include "operator.h" + +class Operator_Cylinder; + +class Operator_Ext_Cylinder : public Operator_Extension +{ + friend class Engine_Ext_Cylinder; +public: + Operator_Ext_Cylinder(Operator_Cylinder* op); + ~Operator_Ext_Cylinder(); + + virtual bool BuildExtension(); + + virtual Engine_Extension* CreateEngineExtention(); + + virtual bool IsCylinderCoordsSave() {return true;} + + virtual std::string GetExtensionName() {return std::string("Extension for the Cylinder-Coords Operator");} + +protected: + Operator_Cylinder* m_Op_Cyl; + + bool CC_closedAlpha; + bool CC_R0_included; + + //special EC operator for R0 + FDTD_FLOAT* vv_R0; //calc new voltage from old voltage + FDTD_FLOAT* vi_R0; //calc new voltage from old current + +}; + +#endif // OPERATOR_EXT_CYLINDER_H diff --git a/openEMS.pro b/openEMS.pro index 1016be6..ffd6601 100644 --- a/openEMS.pro +++ b/openEMS.pro @@ -44,14 +44,15 @@ SOURCES += main.cpp \ openems.cpp \ FDTD/engine_multithread.cpp \ FDTD/operator_cylinder.cpp \ - FDTD/engine_cylinder.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/excitation.cpp + FDTD/excitation.cpp \ + FDTD/operator_ext_cylinder.cpp \ + FDTD/engine_ext_cylinder.cpp HEADERS += tools/ErrorMsg.h \ tools/AdrOp.h \ tools/constants.h \ @@ -68,12 +69,13 @@ HEADERS += tools/ErrorMsg.h \ FDTD/operator_cylinder.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/excitation.h + FDTD/excitation.h \ + FDTD/operator_ext_cylinder.h \ + FDTD/engine_ext_cylinder.h QMAKE_CXXFLAGS_RELEASE = -O3 \ -g \ -march=native diff --git a/openems.cpp b/openems.cpp index 21c45db..193de86 100644 --- a/openems.cpp +++ b/openems.cpp @@ -19,7 +19,7 @@ #include #include "tools/array_ops.h" #include "FDTD/engine.h" -#include "FDTD/engine_cylinder.h" +#include "FDTD/operator_cylinder.h" #include "FDTD/engine_multithread.h" #include "FDTD/engine_sse.h" #include "FDTD/operator_ext_mur_abc.h" @@ -255,7 +255,17 @@ int openEMS::SetupFDTD(const char* file) if (CylinderCoords) { cerr << "openEMS: creating cylinder coordinate FDTD engine..." << endl; - FDTD_Eng = Engine_Cylinder::New((Operator_Cylinder*)FDTD_Op); + switch (m_engine) { +// case EngineType_Multithreaded: +// FDTD_Eng = Engine_Multithread::New(FDTD_Op,m_engine_numThreads); +// break; +// case EngineType_SSE: +// FDTD_Eng = Engine_sse::New(dynamic_cast(FDTD_Op)); +// break; + default: + FDTD_Eng = Engine::New(FDTD_Op); + break; + } } else {