zcy
/
openEMS
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Update to vtk-field dumps + material-dump for debugging

pull/1/head
Thorsten Liebig 2010-03-12 20:39:04 +01:00
parent 1220aef54f
commit b9a3165aac
10 changed files with 171 additions and 23 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

54
FDTD/operator.cpp
View File

@ -17,6 +17,7 @@
#include <fstream> #include <fstream>
#include "operator.h" #include "operator.h"
#include "processfields.h"
#include "tools/array_ops.h" #include "tools/array_ops.h"
Operator::Operator() Operator::Operator()
@ -202,6 +203,59 @@ void Operator::DumpOperator2File(string filename)
file.close(); file.close();
} }
void Operator::DumpMaterial2File(string filename)
{
FDTD_FLOAT*** epsilon;
FDTD_FLOAT*** mue;
FDTD_FLOAT*** kappa;
FDTD_FLOAT*** sigma;
unsigned int pos[3];
double inMat[4];
ofstream file(filename.c_str(),ios_base::out);
// file.open;
if (file.is_open()==false)
{
cerr << "Operator::DumpMaterial2File: Can't open file: " << filename << endl;
return;
}
epsilon = Create3DArray( numLines);
mue = Create3DArray( numLines);
kappa = Create3DArray( numLines);
sigma = Create3DArray( numLines);
for (pos[0]=0;pos[0]<numLines[0];++pos[0])
{
for (pos[1]=0;pos[1]<numLines[1];++pos[1])
{
for (pos[2]=0;pos[2]<numLines[2];++pos[2])
{
for (int n=0;n<3;++n)
{
Calc_EffMatPos(n, pos, inMat);
epsilon[pos[0]][pos[1]][pos[2]]+=inMat[0]/__EPS0__;
mue[pos[0]][pos[1]][pos[2]]+=inMat[2]/__MUE0__;
kappa[pos[0]][pos[1]][pos[2]]+=inMat[1];
sigma[pos[0]][pos[1]][pos[2]]+=inMat[3];
}
epsilon[pos[0]][pos[1]][pos[2]]/=3;
mue[pos[0]][pos[1]][pos[2]]/=3;
kappa[pos[0]][pos[1]][pos[2]]/=3;
sigma[pos[0]][pos[1]][pos[2]]/=3;
}
}
}
string names[] = {"epsilon","mue","kappa","sigma"};
FDTD_FLOAT*** array[] = {epsilon,mue,kappa,sigma};
ProcessFields::DumpMultiScalarArray2VTK(file, names, array, 4, discLines, numLines);
Delete3DArray(epsilon,numLines);
Delete3DArray(mue,numLines);
Delete3DArray(kappa,numLines);
Delete3DArray(sigma,numLines);
file.close();
}
bool Operator::SetGeometryCSX(ContinuousStructure* geo) bool Operator::SetGeometryCSX(ContinuousStructure* geo)
{ {
if (geo==NULL) return false; if (geo==NULL) return false;

1
FDTD/operator.h
View File

@ -53,6 +53,7 @@ public:
void ShowSize(); void ShowSize();
void DumpOperator2File(string filename); void DumpOperator2File(string filename);
void DumpMaterial2File(string filename);
virtual void Reset(); virtual void Reset();

63
FDTD/processfields.cpp
View File

@ -109,7 +109,7 @@ void ProcessFields::DefineStartStopCoord(double* dstart, double* dstop)
// else subSample[dir]=subSampleRate; // else subSample[dir]=subSampleRate;
//} //}
bool ProcessFields::DumpFieldArray2VTK(ofstream &file, string name, FDTD_FLOAT**** array, double** discLines, unsigned int* numLines) void ProcessFields::WriteVTKHeader(ofstream &file, double** discLines, unsigned int* numLines)
{ {
file << "# vtk DataFile Version 2.0" << endl; file << "# vtk DataFile Version 2.0" << endl;
file << "Rectilinear Grid openEMS_ProcessFields" << endl; file << "Rectilinear Grid openEMS_ProcessFields" << endl;
@ -128,8 +128,11 @@ bool ProcessFields::DumpFieldArray2VTK(ofstream &file, string name, FDTD_FLOAT**
for (unsigned int i=0;i<numLines[2];++i) for (unsigned int i=0;i<numLines[2];++i)
file << discLines[2][i] << " "; file << discLines[2][i] << " ";
file << endl << endl; file << endl << endl;
file << "POINT_DATA " << numLines[0]*numLines[1]*numLines[2] << endl; file << "POINT_DATA " << numLines[0]*numLines[1]*numLines[2] << endl;
}
void ProcessFields::WriteVTKVectorArray(ofstream &file, string name, FDTD_FLOAT**** array, unsigned int* numLines)
{
file << "VECTORS " << name << " float " << endl; file << "VECTORS " << name << " float " << endl;
unsigned int pos[3]; unsigned int pos[3];
@ -149,3 +152,59 @@ bool ProcessFields::DumpFieldArray2VTK(ofstream &file, string name, FDTD_FLOAT**
} }
} }
} }
bool ProcessFields::DumpVectorArray2VTK(ofstream &file, string name, FDTD_FLOAT**** array, double** discLines, unsigned int* numLines)
{
WriteVTKHeader(file, discLines, numLines);
WriteVTKVectorArray(file, name, array, numLines);
}
bool ProcessFields::DumpMultiVectorArray2VTK(ofstream &file, string names[], FDTD_FLOAT**** array[], unsigned int numFields, double** discLines, unsigned int* numLines)
{
WriteVTKHeader(file, discLines, numLines);
for (int n=0;n<numFields;++n)
{
WriteVTKVectorArray(file, names[n], array[n], numLines);
file << endl;
}
}
void ProcessFields::WriteVTKScalarArray(ofstream &file, string name, FDTD_FLOAT*** array, unsigned int* numLines)
{
file << "SCALARS " << name << " float " << 1 << endl;
file << "LOOKUP_TABLE default" << endl;
unsigned int pos[3];
int count=0;
for (pos[2]=0;pos[2]<numLines[2];++pos[2])
{
for (pos[1]=0;pos[1]<numLines[1];++pos[1])
{
for (pos[0]=0;pos[0]<numLines[0];++pos[0])
{
file << array[pos[0]][pos[1]][pos[2]] << " ";
++count;
if (count%10==0)
file << endl;
}
}
}
}
bool ProcessFields::DumpScalarArray2VTK(ofstream &file, string name, FDTD_FLOAT*** array, double** discLines, unsigned int* numLines)
{
WriteVTKHeader(file, discLines, numLines);
WriteVTKScalarArray(file, name, array, numLines);
}
bool ProcessFields::DumpMultiScalarArray2VTK(ofstream &file, string names[], FDTD_FLOAT*** array[], unsigned int numFields, double** discLines, unsigned int* numLines)
{
WriteVTKHeader(file, discLines, numLines);
for (int n=0;n<numFields;++n)
{
WriteVTKScalarArray(file, names[n], array[n], numLines);
file << endl;
}
}

