From b34848f32307d7d2f3018d9d6d41bff02cee6eef Mon Sep 17 00:00:00 2001 From: Thorsten Liebig Date: Mon, 7 Nov 2011 12:07:55 +0100 Subject: [PATCH] new engine interface using sse & speed enhancement in energy estimate Note: cylindrical engine interface now based on sse engine interface --> similar to sse engine & cylindrical engine --- Common/engine_interface_base.h | 8 +++ Common/processfields.cpp | 29 +--------- Common/processfields.h | 2 +- FDTD/engine_interface_cylindrical_fdtd.cpp | 2 +- FDTD/engine_interface_cylindrical_fdtd.h | 6 +- FDTD/engine_interface_fdtd.cpp | 47 ++++++++++++++++ FDTD/engine_interface_fdtd.h | 2 + FDTD/engine_interface_sse_fdtd.cpp | 64 ++++++++++++++++++++++ FDTD/engine_interface_sse_fdtd.h | 38 +++++++++++++ FDTD/openems_fdtd_mpi.cpp | 4 +- FDTD/operator_sse.h | 1 + openEMS.pro | 2 + openems.cpp | 12 ++-- 13 files changed, 179 insertions(+), 38 deletions(-) create mode 100644 FDTD/engine_interface_sse_fdtd.cpp create mode 100644 FDTD/engine_interface_sse_fdtd.h diff --git a/Common/engine_interface_base.h b/Common/engine_interface_base.h index 5f1ad5a..d25c8d4 100644 --- a/Common/engine_interface_base.h +++ b/Common/engine_interface_base.h @@ -67,6 +67,14 @@ public: //! Get the current number of timesteps virtual unsigned int GetNumberOfTimesteps() const =0; + //! Calc (roughly) the total energy + /*! + This method only calculates a very rough estimate of the total energy in the simulation domain. + The result may even be roughly proportional to the real system energy only. + Primary goal is speed, not accuracy. + */ + virtual double CalcFastEnergy() const =0; + protected: Engine_Interface_Base(Operator_Base* base_op); diff --git a/Common/processfields.cpp b/Common/processfields.cpp index 1d3c8c2..4fdcdfe 100644 --- a/Common/processfields.cpp +++ b/Common/processfields.cpp @@ -123,34 +123,9 @@ void ProcessFields::DefineStartStopCoord(double* dstart, double* dstop) } } -double ProcessFields::CalcTotalEnergy() const +double ProcessFields::CalcTotalEnergyEstimate() const { - double energy=0.0; - - Engine_Interface_FDTD* EI_FDTD = dynamic_cast(m_Eng_Interface); - - if (EI_FDTD) - { - const Engine* Eng = EI_FDTD->GetFDTDEngine(); - - unsigned int pos[3]; - for (pos[0]=0; pos[0]GetNumberOfLines(0)-1; ++pos[0]) - { - for (pos[1]=0; pos[1]GetNumberOfLines(1)-1; ++pos[1]) - { - for (pos[2]=0; pos[2]GetNumberOfLines(2)-1; ++pos[2]) - { - energy+=fabs(Eng->GetVolt(0,pos[0],pos[1],pos[2]) * Eng->GetCurr(1,pos[0],pos[1],pos[2])); - energy+=fabs(Eng->GetVolt(0,pos[0],pos[1],pos[2]) * Eng->GetCurr(2,pos[0],pos[1],pos[2])); - energy+=fabs(Eng->GetVolt(1,pos[0],pos[1],pos[2]) * Eng->GetCurr(0,pos[0],pos[1],pos[2])); - energy+=fabs(Eng->GetVolt(1,pos[0],pos[1],pos[2]) * Eng->GetCurr(2,pos[0],pos[1],pos[2])); - energy+=fabs(Eng->GetVolt(2,pos[0],pos[1],pos[2]) * Eng->GetCurr(0,pos[0],pos[1],pos[2])); - energy+=fabs(Eng->GetVolt(2,pos[0],pos[1],pos[2]) * Eng->GetCurr(1,pos[0],pos[1],pos[2])); - } - } - } - } - return energy*0.5; + return m_Eng_Interface->CalcFastEnergy(); } void ProcessFields::SetSubSampling(unsigned int subSampleRate, int dir) diff --git a/Common/processfields.h b/Common/processfields.h index 01575db..11a2c5c 100644 --- a/Common/processfields.h +++ b/Common/processfields.h @@ -79,7 +79,7 @@ public: //! Dump a frequency-domain complex-vector dump to an HDF5 file static bool DumpVectorArray2HDF5(string filename, string groupName, string name, std::complex const* const* const* const* array, unsigned int const* numLines, float weight, float frequency); - double CalcTotalEnergy() const; + double CalcTotalEnergyEstimate() const; void SetFileType(FileType fileType) {m_fileType=fileType;} diff --git a/FDTD/engine_interface_cylindrical_fdtd.cpp b/FDTD/engine_interface_cylindrical_fdtd.cpp index d19df87..a140888 100644 --- a/FDTD/engine_interface_cylindrical_fdtd.cpp +++ b/FDTD/engine_interface_cylindrical_fdtd.cpp @@ -17,7 +17,7 @@ #include "engine_interface_cylindrical_fdtd.h" -Engine_Interface_Cylindrical_FDTD::Engine_Interface_Cylindrical_FDTD(Operator* op, Engine* eng) : Engine_Interface_FDTD(op,eng) +Engine_Interface_Cylindrical_FDTD::Engine_Interface_Cylindrical_FDTD(Operator_sse* op, Engine_sse* eng) : Engine_Interface_SSE_FDTD(op,eng) { m_Op_Cyl = dynamic_cast(op); if (m_Op_Cyl==NULL) diff --git a/FDTD/engine_interface_cylindrical_fdtd.h b/FDTD/engine_interface_cylindrical_fdtd.h index 5df12eb..941a8ba 100644 --- a/FDTD/engine_interface_cylindrical_fdtd.h +++ b/FDTD/engine_interface_cylindrical_fdtd.h @@ -15,16 +15,16 @@ * along with this program. If not, see . */ -#include "engine_interface_fdtd.h" +#include "engine_interface_sse_fdtd.h" #include "operator_cylinder.h" #ifndef ENGINE_INTERFACE_CYLINDRICAL_FDTD_H #define ENGINE_INTERFACE_CYLINDRICAL_FDTD_H -class Engine_Interface_Cylindrical_FDTD : public Engine_Interface_FDTD +class Engine_Interface_Cylindrical_FDTD : public Engine_Interface_SSE_FDTD { public: - Engine_Interface_Cylindrical_FDTD(Operator* op, Engine* eng); + Engine_Interface_Cylindrical_FDTD(Operator_sse* op, Engine_sse* eng); virtual ~Engine_Interface_Cylindrical_FDTD(); virtual double* GetHField(const unsigned int* pos, double* out) const; diff --git a/FDTD/engine_interface_fdtd.cpp b/FDTD/engine_interface_fdtd.cpp index ee7e9e0..3768b5b 100644 --- a/FDTD/engine_interface_fdtd.cpp +++ b/FDTD/engine_interface_fdtd.cpp @@ -211,3 +211,50 @@ double Engine_Interface_FDTD::GetRawField(unsigned int n, const unsigned int* po return 0.0; } + +double Engine_Interface_FDTD::CalcFastEnergy() const +{ + double E_energy=0.0; + double H_energy=0.0; + + unsigned int pos[3]; + if (m_Eng->GetType()==Engine::BASIC) + { + for (pos[0]=0; pos[0]GetNumberOfLines(0)-1; ++pos[0]) + { + for (pos[1]=0; pos[1]GetNumberOfLines(1)-1; ++pos[1]) + { + for (pos[2]=0; pos[2]GetNumberOfLines(2)-1; ++pos[2]) + { + E_energy+=m_Eng->Engine::GetVolt(0,pos[0],pos[1],pos[2]) * m_Eng->Engine::GetVolt(0,pos[0],pos[1],pos[2]); + E_energy+=m_Eng->Engine::GetVolt(1,pos[0],pos[1],pos[2]) * m_Eng->Engine::GetVolt(1,pos[0],pos[1],pos[2]); + E_energy+=m_Eng->Engine::GetVolt(2,pos[0],pos[1],pos[2]) * m_Eng->Engine::GetVolt(2,pos[0],pos[1],pos[2]); + + H_energy+=m_Eng->Engine::GetCurr(0,pos[0],pos[1],pos[2]) * m_Eng->Engine::GetCurr(0,pos[0],pos[1],pos[2]); + H_energy+=m_Eng->Engine::GetCurr(1,pos[0],pos[1],pos[2]) * m_Eng->Engine::GetCurr(1,pos[0],pos[1],pos[2]); + H_energy+=m_Eng->Engine::GetCurr(2,pos[0],pos[1],pos[2]) * m_Eng->Engine::GetCurr(2,pos[0],pos[1],pos[2]); + } + } + } + } + else + { + for (pos[0]=0; pos[0]GetNumberOfLines(0)-1; ++pos[0]) + { + for (pos[1]=0; pos[1]GetNumberOfLines(1)-1; ++pos[1]) + { + for (pos[2]=0; pos[2]GetNumberOfLines(2)-1; ++pos[2]) + { + E_energy+=m_Eng->GetVolt(0,pos[0],pos[1],pos[2]) * m_Eng->GetVolt(0,pos[0],pos[1],pos[2]); + E_energy+=m_Eng->GetVolt(1,pos[0],pos[1],pos[2]) * m_Eng->GetVolt(1,pos[0],pos[1],pos[2]); + E_energy+=m_Eng->GetVolt(2,pos[0],pos[1],pos[2]) * m_Eng->GetVolt(2,pos[0],pos[1],pos[2]); + + H_energy+=m_Eng->GetCurr(0,pos[0],pos[1],pos[2]) * m_Eng->GetCurr(0,pos[0],pos[1],pos[2]); + H_energy+=m_Eng->GetCurr(1,pos[0],pos[1],pos[2]) * m_Eng->GetCurr(1,pos[0],pos[1],pos[2]); + H_energy+=m_Eng->GetCurr(2,pos[0],pos[1],pos[2]) * m_Eng->GetCurr(2,pos[0],pos[1],pos[2]); + } + } + } + } + return __EPS0__*E_energy + __MUE0__*H_energy; +} diff --git a/FDTD/engine_interface_fdtd.h b/FDTD/engine_interface_fdtd.h index 8ef20bd..e93f7fd 100644 --- a/FDTD/engine_interface_fdtd.h +++ b/FDTD/engine_interface_fdtd.h @@ -50,6 +50,8 @@ public: virtual double GetTime(bool dualTime=false) const {return ((double)m_Eng->GetNumberOfTimesteps() + (double)dualTime*0.5)*m_Op->GetTimestep();}; virtual unsigned int GetNumberOfTimesteps() const {return m_Eng->GetNumberOfTimesteps();} + virtual double CalcFastEnergy() const; + protected: Operator* m_Op; Engine* m_Eng; diff --git a/FDTD/engine_interface_sse_fdtd.cpp b/FDTD/engine_interface_sse_fdtd.cpp new file mode 100644 index 0000000..8f23547 --- /dev/null +++ b/FDTD/engine_interface_sse_fdtd.cpp @@ -0,0 +1,64 @@ +/* +* Copyright (C) 2011 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_interface_sse_fdtd.h" + +Engine_Interface_SSE_FDTD::Engine_Interface_SSE_FDTD(Operator_sse* op, Engine_sse* eng) : Engine_Interface_FDTD(op, eng) +{ + m_Op_SSE = op; + m_Eng_SSE = eng; +} + +Engine_Interface_SSE_FDTD::~Engine_Interface_SSE_FDTD() +{ + m_Op_SSE=NULL; + m_Eng_SSE=NULL; +} + +double Engine_Interface_SSE_FDTD::CalcFastEnergy() const +{ + f4vector E_energy; + E_energy.f[0]=0; + E_energy.f[1]=0; + E_energy.f[2]=0; + E_energy.f[3]=0; + f4vector H_energy; + H_energy = E_energy; + + if (m_Eng_SSE->GetType()!=Engine::SSE) + return Engine_Interface_FDTD::CalcFastEnergy(); + + unsigned int pos[3]; + for (pos[0]=0; pos[0]GetNumberOfLines(0)-1; ++pos[0]) + { + for (pos[1]=0; pos[1]GetNumberOfLines(1)-1; ++pos[1]) + { + for (pos[2]=0; pos[2]numVectors; ++pos[2]) + { + E_energy.v += m_Eng_SSE->Engine_sse::f4_volt[0][pos[0]][pos[1]][pos[2]].v * m_Eng_SSE->Engine_sse::f4_volt[0][pos[0]][pos[1]][pos[2]].v; + E_energy.v += m_Eng_SSE->Engine_sse::f4_volt[1][pos[0]][pos[1]][pos[2]].v * m_Eng_SSE->Engine_sse::f4_volt[1][pos[0]][pos[1]][pos[2]].v; + E_energy.v += m_Eng_SSE->Engine_sse::f4_volt[2][pos[0]][pos[1]][pos[2]].v * m_Eng_SSE->Engine_sse::f4_volt[2][pos[0]][pos[1]][pos[2]].v; + + H_energy.v += m_Eng_SSE->Engine_sse::f4_curr[0][pos[0]][pos[1]][pos[2]].v * m_Eng_SSE->Engine_sse::f4_curr[0][pos[0]][pos[1]][pos[2]].v; + H_energy.v += m_Eng_SSE->Engine_sse::f4_curr[1][pos[0]][pos[1]][pos[2]].v * m_Eng_SSE->Engine_sse::f4_curr[1][pos[0]][pos[1]][pos[2]].v; + H_energy.v += m_Eng_SSE->Engine_sse::f4_curr[2][pos[0]][pos[1]][pos[2]].v * m_Eng_SSE->Engine_sse::f4_curr[2][pos[0]][pos[1]][pos[2]].v; + } + } + } + + return __EPS0__*(E_energy.f[0]+E_energy.f[1]+E_energy.f[2]+E_energy.f[3]) + __MUE0__*(H_energy.f[0]+H_energy.f[1]+H_energy.f[2]+H_energy.f[3]); +} diff --git a/FDTD/engine_interface_sse_fdtd.h b/FDTD/engine_interface_sse_fdtd.h new file mode 100644 index 0000000..17327e1 --- /dev/null +++ b/FDTD/engine_interface_sse_fdtd.h @@ -0,0 +1,38 @@ +/* +* Copyright (C) 2011 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_INTERFACE_SSE_FDTD_H +#define ENGINE_INTERFACE_SSE_FDTD_H + +#include "engine_interface_fdtd.h" +#include "operator_sse.h" +#include "engine_sse.h" + +class Engine_Interface_SSE_FDTD : public Engine_Interface_FDTD +{ +public: + Engine_Interface_SSE_FDTD(Operator_sse* op, Engine_sse* eng); + virtual ~Engine_Interface_SSE_FDTD(); + + virtual double CalcFastEnergy() const; + +protected: + Operator_sse* m_Op_SSE; + Engine_sse* m_Eng_SSE; +}; + +#endif // ENGINE_INTERFACE_SSE_FDTD_H diff --git a/FDTD/openems_fdtd_mpi.cpp b/FDTD/openems_fdtd_mpi.cpp index bf97c45..7ad24ca 100644 --- a/FDTD/openems_fdtd_mpi.cpp +++ b/FDTD/openems_fdtd_mpi.cpp @@ -346,7 +346,7 @@ bool openEMS_FDTD_MPI::CheckEnergyCalc() double openEMS_FDTD_MPI::CalcEnergy() { double energy = 0; - double loc_energy= m_ProcField->CalcTotalEnergy(); + double loc_energy= m_ProcField->CalcTotalEnergyEstimate(); //calc the sum of all local energies MPI_Reduce(&loc_energy, &energy, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); @@ -429,7 +429,7 @@ void openEMS_FDTD_MPI::RunFDTD() MPI_Bcast(&m_NumberCells, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD); //special handling of a field processing, needed to realize the end criteria... - m_ProcField = new ProcessFields(new Engine_Interface_FDTD(FDTD_Op,FDTD_Eng)); + m_ProcField = new ProcessFields(NewEngineInterface()); PA->AddProcessing(m_ProcField); //init processings diff --git a/FDTD/operator_sse.h b/FDTD/operator_sse.h index 280fa43..317ea72 100644 --- a/FDTD/operator_sse.h +++ b/FDTD/operator_sse.h @@ -23,6 +23,7 @@ class Operator_sse : public Operator { + friend class Engine_Interface_SSE_FDTD; public: //! Create a new operator static Operator_sse* New(); diff --git a/openEMS.pro b/openEMS.pro index 605cf53..df76004 100644 --- a/openEMS.pro +++ b/openEMS.pro @@ -82,6 +82,7 @@ SOURCES += FDTD/engine.cpp \ FDTD/operator_cylindermultigrid.cpp \ FDTD/engine_cylindermultigrid.cpp \ FDTD/engine_interface_fdtd.cpp \ + FDTD/engine_interface_sse_fdtd.cpp \ FDTD/engine_interface_cylindrical_fdtd.cpp # FDTD/extensions source files @@ -143,6 +144,7 @@ HEADERS += FDTD/engine.h \ FDTD/operator_cylindermultigrid.h \ FDTD/engine_cylindermultigrid.h \ FDTD/engine_interface_fdtd.h \ + FDTD/engine_interface_sse_fdtd.h \ FDTD/engine_interface_cylindrical_fdtd.h # FDTD/extensions header files diff --git a/openems.cpp b/openems.cpp index e49bde1..d84800c 100644 --- a/openems.cpp +++ b/openems.cpp @@ -282,8 +282,12 @@ bool openEMS::SetupBoundaryConditions(TiXmlElement* BC) Engine_Interface_FDTD* openEMS::NewEngineInterface() { Operator_Cylinder* op_cyl = dynamic_cast(FDTD_Op); - if (op_cyl) - return new Engine_Interface_Cylindrical_FDTD(FDTD_Op,FDTD_Eng); + Engine_sse* eng_sse = dynamic_cast(FDTD_Eng); + if (op_cyl && eng_sse) + return new Engine_Interface_Cylindrical_FDTD(op_cyl,eng_sse); + Operator_sse* op_sse = dynamic_cast(FDTD_Op); + if (op_sse && eng_sse) + return new Engine_Interface_SSE_FDTD(op_sse,eng_sse); return new Engine_Interface_FDTD(FDTD_Op,FDTD_Eng); } @@ -761,7 +765,7 @@ void openEMS::RunFDTD() if (ProcField->CheckTimestep()) { - currE = ProcField->CalcTotalEnergy(); + currE = ProcField->CalcTotalEnergyEstimate(); if (currE>maxE) maxE=currE; } @@ -775,7 +779,7 @@ void openEMS::RunFDTD() t_diff = CalcDiffTime(currTime,prevTime); if (t_diff>4) { - currE = ProcField->CalcTotalEnergy(); + currE = ProcField->CalcTotalEnergyEstimate(); if (currE>maxE) maxE=currE; cout << "[@" << FormatTime(CalcDiffTime(currTime,startTime)) << "] Timestep: " << setw(12) << currTS << " (" << setw(6) << setprecision(2) << std::fixed << (double)currTS/(double)NrTS*100.0 << "%)" ;