44 #include <boost/math/special_functions/fpclassify.hpp>
53 ModuleKey ProcessEquiSpacedOutput::className =
56 ProcessEquiSpacedOutput::create,
57 "Write data as equi-spaced output using simplices to represent the data for connecting points");
63 f->m_setUpEquiSpacedFields =
true;
66 "Only process tetrahedral elements");
87 cout <<
"Interpolating fields to equispaced" << endl;
90 int coordim =
m_f->m_exp[0]->GetCoordim(0);
91 int shapedim =
m_f->m_exp[0]->GetExp(0)->GetShapeDimension();
92 int npts =
m_f->m_exp[0]->GetTotPoints();
95 int nel =
m_f->m_exp[0]->GetExpSize();
107 map<int,StdRegions::Orientation > face0orient;
108 set<int> prismorient;
112 vector<std::string> fieldNames;
114 vector<Array<OneD, int> > ptsConn;
117 for(
int i = 0; i < nel; ++i)
119 e =
m_f->m_exp[0]->GetExp(i);
123 int fid = e->GetGeom()->GetFid(0);
124 if(face0orient.count(fid))
126 prismorient.insert(i);
146 for(
int i = 0; i < nel; ++i)
148 e =
m_f->m_exp[0]->GetExp(i);
151 int fid = e->GetGeom()->GetFid(2);
153 if(face0orient.count(fid))
168 prismorient.insert(i);
175 prismorient.insert(i);
182 for(
int i = 0; i < nel; ++i)
184 e =
m_f->m_exp[0]->GetExp(i);
187 if(
m_f->m_exp[0]->GetExp(i)->DetShapeType() !=
194 ppe.push_back(newpoints);
195 newtotpoints += newpoints;
197 switch(e->DetShapeType())
201 int npoints0 = e->GetBasis(0)->GetNumPoints();
209 int np0 = e->GetBasis(0)->GetNumPoints();
210 int np1 = e->GetBasis(1)->GetNumPoints();
211 int np = max(np0,np1);
218 int np0 = e->GetBasis(0)->GetNumPoints();
219 int np1 = e->GetBasis(1)->GetNumPoints();
220 int np = max(np0,np1);
228 int np0 = e->GetBasis(0)->GetNumPoints();
229 int np1 = e->GetBasis(1)->GetNumPoints();
230 int np2 = e->GetBasis(2)->GetNumPoints();
231 int np = max(np0,max(np1,np2));
239 int np0 = e->GetBasis(0)->GetNumPoints();
240 int np1 = e->GetBasis(1)->GetNumPoints();
241 int np2 = e->GetBasis(2)->GetNumPoints();
242 int np = max(np0,max(np1,np2));
250 int np0 = e->GetBasis(0)->GetNumPoints();
251 int np1 = e->GetBasis(1)->GetNumPoints();
252 int np2 = e->GetBasis(2)->GetNumPoints();
253 int np = max(np0,max(np1,np2));
261 int np0 = e->GetBasis(0)->GetNumPoints();
262 int np1 = e->GetBasis(1)->GetNumPoints();
263 int np2 = e->GetBasis(2)->GetNumPoints();
264 int np = max(np0,max(np1,np2));
278 bool standard =
true;
280 if(prismorient.count(i))
285 e->GetSimplexEquiSpacedConnectivity(conn,standard);
290 if((prevNcoeffs != e->GetNcoeffs()) ||
291 (prevNpoints != e->GetTotPoints()))
293 prevNcoeffs = e->GetNcoeffs();
294 prevNpoints = e->GetTotPoints();
296 e->GetSimplexEquiSpacedConnectivity(conn);
300 for(
int j = 0; j < conn.num_elements(); ++j)
302 newconn[j] = conn[j] + cnt;
305 ptsConn.push_back(newconn);
309 if(
m_f->m_fielddef.size())
311 nfields =
m_f->m_exp.size();
320 for(
int i = 0; i < nfields + coordim; ++i)
326 for(
int i = 0; i < coordim; ++i)
331 for(
int i = coordim; i < 3; ++i)
336 m_f->m_exp[0]->GetCoords(coords[0],coords[1],coords[2]);
338 int nq1 =
m_f->m_exp[0]->GetTotPoints();
344 m_f->m_exp[0]->GetCoords(x1, y1, z1);
349 for(
int n = 0; n < coordim; ++n)
353 for(
int i = 0; i < nel; ++i)
355 m_f->m_exp[0]->GetExp(i)->PhysInterpToSimplexEquiSpaced(
357 tmp = pts[n] + cnt1);
359 cnt +=
m_f->m_exp[0]->GetExp(i)->GetTotPoints();
363 if(
m_f->m_fielddef.size())
365 ASSERTL0(
m_f->m_fielddef[0]->m_fields.size() ==
m_f->m_exp.size(),
366 "More expansion defined than fields");
368 for(
int n = 0; n <
m_f->m_exp.size(); ++n)
373 if(
m_config[
"modalenergy"].m_beenSet)
376 for(
int i = 0; i < nel; ++i)
380 cnt +=
m_f->m_exp[0]->GetExp(i)->GetTotPoints();
386 for(
int i = 0; i < nel; ++i)
388 m_f->m_exp[0]->GetExp(i)->PhysInterpToSimplexEquiSpaced(
390 tmp = pts[coordim + n] + cnt1);
392 cnt +=
m_f->m_exp[0]->GetExp(i)->GetTotPoints();
397 fieldNames.push_back(
m_f->m_fielddef[0]->m_fields[n]);
406 else if (shapedim == 3)
410 m_f->m_fieldPts->SetConnectivity(ptsConn);
420 e =
m_f->m_exp[0]->GetExp(n);
422 switch(e->DetShapeType())
426 int np0 = e->GetBasis(0)->GetNumPoints();
427 int np1 = e->GetBasis(1)->GetNumPoints();
428 int np = max(np0,np1);
442 e->GetBasis(1)->GetBasisKey(),
450 int np0 = e->GetBasis(0)->GetNumPoints();
451 int np1 = e->GetBasis(1)->GetNumPoints();
452 int np = max(np0,np1);
464 e->GetBasis(1)->GetBasisKey(),
471 ASSERTL0(
false,
"Shape needs setting up");
#define ASSERTL0(condition, msg)
int getNumberOfCoefficients(int Na, int Nb, int Nc)
pair< ModuleType, string > ModuleKey
static Array< OneD, NekDouble > NullNekDouble1DArray
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
void SetupEquiSpacedField(void)
map< string, ConfigOption > m_config
List of configuration values.
void GenOrthoModes(int n, const Array< OneD, const NekDouble > &phys, Array< OneD, NekDouble > &coeffs)
int getNumberOfCoefficients(int Na, int Nb, int Nc)
FieldSharedPtr m_f
Field object.
Principle Orthogonal Functions .
void Interp2D(const BasisKey &fbasis0, const BasisKey &fbasis1, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, Array< OneD, NekDouble > &to)
this function interpolates a 2D function evaluated at the quadrature points of the 2D basis...
int getNumberOfCoefficients(int Na)
int getNumberOfCoefficients(int Na, int Nb)
Principle Orthogonal Functions .
Defines a specification for a set of points.
boost::shared_ptr< Field > FieldSharedPtr
Represents a command-line configuration option.
int getNumberOfCoefficients(int Na, int Nb, int Nc)
virtual ~ProcessEquiSpacedOutput()
int getNumberOfCoefficients(int Na, int Nb)
boost::shared_ptr< Expansion > ExpansionSharedPtr
Describes the specification for a Basis.
ModuleFactory & GetModuleFactory()
Abstract base class for processing modules.
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, tDescription pDesc="")
Register a class with the factory.
void FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function performs the Forward transformation from physical space to coefficient space...