9
FDTD/processfields.h
View File

@ -41,11 +41,18 @@ public:
//! Set dump type: 0 for E-fields, 1 for H-fields, 2 for D-fields, 3 for B-fields, 4 for J-fields, etc... //! Set dump type: 0 for E-fields, 1 for H-fields, 2 for D-fields, 3 for B-fields, 4 for J-fields, etc...
void SetDumpType(int type) {DumpType=type;} void SetDumpType(int type) {DumpType=type;}
static bool DumpVectorArray2VTK(ofstream &file, string name, FDTD_FLOAT**** array, double** discLines, unsigned int* numLines);
static bool DumpMultiVectorArray2VTK(ofstream &file, string names[], FDTD_FLOAT**** array[], unsigned int numFields, double** discLines, unsigned int* numLines);
static bool DumpScalarArray2VTK(ofstream &file, string name, FDTD_FLOAT*** array, double** discLines, unsigned int* numLines);
static bool DumpMultiScalarArray2VTK(ofstream &file, string names[], FDTD_FLOAT*** array[], unsigned int numFields, double** discLines, unsigned int* numLines);
// virtual void Process(); // virtual void Process();
protected: protected:
ProcessFields(Operator* op, Engine* eng); ProcessFields(Operator* op, Engine* eng);
static bool DumpFieldArray2VTK(ofstream &file, string name, FDTD_FLOAT**** array, double** discLines, unsigned int* numLines); static void WriteVTKHeader(ofstream &file, double** discLines, unsigned int* numLines);
static void WriteVTKVectorArray(ofstream &file, string name, FDTD_FLOAT**** array, unsigned int* numLines);
static void WriteVTKScalarArray(ofstream &file, string name, FDTD_FLOAT*** array, unsigned int* numLines);
int DumpMode; int DumpMode;
int DumpType; int DumpType;

8
FDTD/processfields_td.cpp
View File

