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
This commit is contained in:
Thorsten Liebig 2011-11-07 12:07:55 +01:00
parent 0f27f2c886
commit b34848f323
13 changed files with 179 additions and 38 deletions

View File

@ -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);

View File

@ -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<Engine_Interface_FDTD*>(m_Eng_Interface);
if (EI_FDTD)
{
const Engine* Eng = EI_FDTD->GetFDTDEngine();
unsigned int pos[3];
for (pos[0]=0; pos[0]<Op->GetNumberOfLines(0)-1; ++pos[0])
{
for (pos[1]=0; pos[1]<Op->GetNumberOfLines(1)-1; ++pos[1])
{
for (pos[2]=0; pos[2]<Op->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)

View File

@ -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<float> 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;}

View File

@ -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<Operator_Cylinder*>(op);
if (m_Op_Cyl==NULL)

View File

@ -15,16 +15,16 @@
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
#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;

View File

@ -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]<m_Op->GetNumberOfLines(0)-1; ++pos[0])
{
for (pos[1]=0; pos[1]<m_Op->GetNumberOfLines(1)-1; ++pos[1])
{
for (pos[2]=0; pos[2]<m_Op->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]<m_Op->GetNumberOfLines(0)-1; ++pos[0])
{
for (pos[1]=0; pos[1]<m_Op->GetNumberOfLines(1)-1; ++pos[1])
{
for (pos[2]=0; pos[2]<m_Op->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;
}

View File

@ -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;

View File

@ -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 <http://www.gnu.org/licenses/>.
*/
#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]<m_Op_SSE->GetNumberOfLines(0)-1; ++pos[0])
{
for (pos[1]=0; pos[1]<m_Op_SSE->GetNumberOfLines(1)-1; ++pos[1])
{
for (pos[2]=0; pos[2]<m_Op_SSE->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]);
}

View File

@ -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 <http://www.gnu.org/licenses/>.
*/
#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

View File

@ -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

View File

@ -23,6 +23,7 @@
class Operator_sse : public Operator
{
friend class Engine_Interface_SSE_FDTD;
public:
//! Create a new operator
static Operator_sse* New();

View File

@ -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

View File

@ -282,8 +282,12 @@ bool openEMS::SetupBoundaryConditions(TiXmlElement* BC)
Engine_Interface_FDTD* openEMS::NewEngineInterface()
{
Operator_Cylinder* op_cyl = dynamic_cast<Operator_Cylinder*>(FDTD_Op);
if (op_cyl)
return new Engine_Interface_Cylindrical_FDTD(FDTD_Op,FDTD_Eng);
Engine_sse* eng_sse = dynamic_cast<Engine_sse*>(FDTD_Eng);
if (op_cyl && eng_sse)
return new Engine_Interface_Cylindrical_FDTD(op_cyl,eng_sse);
Operator_sse* op_sse = dynamic_cast<Operator_sse*>(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 << "%)" ;