From b296c441f9764edcec888579090cfd8d5a983b7e Mon Sep 17 00:00:00 2001 From: Thorsten Liebig Date: Tue, 27 Apr 2010 23:06:42 +0200 Subject: [PATCH] extention updates & new extention: Mur's absorbing boundary condition --- FDTD/engine.cpp | 21 +++++-- FDTD/engine.h | 3 +- FDTD/engine_ext_mur_abc.cpp | 105 ++++++++++++++++++++++++++++++++++ FDTD/engine_ext_mur_abc.h | 50 ++++++++++++++++ FDTD/engine_extension.cpp | 21 ++++++- FDTD/engine_extension.h | 8 ++- FDTD/operator.cpp | 3 + FDTD/operator.h | 2 + FDTD/operator_ext_mur_abc.cpp | 80 ++++++++++++++++++++++++++ FDTD/operator_ext_mur_abc.h | 49 ++++++++++++++++ FDTD/operator_extension.cpp | 17 ++++++ FDTD/operator_extension.h | 3 + openEMS.pro | 8 ++- openems.cpp | 13 +++++ tools/array_ops.cpp | 27 +++++++++ tools/array_ops.h | 3 + 16 files changed, 400 insertions(+), 13 deletions(-) create mode 100644 FDTD/engine_ext_mur_abc.cpp create mode 100644 FDTD/engine_ext_mur_abc.h create mode 100644 FDTD/operator_ext_mur_abc.cpp create mode 100644 FDTD/operator_ext_mur_abc.h diff --git a/FDTD/engine.cpp b/FDTD/engine.cpp index 9a29692..c2f0411 100644 --- a/FDTD/engine.cpp +++ b/FDTD/engine.cpp @@ -17,6 +17,7 @@ #include "engine.h" #include "engine_extension.h" +#include "operator_extension.h" #include "tools/array_ops.h" //! \brief construct an Engine instance @@ -35,10 +36,25 @@ Engine::Engine(const Operator* op) { numLines[n] = Op->GetNumberOfLines(n); } + + for (size_t n=0;nGetNumberOfExtentions();++n) + { + Operator_Extension* op_ext = Op->GetExtension(n); + Engine_Extension* eng_ext = op_ext->CreateEngineExtention(); + if (eng_ext) + { + eng_ext->SetEngine(this); + m_Eng_exts.push_back(eng_ext); + } + } } Engine::~Engine() { + for (size_t n=0;nReset(); } @@ -143,11 +159,6 @@ void Engine::ApplyCurrentExcite() { } -void Engine::AddExtension(Engine_Extension* eng_ext) -{ - m_Eng_exts.push_back(eng_ext); -} - bool Engine::IterateTS(unsigned int iterTS) { for (unsigned int iter=0;iter. +*/ + +#include "engine_ext_mur_abc.h" +#include "operator_ext_mur_abc.h" +#include "engine.h" +#include "tools/array_ops.h" + +Engine_Ext_Mur_ABC::Engine_Ext_Mur_ABC(Operator_Ext_Mur_ABC* op_ext) : Engine_Extension(op_ext) +{ + m_Op_mur = op_ext; + m_numLines[0] = m_Op_mur->m_numLines[0]; + m_numLines[1] = m_Op_mur->m_numLines[1]; + m_ny = m_Op_mur->m_ny; + m_nyP = m_Op_mur->m_nyP; + m_nyPP = m_Op_mur->m_nyPP; + m_LineNr = m_Op_mur->m_LineNr; + m_LineNr_Shift = m_Op_mur->m_LineNr_Shift; + + m_Mur_Coeff = m_Op_mur->m_Mur_Coeff; + + m_volt_nyP = Create2DArray(m_numLines); + m_volt_nyPP = Create2DArray(m_numLines); +} + +Engine_Ext_Mur_ABC::~Engine_Ext_Mur_ABC() +{ + Delete2DArray(m_volt_nyP,m_numLines); + m_volt_nyP = NULL; + Delete2DArray(m_volt_nyPP,m_numLines); + m_volt_nyPP = NULL; +} + +void Engine_Ext_Mur_ABC::DoPreVoltageUpdates() +{ + if (m_Eng==NULL) return; + if (m_Mur_Coeff==0) return; + unsigned int pos[] = {0,0,0}; + unsigned int pos_shift[] = {0,0,0}; + pos[m_ny] = m_LineNr; + pos_shift[m_ny] = m_LineNr_Shift; + + for (pos[m_nyP]=0;pos[m_nyP]GetVolt(m_nyP,pos_shift) - m_Mur_Coeff * m_Eng->GetVolt(m_nyP,pos); + m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]] = m_Eng->GetVolt(m_nyPP,pos_shift) - m_Mur_Coeff * m_Eng->GetVolt(m_nyPP,pos); + } + } +} + +void Engine_Ext_Mur_ABC::DoPostVoltageUpdates() +{ + if (m_Eng==NULL) return; + if (m_Mur_Coeff==0) return; + unsigned int pos[] = {0,0,0}; + unsigned int pos_shift[] = {0,0,0}; + pos[m_ny] = m_LineNr; + pos_shift[m_ny] = m_LineNr_Shift; + + for (pos[m_nyP]=0;pos[m_nyP]GetVolt(m_nyP,pos_shift); + m_volt_nyPP[pos[m_nyP]][pos[m_nyPP]] += m_Mur_Coeff * m_Eng->GetVolt(m_nyPP,pos_shift); + } + } +} + +void Engine_Ext_Mur_ABC::Apply2Voltages() +{ + if (m_Eng==NULL) return; + if (m_Mur_Coeff==0) return; + unsigned int pos[] = {0,0,0}; + pos[m_ny] = m_LineNr; + + 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]]; + } + } +} diff --git a/FDTD/engine_ext_mur_abc.h b/FDTD/engine_ext_mur_abc.h new file mode 100644 index 0000000..98f38cf --- /dev/null +++ b/FDTD/engine_ext_mur_abc.h @@ -0,0 +1,50 @@ +/* +* 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 . +*/ + +#ifndef ENGINE_EXT_MUR_ABC_H +#define ENGINE_EXT_MUR_ABC_H + +#include "engine_extension.h" +#include "operator.h" + +class Operator_Ext_Mur_ABC; + +class Engine_Ext_Mur_ABC : public Engine_Extension +{ +public: + Engine_Ext_Mur_ABC(Operator_Ext_Mur_ABC* op_ext); + virtual ~Engine_Ext_Mur_ABC(); + + virtual void DoPreVoltageUpdates(); + virtual void DoPostVoltageUpdates(); + virtual void Apply2Voltages(); + +protected: + Operator_Ext_Mur_ABC* m_Op_mur; + + int m_ny; + int m_nyP,m_nyPP; + unsigned int m_LineNr; + int m_LineNr_Shift; + unsigned int m_numLines[2]; + + FDTD_FLOAT m_Mur_Coeff; + FDTD_FLOAT** m_volt_nyP; //n+1 direction + FDTD_FLOAT** m_volt_nyPP; //n+2 direction +}; + +#endif // ENGINE_EXT_MUR_ABC_H diff --git a/FDTD/engine_extension.cpp b/FDTD/engine_extension.cpp index 653d99b..ec4b2b9 100644 --- a/FDTD/engine_extension.cpp +++ b/FDTD/engine_extension.cpp @@ -1,9 +1,26 @@ +/* +* 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_extension.h" #include "engine.h" -Engine_Extension::Engine_Extension(Operator_Extension* op_ext, Engine* eng) +Engine_Extension::Engine_Extension(Operator_Extension* op_ext) { m_Op_ext = op_ext; - m_Eng = eng; + m_Eng = NULL; } diff --git a/FDTD/engine_extension.h b/FDTD/engine_extension.h index 756157a..8bbf830 100644 --- a/FDTD/engine_extension.h +++ b/FDTD/engine_extension.h @@ -29,18 +29,20 @@ public: virtual void DoPreVoltageUpdates() {} //! This methode will be called __after__ the main engine does the usual voltage updates. This methode may __not__ change the engine voltages!!! virtual void DoPostVoltageUpdates() {} - //! This methode will be called __after__ all updates to the voltages and extensions and may add its results to the engine voltages, but may __not__ rely on the current value of the engine voltages!!! + //! This methode will be called __after__ all updates to the voltages and extensions and may add/set its results to the engine voltages, but may __not__ rely on the current value of the engine voltages!!! virtual void Apply2Voltages() {} //! This methode will be called __before__ the main engine does the usual current updates. This methode may __not__ change the engine current!!! virtual void DoPreCurrentUpdates() {} //! This methode will be called __after__ the main engine does the usual current updates. This methode may __not__ change the engine current!!! virtual void DoPostCurrentUpdates() {} - //! This methode will be called __after__ all updates to the current and extensions and may add its results to the engine current, but may __not__ rely on the current value of the engine current!!! + //! This methode will be called __after__ all updates to the current and extensions and may add/set its results to the engine current, but may __not__ rely on the current value of the engine current!!! virtual void Apply2Current() {} + //! Set the Engine to this extention. This will usually done automatically by Engine::AddExtension + virtual void SetEngine(Engine* eng) {m_Eng=eng;} protected: - Engine_Extension(Operator_Extension* op_ext, Engine* eng); + Engine_Extension(Operator_Extension* op_ext); Operator_Extension* m_Op_ext; Engine* m_Eng; diff --git a/FDTD/operator.cpp b/FDTD/operator.cpp index c4caa47..80882b5 100644 --- a/FDTD/operator.cpp +++ b/FDTD/operator.cpp @@ -35,6 +35,9 @@ Operator::Operator() Operator::~Operator() { + for (size_t n=0;n. +*/ + +#include "operator_ext_mur_abc.h" +#include "engine_ext_mur_abc.h" + +#include "tools/array_ops.h" + +Operator_Ext_Mur_ABC::Operator_Ext_Mur_ABC(Operator* op) : Operator_Extension(op) +{ + m_Mur_Coeff = NULL; + m_ny = -1; + m_nyP = -1; + m_nyPP = -1; + m_LineNr = 0; + m_LineNr_Shift = 0; + m_Mur_Coeff = 0; +} + +Operator_Ext_Mur_ABC::~Operator_Ext_Mur_ABC() +{ +} + +void Operator_Ext_Mur_ABC::SetDirection(int ny, bool top_ny) +{ + if ((ny<0) || (ny>2)) + return; + m_ny = ny; + m_nyP = (ny+1)%3; + m_nyPP = (ny+2)%3; + if (!top_ny) + { + m_LineNr = 0; + m_LineNr_Shift = 1; + } + else + { + m_LineNr = m_Op->GetNumberOfLines(m_ny)-1; + m_LineNr_Shift = m_Op->GetNumberOfLines(m_ny) - 2; + } + + m_numLines[0] = m_Op->GetNumberOfLines(m_nyP); + m_numLines[1] = m_Op->GetNumberOfLines(m_nyPP); +} + +bool Operator_Ext_Mur_ABC::BuildExtension() +{ + if (m_ny<0) + { + cerr << "Operator_Ext_Mur_ABC::BuildExtension: Warning, Extension not initialized! Use SetDirection!! Abort build!!" << endl; + return false; + } + double dT = m_Op->GetTimestep(); + unsigned int pos[] = {0,0,0}; + pos[m_ny] = m_LineNr; + double delta = fabs(m_Op->GetMeshDelta(m_ny,pos)); + m_Mur_Coeff = (__C0__ * dT - delta) / (__C0__ * dT + delta); +// cerr << "Operator_Ext_Mur_ABC::BuildExtension(): " << m_Mur_Coeff << endl; + return true; +} + +Engine_Extension* Operator_Ext_Mur_ABC::CreateEngineExtention() +{ + Engine_Ext_Mur_ABC* eng_ext = new Engine_Ext_Mur_ABC(this); + return eng_ext; +} diff --git a/FDTD/operator_ext_mur_abc.h b/FDTD/operator_ext_mur_abc.h new file mode 100644 index 0000000..52b15d1 --- /dev/null +++ b/FDTD/operator_ext_mur_abc.h @@ -0,0 +1,49 @@ +/* +* 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 . +*/ + +#ifndef OPERATOR_EXT_MUR_ABC_H +#define OPERATOR_EXT_MUR_ABC_H + +#include "operator.h" +#include "operator_extension.h" + +class Operator_Ext_Mur_ABC : public Operator_Extension +{ + friend class Engine_Ext_Mur_ABC; +public: + Operator_Ext_Mur_ABC(Operator* op); + ~Operator_Ext_Mur_ABC(); + + //! Define the direction of this ABC: ny=0,1,2 -> x,y,z and if at bottom_ny -> e.g. x=0 or x=end + void SetDirection(int ny, bool top_ny); + + virtual bool BuildExtension(); + + virtual Engine_Extension* CreateEngineExtention(); + +protected: + int m_ny; + int m_nyP,m_nyPP; + unsigned int m_LineNr; + int m_LineNr_Shift; + + unsigned int m_numLines[2]; + + FDTD_FLOAT m_Mur_Coeff; +}; + +#endif // OPERATOR_EXT_MUR_ABC_H diff --git a/FDTD/operator_extension.cpp b/FDTD/operator_extension.cpp index 113890a..a94b2a1 100644 --- a/FDTD/operator_extension.cpp +++ b/FDTD/operator_extension.cpp @@ -1,3 +1,20 @@ +/* +* 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 "operator_extension.h" #include "operator.h" diff --git a/FDTD/operator_extension.h b/FDTD/operator_extension.h index 68f0a5e..986ba3b 100644 --- a/FDTD/operator_extension.h +++ b/FDTD/operator_extension.h @@ -19,6 +19,7 @@ #define OPERATOR_EXTENSION_H class Operator; +class Engine_Extension; //! Abstract base-class for all operator extensions class Operator_Extension @@ -26,6 +27,8 @@ class Operator_Extension public: virtual bool BuildExtension() {return true;} + virtual Engine_Extension* CreateEngineExtention() {return 0;} + protected: Operator_Extension(Operator* op); Operator* m_Op; diff --git a/openEMS.pro b/openEMS.pro index 156d065..cb5f43a 100644 --- a/openEMS.pro +++ b/openEMS.pro @@ -40,7 +40,9 @@ SOURCES += main.cpp \ FDTD/operator_cylinder.cpp \ FDTD/engine_cylinder.cpp \ FDTD/operator_extension.cpp \ - FDTD/engine_extension.cpp + FDTD/engine_extension.cpp \ + FDTD/engine_ext_mur_abc.cpp \ + FDTD/operator_ext_mur_abc.cpp HEADERS += tools/ErrorMsg.h \ tools/AdrOp.h \ tools/constants.h \ @@ -58,7 +60,9 @@ HEADERS += tools/ErrorMsg.h \ FDTD/operator_cylinder.h \ FDTD/engine_cylinder.h \ FDTD/operator_extension.h \ - FDTD/engine_extension.h + FDTD/engine_extension.h \ + FDTD/engine_ext_mur_abc.h \ + FDTD/operator_ext_mur_abc.h QMAKE_CXXFLAGS_RELEASE = -O2 \ -g \ -march=native diff --git a/openems.cpp b/openems.cpp index a918e6d..8a0f6ce 100644 --- a/openems.cpp +++ b/openems.cpp @@ -21,6 +21,7 @@ #include "FDTD/engine.h" #include "FDTD/engine_cylinder.h" #include "FDTD/engine_multithread.h" +#include "FDTD/operator_ext_mur_abc.h" #include "FDTD/processvoltage.h" #include "FDTD/processcurrent.h" #include "FDTD/processfields_td.h" @@ -254,6 +255,18 @@ int openEMS::SetupFDTD(const char* file) if (FDTD_Op->SetGeometryCSX(&CSX)==false) return(2); + /**************************** create all operator/engine extensions here !!!! **********************************/ + //Mur-ABC + for (int n=0;n<6;++n) + { + if (bounds[n]==2) + { + Operator_Ext_Mur_ABC* op_ext_mur = new Operator_Ext_Mur_ABC(FDTD_Op); + op_ext_mur->SetDirection(n/2,n%2); + FDTD_Op->AddExtension(op_ext_mur); + } + } + FDTD_Op->CalcECOperator(); SetupExcitation(FDTD_Opts->FirstChildElement("Excitation")); diff --git a/tools/array_ops.cpp b/tools/array_ops.cpp index 0393466..ef08217 100644 --- a/tools/array_ops.cpp +++ b/tools/array_ops.cpp @@ -18,6 +18,33 @@ #include "array_ops.h" #include +FDTD_FLOAT** Create2DArray(const unsigned int* numLines) +{ + FDTD_FLOAT** array=NULL; + unsigned int pos[3]; + array = new FDTD_FLOAT*[numLines[0]]; + for (pos[0]=0;pos[0]