@ -63,7 +63,7 @@ void ProcessFieldsTD::DumpCellInterpol(ofstream &file)
} }
} }
} }
DumpFieldArray2VTK(file,string("E-Field"),E_T,discDLines,numDLines); DumpVectorArray2VTK(file,string("E-Field"),E_T,discDLines,numDLines);
Delete_N_3DArray(E_T,numDLines); Delete_N_3DArray(E_T,numDLines);
E_T = NULL; E_T = NULL;
} }
@ -101,7 +101,7 @@ void ProcessFieldsTD::DumpCellInterpol(ofstream &file)
} }
} }
} }
DumpFieldArray2VTK(file,string("H-Field"),H_T,discDLines,numDLines); DumpVectorArray2VTK(file,string("H-Field"),H_T,discDLines,numDLines);
Delete_N_3DArray(H_T,numDLines); Delete_N_3DArray(H_T,numDLines);
H_T = NULL; H_T = NULL;
} }
@ -126,7 +126,7 @@ void ProcessFieldsTD::DumpNoInterpol(ofstream &file)
} }
} }
} }
DumpFieldArray2VTK(file,string("E-Field"),E_T,discLines,numLines); DumpVectorArray2VTK(file,string("E-Field"),E_T,discLines,numLines);
Delete_N_3DArray(E_T,numLines); Delete_N_3DArray(E_T,numLines);
E_T = NULL; E_T = NULL;
} }
@ -149,7 +149,7 @@ void ProcessFieldsTD::DumpNoInterpol(ofstream &file)
} }
} }
} }
DumpFieldArray2VTK(file,string("H-Field"),H_T,discLines,numLines); DumpVectorArray2VTK(file,string("H-Field"),H_T,discLines,numLines);
Delete_N_3DArray(H_T,numLines); Delete_N_3DArray(H_T,numLines);
H_T = NULL; H_T = NULL;
} }

8
examples/FDTD_examples.cpp
View File

@ -285,10 +285,10 @@ void BuildMSL(const char* filename)
//MSL //MSL
CSProperties* MSL = NULL; CSProperties* MSL = NULL;
// CSPropMaterial* MSL_mat = new CSPropMaterial(CSX.GetParameterSet()); CSPropMaterial* MSL_mat = new CSPropMaterial(CSX.GetParameterSet());
// MSL_mat->SetKappa(56e6); MSL_mat->SetKappa(56e6);
// MSL = MSL_mat; MSL = MSL_mat;
MSL = new CSPropMetal(CSX.GetParameterSet()); // MSL = new CSPropMetal(CSX.GetParameterSet());
CSX.AddProperty(MSL); CSX.AddProperty(MSL);
box = new CSPrimBox(CSX.GetParameterSet(),MSL); box = new CSPrimBox(CSX.GetParameterSet(),MSL);

10
examples/MSL.xml
View File

@ -21,10 +21,16 @@
<WeightY Epsilon="1.000000e+00" Mue="1.000000e+00" Kappa="term:pow(abs(z)-800,4)" Sigma="term:pow(abs(z)-800,4)" /> <WeightY Epsilon="1.000000e+00" Mue="1.000000e+00" Kappa="term:pow(abs(z)-800,4)" Sigma="term:pow(abs(z)-800,4)" />
<WeightZ Epsilon="1.000000e+00" Mue="1.000000e+00" Kappa="term:pow(abs(z)-800,4)" Sigma="term:pow(abs(z)-800,4)" /> <WeightZ Epsilon="1.000000e+00" Mue="1.000000e+00" Kappa="term:pow(abs(z)-800,4)" Sigma="term:pow(abs(z)-800,4)" />
</Material> </Material>
<Metal ID="1" Name=""> <Material ID="1" Name="" Isotropy="1">
<FillColor R="118" G="90" B="46" a="255" /> <FillColor R="118" G="90" B="46" a="255" />
<EdgeColor R="118" G="90" B="46" a="255" /> <EdgeColor R="118" G="90" B="46" a="255" />
</Metal> <Property Epsilon="1.000000e+00" Mue="1.000000e+00" Kappa="5.600000e+07" Sigma="0.000000e+00" />
<PropertyY Epsilon="1.000000e+00" Mue="1.000000e+00" Kappa="0.000000e+00" Sigma="0.000000e+00" />
<PropertyZ Epsilon="1.000000e+00" Mue="1.000000e+00" Kappa="0.000000e+00" Sigma="0.000000e+00" />
<WeightX Epsilon="1.000000e+00" Mue="1.000000e+00" Kappa="1.000000e+00" Sigma="1.000000e+00" />
<WeightY Epsilon="1.000000e+00" Mue="1.000000e+00" Kappa="1.000000e+00" Sigma="1.000000e+00" />
<WeightZ Epsilon="1.000000e+00" Mue="1.000000e+00" Kappa="1.000000e+00" Sigma="1.000000e+00" />
</Material>
<Electrode ID="2" Name="" Number="0" Delay="0.000000e+00"> <Electrode ID="2" Name="" Number="0" Delay="0.000000e+00">
<FillColor R="99" G="51" B="159" a="255" /> <FillColor R="99" G="51" B="159" a="255" />
<EdgeColor R="99" G="51" B="159" a="255" /> <EdgeColor R="99" G="51" B="159" a="255" />

