From 5e5a278ac7a9eb348876f8fcd8f7da3681509127 Mon Sep 17 00:00:00 2001 From: Thorsten Liebig Date: Tue, 11 May 2010 20:38:58 +0200 Subject: [PATCH] new: node interpolated dump in processfields --- FDTD/processfields.cpp | 31 +++++++++- FDTD/processfields.h | 4 +- FDTD/processfields_td.cpp | 122 ++++++++++++++++++++++++++++++++++++++ FDTD/processfields_td.h | 1 + 4 files changed, 156 insertions(+), 2 deletions(-) diff --git a/FDTD/processfields.cpp b/FDTD/processfields.cpp index dae9c38..821c090 100644 --- a/FDTD/processfields.cpp +++ b/FDTD/processfields.cpp @@ -138,6 +138,35 @@ void ProcessFields::DefineStartStopCoord(double* dstart, double* dstop) discLines[n][i] = lines.at(i); } } + else if (m_DumpMode==NODE_INTERPOLATE) + { + if (Op->SnapToMesh(dstart,start)==false) cerr << "ProcessFields::DefineStartStopCoord: Warning: Snapping problem, check start value!!" << endl; + if (Op->SnapToMesh(dstop,stop)==false) cerr << "ProcessFields::DefineStartStopCoord: Warning: Snapping problem, check stop value!!" << endl; + + //create mesh + for (int n=0;n<3;++n) + { + if (start[n]>stop[n]) + { + unsigned int help = start[n]; + start[n]=stop[n]; + stop[n]=help; + } + if (stop[n]==Op->GetNumberOfLines(n)-1) + --stop[n]; +// cerr << "start " << start[n] << "stop " << stop[n]; + lines.clear(); + for (unsigned int i=start[n];i<=stop[n];i+=subSample[n]) + { + lines.push_back(Op->GetDiscLine(n,i));//0.5*(Op->discLines[n][i+1] + Op->discLines[n][i])); + } + numLines[n] = lines.size(); + delete[] discLines[n]; + discLines[n] = new double[numLines[n]]; + for (unsigned int i=0;iSnapToMesh(dstart,start,true)==false) cerr << "ProcessFields::DefineStartStopCoord: Warning: Snapping problem, check start value!!" << endl; @@ -158,7 +187,7 @@ void ProcessFields::DefineStartStopCoord(double* dstart, double* dstop) for (unsigned int i=start[n];iGetDiscLine(n,i,true));//0.5*(Op->discLines[n][i+1] + Op->discLines[n][i])); - } + } numDLines[n] = lines.size(); delete[] discDLines[n]; discDLines[n] = new double[numDLines[n]]; diff --git a/FDTD/processfields.h b/FDTD/processfields.h index 5aef9c9..35c0c8e 100644 --- a/FDTD/processfields.h +++ b/FDTD/processfields.h @@ -46,7 +46,9 @@ public: //! Define the Dump-Mode void SetDumpMode(DumpMode mode) {m_DumpMode=mode;} - //! This methode will dump all fields in the center of a main cell (dual-node) using 4 E-field and 2 H-fields per direction. (default) + //! This methode will dump all fields on a main cell node using 2 E-field and 4 H-fields per direction. + void SetDumpMode2Node() {m_DumpMode=NODE_INTERPOLATE;} + //! This methode will dump all fields in the center of a main cell (dual-node) using 4 E-field and 2 H-fields per direction. void SetDumpMode2Cell() {m_DumpMode=CELL_INTERPOLATE;} //! Set dump type: 0 for E-fields, 1 for H-fields, 2 for D-fields, 3 for B-fields, 4 for J-fields, etc... diff --git a/FDTD/processfields_td.cpp b/FDTD/processfields_td.cpp index 23704a9..e67c635 100644 --- a/FDTD/processfields_td.cpp +++ b/FDTD/processfields_td.cpp @@ -29,6 +29,126 @@ ProcessFieldsTD::~ProcessFieldsTD() { } +void ProcessFieldsTD::DumpNodeInterpol(string filename) +{ + if (m_DumpType==H_FIELD_DUMP) + { + //create array + FDTD_FLOAT**** H_T = Create_N_3DArray(numLines); + unsigned int pos[3] = {start[0],start[1],start[2]}; + unsigned int OpPos[3]; + double delta; +// cerr << "processing e-fields... " << endl; + for (pos[0]=0;pos[0]GetMeshDelta(0,OpPos,true); + if (delta) + { + H_T[0][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(0,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetCurr(0,OpPos[0],OpPos[1]+1,OpPos[2]) + Eng->GetCurr(0,OpPos[0],OpPos[1],OpPos[2]+1) + Eng->GetCurr(0,OpPos[0],OpPos[1]+1,OpPos[2]+1); + H_T[0][pos[0]][pos[1]][pos[2]] /= (4*delta); + } + //in y + delta = Op->GetMeshDelta(1,OpPos,true); + if (delta) + { + H_T[1][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(1,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetCurr(1,OpPos[0]+1,OpPos[1],OpPos[2]) + Eng->GetCurr(1,OpPos[0],OpPos[1],OpPos[2]+1) + Eng->GetCurr(1,OpPos[0]+1,OpPos[1],OpPos[2]+1); + H_T[1][pos[0]][pos[1]][pos[2]] /= (4*delta);//*Op->gridDelta); + } + //in z + delta = Op->GetMeshDelta(2,OpPos,true); + if (delta) + { + H_T[2][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(2,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetCurr(2,OpPos[0],OpPos[1]+1,OpPos[2]) + Eng->GetCurr(2,OpPos[0]+1,OpPos[1],OpPos[2]) + Eng->GetCurr(2,OpPos[0]+1,OpPos[1]+1,OpPos[2]); + H_T[2][pos[0]][pos[1]][pos[2]] /= (4*delta); + } + } + } + } + + if (m_fileType==VTK_FILETYPE) + { + ofstream file(filename.c_str()); + if (file.is_open()==false) { cerr << "ProcessFieldsTD::Process: can't open file '" << filename << "' for writing... abort! " << endl;}; + DumpVectorArray2VTK(file,string("H-Field"),H_T,discLines,numLines,m_precision); + file.close(); + } + else if (m_fileType==HDF5_FILETYPE) + { + stringstream ss; + ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); + DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),H_T,numLines,Eng->GetNumberOfTimesteps()*Op->GetTimestep()); + } + else + cerr << "ProcessFieldsTD::DumpNodeInterpol: unknown File-Type" << endl; + Delete_N_3DArray(H_T,numLines); + H_T = NULL; + } + + if (m_DumpType==E_FIELD_DUMP) + { + //create array + FDTD_FLOAT**** E_T = Create_N_3DArray(numLines); + unsigned int pos[3] = {start[0],start[1],start[2]}; + unsigned int OpPos[3]; + unsigned int OpPosUp[3]; + double delta, deltaUp, deltaRel; +// cerr << "processing h-fields... " << endl; + for (pos[0]=0;pos[0]GetMeshDelta(n,OpPos); + ++OpPosUp[n]; + deltaUp = Op->GetMeshDelta(n,OpPos); + deltaRel = delta / (delta+deltaUp); + if (delta*deltaUp) + { + E_T[n][pos[0]][pos[1]][pos[2]] = Eng->GetVolt(n,OpPos)*(1-deltaRel)/delta + Eng->GetVolt(n,OpPosUp)/deltaUp*deltaRel; + } + --OpPosUp[n]; + } + } + } + } + if (m_fileType==VTK_FILETYPE) + { + ofstream file(filename.c_str()); + if (file.is_open()==false) { cerr << "ProcessFieldsTD::Process: can't open file '" << filename << "' for writing... abort! " << endl;}; + DumpVectorArray2VTK(file,string("E-Field"),E_T,discLines,numLines,m_precision); + file.close(); + } + else if (m_fileType==HDF5_FILETYPE) + { + stringstream ss; + ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); + DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),E_T,numLines,(0.5+Eng->GetNumberOfTimesteps())*Op->GetTimestep()); + } + else + cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl; + Delete_N_3DArray(E_T,numLines); + E_T = NULL; + } +} + void ProcessFieldsTD::DumpCellInterpol(string filename) { if (m_DumpType==E_FIELD_DUMP) @@ -259,6 +379,8 @@ int ProcessFieldsTD::Process() ss << ".vtk"; if (m_DumpMode==NO_INTERPOLATION) DumpNoInterpol(ss.str()); + if (m_DumpMode==NODE_INTERPOLATE) + DumpNodeInterpol(ss.str()); if (m_DumpMode==CELL_INTERPOLATE) DumpCellInterpol(ss.str()); } diff --git a/FDTD/processfields_td.h b/FDTD/processfields_td.h index 84bfc63..6233762 100644 --- a/FDTD/processfields_td.h +++ b/FDTD/processfields_td.h @@ -35,6 +35,7 @@ protected: int pad_length; void DumpNoInterpol(string filename); + void DumpNodeInterpol(string filename); void DumpCellInterpol(string filename); };