14
main.cpp
View File

@ -34,7 +34,7 @@ int main(int argc, char *argv[])
#ifdef STANDALONE #ifdef STANDALONE
if (argc<=1) if (argc<=1)
{ {
cerr << " usage: openEMS FDTD_XML_FILE [--disable-dumps]" << endl; cout << " usage: openEMS FDTD_XML_FILE [--disable-dumps] [--debug-material]" << endl;
exit(-1); exit(-1);
} }
@ -43,7 +43,17 @@ int main(int argc, char *argv[])
for (int n=2;n<argc;++n) for (int n=2;n<argc;++n)
{ {
if (strcmp(argv[n],"--disable-dumps")==0) if (strcmp(argv[n],"--disable-dumps")==0)
{
cout << "openEMS - disabling all field dumps" << endl;
FDTD.SetEnableDumps(false); FDTD.SetEnableDumps(false);
}
else if (strcmp(argv[n],"--debug-material")==0)
{
cout << "openEMS - dumping material to 'material_dump.vtk'" << endl;
FDTD.DebugMaterial();
}
else
cout << "openEMS - unknown argument: " << argv[n] << endl;
} }
} }
@ -63,7 +73,7 @@ int main(int argc, char *argv[])
const char* fileCoax="../examples/Coax_Cart.xml"; const char* fileCoax="../examples/Coax_Cart.xml";
BuildCoaxial_Cartesian(fileCoax); BuildCoaxial_Cartesian(fileCoax);
const char* file=fileCoax; const char* file=fileMSL;
// cerr << CSX.ReadFromXML("examples/PlaneWave.xml") << endl; // cerr << CSX.ReadFromXML("examples/PlaneWave.xml") << endl;
#endif #endif

24
openems.cpp
View File

@ -33,6 +33,7 @@ openEMS::openEMS()
FDTD_Eng=NULL; FDTD_Eng=NULL;
PA=NULL; PA=NULL;
Enable_Dumps = true; Enable_Dumps = true;
DebugMat = false;
} }
openEMS::~openEMS() openEMS::~openEMS()
@ -74,7 +75,7 @@ int openEMS::SetupFDTD(const char* file)
exit(-1); exit(-1);
} }
cerr << "Read openEMS Settings..." << endl; cout << "Read openEMS Settings..." << endl;
TiXmlElement* FDTD_Opts = doc.FirstChildElement("openEMS-Parameter"); TiXmlElement* FDTD_Opts = doc.FirstChildElement("openEMS-Parameter");
if (FDTD_Opts==NULL) if (FDTD_Opts==NULL)
{ {
@ -109,7 +110,7 @@ int openEMS::SetupFDTD(const char* file)
BC->QueryIntAttribute("zmin",&bounds[4]); BC->QueryIntAttribute("zmin",&bounds[4]);
BC->QueryIntAttribute("zmax",&bounds[5]); BC->QueryIntAttribute("zmax",&bounds[5]);
cerr << "Read Geometry..." << endl; cout << "Read Geometry..." << endl;
ContinuousStructure CSX; ContinuousStructure CSX;
string EC(CSX.ReadFromXML(&doc)); string EC(CSX.ReadFromXML(&doc));
if (EC.empty()==false) if (EC.empty()==false)
@ -123,10 +124,15 @@ int openEMS::SetupFDTD(const char* file)
PMC[n]=(bounds[n]==1); PMC[n]=(bounds[n]==1);
//*************** setup operator ************// //*************** setup operator ************//
cerr << "Create Operator..." << endl; cout << "Create Operator..." << endl;
FDTD_Op = new Operator(); FDTD_Op = new Operator();
if (FDTD_Op->SetGeometryCSX(&CSX)==false) return(-1); if (FDTD_Op->SetGeometryCSX(&CSX)==false) return(-1);
if (DebugMat)
{
FDTD_Op->DumpMaterial2File("material_dump.vtk");
}
FDTD_Op->CalcECOperator(); FDTD_Op->CalcECOperator();
if (Excit_Type==0) if (Excit_Type==0)
@ -143,10 +149,10 @@ int openEMS::SetupFDTD(const char* file)
FDTD_Op->ApplyMagneticBC(PMC); FDTD_Op->ApplyMagneticBC(PMC);
cerr << "Nyquist number of timesteps: " << FDTD_Op->GetNyquistNum(f0+fc) << endl; cout << "Nyquist number of timesteps: " << FDTD_Op->GetNyquistNum(f0+fc) << endl;
unsigned int Nyquist = FDTD_Op->GetNyquistNum(f0+fc); unsigned int Nyquist = FDTD_Op->GetNyquistNum(f0+fc);
cerr << "Time for operator: " << difftime(OpDoneTime,startTime) << endl; cout << "Creation time for operator: " << difftime(OpDoneTime,startTime) << " s" << endl;
//create FDTD engine //create FDTD engine
FDTD_Eng = new Engine(FDTD_Op); FDTD_Eng = new Engine(FDTD_Op);
@ -154,6 +160,7 @@ int openEMS::SetupFDTD(const char* file)
time_t currTime = time(NULL); time_t currTime = time(NULL);
//*************** setup processing ************// //*************** setup processing ************//
cout << "Setting up processing..." << endl;
PA = new ProcessingArray(); PA = new ProcessingArray();
double start[3]; double start[3];
@ -229,6 +236,7 @@ int openEMS::SetupFDTD(const char* file)
void openEMS::RunFDTD() void openEMS::RunFDTD()
{ {
cout << "Running FDTD engine... this may take a while... grab a coup of coffee?!?" << endl;
time_t currTime = time(NULL); time_t currTime = time(NULL);
//*************** simulate ************// //*************** simulate ************//
int step=PA->Process(); int step=PA->Process();
@ -237,7 +245,7 @@ void openEMS::RunFDTD()
{ {
FDTD_Eng->IterateTS(step); FDTD_Eng->IterateTS(step);
step=PA->Process(); step=PA->Process();
// cerr << " do " << step << " steps; current: " << eng.GetNumberOfTimesteps() << endl; // cout << " do " << step << " steps; current: " << eng.GetNumberOfTimesteps() << endl;
if ((step<0) || (step>NrTS - FDTD_Eng->GetNumberOfTimesteps())) step=NrTS - FDTD_Eng->GetNumberOfTimesteps(); if ((step<0) || (step>NrTS - FDTD_Eng->GetNumberOfTimesteps())) step=NrTS - FDTD_Eng->GetNumberOfTimesteps();
} }
@ -247,6 +255,6 @@ void openEMS::RunFDTD()
double t_diff = difftime(currTime,prevTime); double t_diff = difftime(currTime,prevTime);
cerr << "Time for " << FDTD_Eng->GetNumberOfTimesteps() << " iterations with " << FDTD_Op->GetNumberCells() << " cells : " << t_diff << " sec" << endl; cout << "Time for " << FDTD_Eng->GetNumberOfTimesteps() << " iterations with " << FDTD_Op->GetNumberCells() << " cells : " << t_diff << " sec" << endl;
cerr << "Speed: " << (double)FDTD_Op->GetNumberCells()*(double)FDTD_Eng->GetNumberOfTimesteps()/t_diff/1e6 << " MCells/s " << endl; cout << "Speed: " << (double)FDTD_Op->GetNumberCells()*(double)FDTD_Eng->GetNumberOfTimesteps()/t_diff/1e6 << " MCells/s " << endl;
} }

3
openems.h
View File

@ -36,10 +36,13 @@ public:
void SetEnableDumps(bool val) {Enable_Dumps=val;} void SetEnableDumps(bool val) {Enable_Dumps=val;}
void DebugMaterial() {DebugMat=true;}
protected: protected:
//! Number of Timesteps //! Number of Timesteps
int NrTS; int NrTS;
bool Enable_Dumps; bool Enable_Dumps;
bool DebugMat;
Operator* FDTD_Op; Operator* FDTD_Op;
Engine* FDTD_Eng; Engine* FDTD_Eng;
ProcessingArray* PA; ProcessingArray* PA;