Nektar++
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Classes | Typedefs | Enumerations | Functions | Variables
Nektar::FieldUtils Namespace Reference

Classes

struct  Field
 
class  FieldConvertComm
 
class  InputDat
 Input module for Xml files. More...
 
class  InputFld
 
class  InputNek5000
 
class  InputPts
 
class  InputSemtex
 
class  InputXml
 
class  Interpolator
 A class that contains algorithms for interpolation between pts fields, expansions and different meshes. More...
 
struct  ConfigOption
 Represents a command-line configuration option. More...
 
class  Module
 
class  InputModule
 Abstract base class for input modules. More...
 
class  ProcessModule
 Abstract base class for processing modules. More...
 
class  OutputModule
 Abstract base class for output modules. More...
 
class  Octree
 
class  OutputFileBase
 Converter from fld to vtk. More...
 
class  OutputFld
 Output to fld format. More...
 
class  OutputInfo
 
class  OutputPts
 Converter from fld to pts. More...
 
class  OutputStdOut
 
class  OutputTecplot
 Tecplot output class. More...
 
class  OutputTecplotBinary
 Tecplot output class, specifically for binary field output. More...
 
class  OutputVtk
 Converter from fld to vtk. More...
 
class  OutputXml
 Converter from fld to vtk. More...
 
class  ProcessAddCompositeID
 This processing module adds a fld with the composite ID. More...
 
class  ProcessAddFld
 This processing module scales the input fld file. More...
 
class  ProcessBoundaryExtract
 This processing module sets up for the boundary field to be extracted. More...
 
class  ProcessC0Projection
 This processing module calculates the Q Criterion and adds it as an extra-field to the output file. More...
 
class  ProcessCombineAvg
 This processing module combines two fld files containing average fields. More...
 
class  ProcessConcatenateFld
 This processing module sets up for the boundary field to be extracted. More...
 
class  ProcessCreateExp
 This processing module scales the input fld file. More...
 
class  ProcessDeform
 
struct  TriFaceIDs
 
struct  TriFaceHash
 
class  ProcessDisplacement
 
class  ProcessDOF
 This processing module calculates the number of DOF. More...
 
class  ProcessEquiSpacedOutput
 This processing module interpolates one field to another. More...
 
class  ProcessFieldFromString
 This processing module adds a new field from a string definition. More...
 
class  ProcessGrad
 This processing module calculates the vorticity and adds it as an extra-field to the output file. More...
 
class  ProcessHalfModeToFourier
 This processing modifies Fourier Half Mode to a Fourier expansion. More...
 
class  ProcessHomogeneousPlane
 This processing module replaces all expansions by a single plane from 3DH1D fields, defined by the parameter planeid. More...
 
class  ProcessHomogeneousStretch
 This processing module stretches the homogeneous direction of a 3DH1D expansion by an integer factor. More...
 
class  ProcessInnerProduct
 This processing module computes the inner product between two fields. More...
 
class  ProcessInterpField
 This processing module interpolates one field to another. More...
 
class  ProcessInterpPointDataToFld
 This processing module interpolates one field to another. More...
 
class  ProcessInterpPoints
 This processing module interpolates one field to another. More...
 
class  ProcessInterpPtsToPts
 This processing module interpolates one field to another. More...
 
class  Iso
 
class  IsoVertex
 
class  ProcessIsoContour
 This processing module extracts an isocontour. More...
 
class  ProcessJacobianEnergy
 This processing module scales the input fld file. More...
 
class  ProcessL2Criterion
 This processing module calculates the Lambda 2 Criterion and adds it as an extra-field to the output file. More...
 
class  ProcessMapping
 This processing module scales the input fld file. More...
 
class  ProcessMean
 This processing module computes the mean of each field. More...
 
class  ProcessMeanMode
 This processing module replaces all expansions by the mean mode from 3DH1D fields. More...
 
class  ProcessMultiShear
 This processing module calculates the shear stress metrics and writes it to a surface output file. More...
 
class  ProcessNumModes
 This processing module calculates the vorticity and adds it as an extra-field to the output file. More...
 
class  ProcessPhiFromFile
 
class  ProcessPointDataToFld
 This processing module interpolates one field to another. More...
 
class  ProcessPrintFldNorms
 This processing module prints the L2 and LInf norms of the variables in the field. More...
 
class  ProcessQCriterion
 This processing module calculates the Q Criterion and adds it as an extra-field to the output file. More...
 
class  ProcessQualityMetric
 This processing module scales the input fld file. More...
 
class  ProcessRemoveField
 This processing module adds a new field from a string definition. More...
 
class  ProcessScaleInFld
 This processing module scales the input fld file. More...
 
class  ProcessScalGrad
 This processing module calculates the scalar gradient field and writes it to a surface output file. More...
 
class  ProcessStreamFunction
 This processing module calculates the stream function of a 2D field and adds it as an extra-field to the output file. More...
 
class  ProcessSurfDistance
 This processing module calculates the height of an element connected to a surface and adds it as an extra-field to the output file. More...
 
class  ProcessVorticity
 This processing module calculates the vorticity and adds it as an extra-field to the output file. More...
 
class  ProcessWallNormalData
 This processing module calculates the wall shear stress and adds it as an extra-field to the output file, and writes it to a surface output file. More...
 
class  ProcessWSS
 This processing module calculates the wall shear stress and adds it as an extra-field to the output file, and writes it to a surface output file. More...
 

Typedefs

typedef std::shared_ptr< FieldFieldSharedPtr
 
typedef std::shared_ptr< InterpolatorInterpolatorSharedPtr
 
typedef std::shared_ptr< InputModuleInputModuleSharedPtr
 
typedef std::pair< ModuleType, std::string > ModuleKey
 
typedef std::shared_ptr< ModuleModuleSharedPtr
 
typedef LibUtilities::NekFactory< ModuleKey, Module, FieldSharedPtrModuleFactory
 
typedef std::unordered_map< TriFaceIDs, int, TriFaceHashTriFaceMap
 
typedef std::shared_ptr< IsoIsoSharedPtr
 

Enumerations

enum  ModuleType { eInputModule , eProcessModule , eOutputModule , SIZE_ModuleType }
 
enum  ModulePriority {
  eCreateGraph , eCreateFieldData , eModifyFieldData , eCreateExp ,
  eFillExp , eModifyExp , eBndExtraction , eCreatePts ,
  eConvertExpToPts , eModifyPts , eOutput , SIZE_ModulePriority
}
 
enum  TecplotZoneType {
  eOrdered = 0 , eFELineSeg , eFETriangle , eFEQuadrilateral ,
  eFETetrahedron , eFEBrick , eFEPolygon , eFEPolyhedron
}
 

Functions

ModuleFactoryGetModuleFactory ()
 
std::ostream & operator<< (std::ostream &os, const ModuleKey &rhs)
 
template<typename T >
void swap_endian (T &u)
 Swap endian ordering of the input variable. More...
 
template<typename T >
void swap_endian (std::vector< T > &u)
 
template<typename T >
void WriteStream (std::ostream &outfile, T data)
 Helper function to write binary data to stream. More...
 
template<>
void WriteStream (std::ostream &outfile, std::string data)
 Specialisation of WriteStream to support writing std::string. More...
 
template<typename T >
void WriteStream (std::ostream &outfile, Array< OneD, T > data)
 Specialisation of WriteStream to support writing Nektar::Array datatype. More...
 
template<typename T >
void WriteStream (std::ostream &outfile, std::vector< T > data)
 Specialisation of WriteStream to support writing std::vector datatype. More...
 
bool operator== (TriFaceIDs const &p1, TriFaceIDs const &p2)
 
void TwoPairs (Array< OneD, NekDouble > &cx, Array< OneD, NekDouble > &cy, Array< OneD, NekDouble > &cz, int &pr)
 
void ThreeSimilar (const int i, const int j, const int k, int &pr, int &ps)
 
bool operator== (const IsoVertex &x, const IsoVertex &y)
 
bool operator!= (const IsoVertex &x, const IsoVertex &y)
 
void MatSymEVals (NekDouble d1, NekDouble d2, NekDouble d3, NekDouble a, NekDouble b, NekDouble c, NekDouble &l1, NekDouble &l2, NekDouble &l3)
 Calculates eigenvalues of a 3x3 Symmetric matrix. More...
 
vector< DNekMatMappingIdealToRef (SpatialDomains::GeometrySharedPtr geom, StdRegions::StdExpansionSharedPtr chi)
 

Variables

const std::string ModuleTypeMap [] = {"Input", "Process", "Output"}
 
const char *const ModulePriorityMap []
 
std::string TecplotZoneTypeMap []
 

Typedef Documentation

◆ FieldSharedPtr

typedef std::shared_ptr<Field> Nektar::FieldUtils::FieldSharedPtr

Definition at line 989 of file Field.hpp.

◆ InputModuleSharedPtr

Definition at line 258 of file Module.h.

◆ InterpolatorSharedPtr

Definition at line 108 of file FieldUtils/Interpolator.h.

◆ IsoSharedPtr

typedef std::shared_ptr<Iso> Nektar::FieldUtils::IsoSharedPtr

Definition at line 186 of file ProcessIsoContour.h.

◆ ModuleFactory

Definition at line 296 of file Module.h.

◆ ModuleKey

typedef std::pair<ModuleType, std::string> Nektar::FieldUtils::ModuleKey

Definition at line 290 of file Module.h.

◆ ModuleSharedPtr

typedef std::shared_ptr<Module> Nektar::FieldUtils::ModuleSharedPtr

Definition at line 294 of file Module.h.

◆ TriFaceMap

typedef std::unordered_map<TriFaceIDs, int, TriFaceHash> Nektar::FieldUtils::TriFaceMap

Definition at line 95 of file ProcessDisplacement.cpp.

Enumeration Type Documentation

◆ ModulePriority

Enumerator
eCreateGraph 
eCreateFieldData 
eModifyFieldData 
eCreateExp 
eFillExp 
eModifyExp 
eBndExtraction 
eCreatePts 
eConvertExpToPts 
eModifyPts 
eOutput 
SIZE_ModulePriority 

Definition at line 76 of file Module.h.

◆ ModuleType

Denotes different types of mesh converter modules: so far only input, output and process modules are defined.

Enumerator
eInputModule 
eProcessModule 
eOutputModule 
SIZE_ModuleType 

Definition at line 66 of file Module.h.

67 {
72 };

◆ TecplotZoneType

Enumerator
eOrdered 
eFELineSeg 
eFETriangle 
eFEQuadrilateral 
eFETetrahedron 
eFEBrick 
eFEPolygon 
eFEPolyhedron 

Definition at line 46 of file OutputTecplot.h.

Function Documentation

◆ GetModuleFactory()

FIELD_UTILS_EXPORT ModuleFactory & Nektar::FieldUtils::GetModuleFactory ( )

Returns an instance of the module factory, held as a singleton.

Definition at line 49 of file Module.cpp.

50 {
51  static ModuleFactory instance;
52  return instance;
53 }
LibUtilities::NekFactory< ModuleKey, Module, FieldSharedPtr > ModuleFactory
Definition: Module.h:296

Referenced by Nektar::SolverUtils::FilterFieldConvert::CreateModules(), main(), Module_Create(), and Module_Register().

◆ MappingIdealToRef()

vector<DNekMat> Nektar::FieldUtils::MappingIdealToRef ( SpatialDomains::GeometrySharedPtr  geom,
StdRegions::StdExpansionSharedPtr  chi 
)
inline

Definition at line 118 of file ProcessQualityMetric.cpp.

120 {
121  vector<DNekMat> ret;
122 
123  if (geom->GetShapeType() == LibUtilities::eQuadrilateral)
124  {
125  vector<Array<OneD, NekDouble> > xy;
126  for (int i = 0; i < geom->GetNumVerts(); i++)
127  {
128  Array<OneD, NekDouble> loc(2);
129  SpatialDomains::PointGeomSharedPtr p = geom->GetVertex(i);
130  p->GetCoords(loc);
131  xy.push_back(loc);
132  }
133 
134  Array<OneD, const LibUtilities::BasisSharedPtr> b = chi->GetBase();
135  Array<OneD, NekDouble> u = b[0]->GetZ();
136  Array<OneD, NekDouble> v = b[1]->GetZ();
137 
138  for (int j = 0; j < b[1]->GetNumPoints(); j++)
139  {
140  for (int i = 0; i < b[0]->GetNumPoints(); i++)
141  {
142  NekDouble a1 = 0.5 * (1.0 - u[i]), a2 = 0.5 * (1.0 + u[i]);
143  NekDouble b1 = 0.5 * (1.0 - v[j]), b2 = 0.5 * (1.0 + v[j]);
144  DNekMat dxdz(2, 2, 1.0, eFULL);
145 
146  dxdz(0, 0) = 0.5 * (-b1 * xy[0][0] + b1 * xy[1][0] +
147  b2 * xy[2][0] - b2 * xy[3][0]);
148  dxdz(1, 0) = 0.5 * (-b1 * xy[0][1] + b1 * xy[1][1] +
149  b2 * xy[2][1] - b2 * xy[3][1]);
150 
151  dxdz(0, 1) = 0.5 * (-a1 * xy[0][0] - a2 * xy[1][0] +
152  a2 * xy[2][0] + a1 * xy[3][0]);
153  dxdz(1, 1) = 0.5 * (-a1 * xy[0][1] - a2 * xy[1][1] +
154  a2 * xy[2][1] + a1 * xy[3][1]);
155 
156  dxdz.Invert();
157  ret.push_back(dxdz);
158  }
159  }
160  }
161  else if (geom->GetShapeType() == LibUtilities::eTriangle)
162  {
163  vector<Array<OneD, NekDouble> > xy;
164  for (int i = 0; i < geom->GetNumVerts(); i++)
165  {
166  Array<OneD, NekDouble> loc(2);
167  SpatialDomains::PointGeomSharedPtr p = geom->GetVertex(i);
168  p->GetCoords(loc);
169  xy.push_back(loc);
170  }
171 
172  Array<OneD, const LibUtilities::BasisSharedPtr> b = chi->GetBase();
173  Array<OneD, NekDouble> u = b[0]->GetZ();
174  Array<OneD, NekDouble> v = b[1]->GetZ();
175 
176  for (int i = 0; i < b[0]->GetNumPoints(); i++)
177  {
178  for (int j = 0; j < b[1]->GetNumPoints(); j++)
179  {
180  DNekMat dxdz(2, 2, 1.0, eFULL);
181  dxdz(0, 0) = -xy[0][0] / 2.0 + xy[1][0] / 2.0;
182 
183  dxdz(0, 1) = -xy[0][0] / 2.0 + xy[2][0] / 2.0;
184 
185  dxdz(1, 0) = -xy[0][1] / 2.0 + xy[1][1] / 2.0;
186 
187  dxdz(1, 1) = -xy[0][1] / 2.0 + xy[2][1] / 2.0;
188 
189  dxdz.Invert();
190  ret.push_back(dxdz);
191  }
192  }
193  }
194  else if (geom->GetShapeType() == LibUtilities::eTetrahedron)
195  {
196  vector<Array<OneD, NekDouble> > xyz;
197  for (int i = 0; i < geom->GetNumVerts(); i++)
198  {
199  Array<OneD, NekDouble> loc(3);
200  SpatialDomains::PointGeomSharedPtr p = geom->GetVertex(i);
201  p->GetCoords(loc);
202  xyz.push_back(loc);
203  }
204 
205  Array<OneD, const LibUtilities::BasisSharedPtr> b = chi->GetBase();
206  Array<OneD, NekDouble> u = b[0]->GetZ();
207  Array<OneD, NekDouble> v = b[1]->GetZ();
208  Array<OneD, NekDouble> z = b[2]->GetZ();
209 
210  for (int i = 0; i < b[0]->GetNumPoints(); i++)
211  {
212  for (int j = 0; j < b[1]->GetNumPoints(); j++)
213  {
214  for (int k = 0; k < b[2]->GetNumPoints(); k++)
215  {
216  DNekMat dxdz(3, 3, 1.0, eFULL);
217  dxdz(0, 0) = -xyz[0][0] / 2.0 + xyz[1][0] / 2.0;
218 
219  dxdz(0, 1) = -xyz[0][0] / 2.0 + xyz[2][0] / 2.0;
220 
221  dxdz(0, 2) = -xyz[0][0] / 2.0 + xyz[3][0] / 2.0;
222 
223  dxdz(1, 0) = -xyz[0][1] / 2.0 + xyz[1][1] / 2.0;
224 
225  dxdz(1, 1) = -xyz[0][1] / 2.0 + xyz[2][1] / 2.0;
226 
227  dxdz(1, 2) = -xyz[0][1] / 2.0 + xyz[3][1] / 2.0;
228 
229  dxdz(2, 0) = -xyz[0][2] / 2.0 + xyz[1][2] / 2.0;
230 
231  dxdz(2, 1) = -xyz[0][2] / 2.0 + xyz[2][2] / 2.0;
232 
233  dxdz(2, 2) = -xyz[0][2] / 2.0 + xyz[3][2] / 2.0;
234 
235  dxdz.Invert();
236  ret.push_back(dxdz);
237  }
238  }
239  }
240  }
241  else if (geom->GetShapeType() == LibUtilities::ePrism)
242  {
243  vector<Array<OneD, NekDouble> > xyz;
244  for (int i = 0; i < geom->GetNumVerts(); i++)
245  {
246  Array<OneD, NekDouble> loc(3);
247  SpatialDomains::PointGeomSharedPtr p = geom->GetVertex(i);
248  p->GetCoords(loc);
249  xyz.push_back(loc);
250  }
251 
252  Array<OneD, const LibUtilities::BasisSharedPtr> b = chi->GetBase();
253  Array<OneD, NekDouble> eta1 = b[0]->GetZ();
254  Array<OneD, NekDouble> eta2 = b[1]->GetZ();
255  Array<OneD, NekDouble> eta3 = b[2]->GetZ();
256 
257  for (int k = 0; k < b[2]->GetNumPoints(); k++)
258  {
259  for (int j = 0; j < b[1]->GetNumPoints(); j++)
260  {
261  for (int i = 0; i < b[0]->GetNumPoints(); i++)
262  {
263  NekDouble xi1 = 0.5 * (1 + eta1[i]) * (1 - eta3[k]) - 1.0;
264  NekDouble a2 = 0.5 * (1 + xi1);
265  NekDouble b1 = 0.5 * (1 - eta2[j]),
266  b2 = 0.5 * (1 + eta2[j]);
267  NekDouble c1 = 0.5 * (1 - eta3[k]),
268  c2 = 0.5 * (1 + eta3[k]);
269 
270  DNekMat dxdz(3, 3, 1.0, eFULL);
271 
272  dxdz(0, 0) = 0.5 * (-b1 * xyz[0][0] + b1 * xyz[1][0] +
273  b2 * xyz[2][0] - b2 * xyz[3][0]);
274  dxdz(1, 0) = 0.5 * (-b1 * xyz[0][1] + b1 * xyz[1][1] +
275  b2 * xyz[2][1] - b2 * xyz[3][1]);
276  dxdz(2, 0) = 0.5 * (-b1 * xyz[0][2] + b1 * xyz[1][2] +
277  b2 * xyz[2][2] - b2 * xyz[3][2]);
278 
279  dxdz(0, 1) = 0.5 * ((a2 - c1) * xyz[0][0] - a2 * xyz[1][0] +
280  a2 * xyz[2][0] + (c1 - a2) * xyz[3][0] -
281  c2 * xyz[4][0] + c2 * xyz[5][0]);
282  dxdz(1, 1) = 0.5 * ((a2 - c1) * xyz[0][1] - a2 * xyz[1][1] +
283  a2 * xyz[2][1] + (c1 - a2) * xyz[3][1] -
284  c2 * xyz[4][1] + c2 * xyz[5][1]);
285  dxdz(2, 1) = 0.5 * ((a2 - c1) * xyz[0][2] - a2 * xyz[1][2] +
286  a2 * xyz[2][2] + (c1 - a2) * xyz[3][2] -
287  c2 * xyz[4][2] + c2 * xyz[5][2]);
288 
289  dxdz(0, 2) = 0.5 * (-b1 * xyz[0][0] - b2 * xyz[3][0] +
290  b1 * xyz[4][0] + b2 * xyz[5][0]);
291  dxdz(1, 2) = 0.5 * (-b1 * xyz[0][1] - b2 * xyz[3][1] +
292  b1 * xyz[4][1] + b2 * xyz[5][1]);
293  dxdz(2, 2) = 0.5 * (-b1 * xyz[0][2] - b2 * xyz[3][2] +
294  b1 * xyz[4][2] + b2 * xyz[5][2]);
295 
296  dxdz.Invert();
297  ret.push_back(dxdz);
298  }
299  }
300  }
301  }
302  else if (geom->GetShapeType() == LibUtilities::eHexahedron)
303  {
304  vector<Array<OneD, NekDouble> > xyz;
305  for (int i = 0; i < geom->GetNumVerts(); i++)
306  {
307  Array<OneD, NekDouble> loc(3);
308  SpatialDomains::PointGeomSharedPtr p = geom->GetVertex(i);
309  p->GetCoords(loc);
310  xyz.push_back(loc);
311  }
312 
313  Array<OneD, const LibUtilities::BasisSharedPtr> b = chi->GetBase();
314  Array<OneD, NekDouble> eta1 = b[0]->GetZ();
315  Array<OneD, NekDouble> eta2 = b[1]->GetZ();
316  Array<OneD, NekDouble> eta3 = b[2]->GetZ();
317 
318  for (int k = 0; k < b[2]->GetNumPoints(); k++)
319  {
320  for (int j = 0; j < b[1]->GetNumPoints(); j++)
321  {
322  for (int i = 0; i < b[0]->GetNumPoints(); i++)
323  {
324  NekDouble a1 = 0.5 * (1 - eta1[i]);
325  NekDouble a2 = 0.5 * (1 + eta1[i]);
326  NekDouble b1 = 0.5 * (1 - eta2[j]),
327  b2 = 0.5 * (1 + eta2[j]);
328  NekDouble c1 = 0.5 * (1 - eta3[k]),
329  c2 = 0.5 * (1 + eta3[k]);
330 
331  DNekMat dxdz(3, 3, 1.0, eFULL);
332 
333  dxdz(0, 0) =
334  -0.5 * b1 * c1 * xyz[0][0] + 0.5 * b1 * c1 * xyz[1][0] +
335  0.5 * b2 * c1 * xyz[2][0] - 0.5 * b2 * c1 * xyz[3][0] -
336  0.5 * b1 * c2 * xyz[5][0] + 0.5 * b1 * c2 * xyz[5][0] +
337  0.5 * b2 * c2 * xyz[6][0] - 0.5 * b2 * c2 * xyz[7][0];
338  dxdz(1, 0) =
339  -0.5 * b1 * c1 * xyz[0][1] + 0.5 * b1 * c1 * xyz[1][1] +
340  0.5 * b2 * c1 * xyz[2][1] - 0.5 * b2 * c1 * xyz[3][1] -
341  0.5 * b1 * c2 * xyz[5][1] + 0.5 * b1 * c2 * xyz[5][1] +
342  0.5 * b2 * c2 * xyz[6][1] - 0.5 * b2 * c2 * xyz[7][1];
343  dxdz(2, 0) =
344  -0.5 * b1 * c1 * xyz[0][2] + 0.5 * b1 * c1 * xyz[1][2] +
345  0.5 * b2 * c1 * xyz[2][2] - 0.5 * b2 * c1 * xyz[3][2] -
346  0.5 * b1 * c2 * xyz[5][2] + 0.5 * b1 * c2 * xyz[5][2] +
347  0.5 * b2 * c2 * xyz[6][2] - 0.5 * b2 * c2 * xyz[7][2];
348 
349  dxdz(0, 1) =
350  -0.5 * a1 * c1 * xyz[0][0] - 0.5 * a2 * c1 * xyz[1][0] +
351  0.5 * a2 * c1 * xyz[2][0] + 0.5 * a1 * c1 * xyz[3][0] -
352  0.5 * a1 * c2 * xyz[5][0] - 0.5 * a2 * c2 * xyz[5][0] +
353  0.5 * a2 * c2 * xyz[6][0] + 0.5 * a1 * c2 * xyz[7][0];
354  dxdz(1, 1) =
355  -0.5 * a1 * c1 * xyz[0][1] - 0.5 * a2 * c1 * xyz[1][1] +
356  0.5 * a2 * c1 * xyz[2][1] + 0.5 * a1 * c1 * xyz[3][1] -
357  0.5 * a1 * c2 * xyz[5][1] - 0.5 * a2 * c2 * xyz[5][1] +
358  0.5 * a2 * c2 * xyz[6][1] + 0.5 * a1 * c2 * xyz[7][1];
359  dxdz(2, 1) =
360  -0.5 * a1 * c1 * xyz[0][2] - 0.5 * a2 * c1 * xyz[1][2] +
361  0.5 * a2 * c1 * xyz[2][2] + 0.5 * a1 * c1 * xyz[3][2] -
362  0.5 * a1 * c2 * xyz[5][2] - 0.5 * a2 * c2 * xyz[5][2] +
363  0.5 * a2 * c2 * xyz[6][2] + 0.5 * a1 * c2 * xyz[7][2];
364 
365  dxdz(0, 0) =
366  -0.5 * b1 * a1 * xyz[0][0] - 0.5 * b1 * a2 * xyz[1][0] -
367  0.5 * b2 * a2 * xyz[2][0] - 0.5 * b2 * a1 * xyz[3][0] +
368  0.5 * b1 * a1 * xyz[5][0] + 0.5 * b1 * a2 * xyz[5][0] +
369  0.5 * b2 * a2 * xyz[6][0] + 0.5 * b2 * a1 * xyz[7][0];
370  dxdz(1, 0) =
371  -0.5 * b1 * a1 * xyz[0][1] - 0.5 * b1 * a2 * xyz[1][1] -
372  0.5 * b2 * a2 * xyz[2][1] - 0.5 * b2 * a1 * xyz[3][1] +
373  0.5 * b1 * a1 * xyz[5][1] + 0.5 * b1 * a2 * xyz[5][1] +
374  0.5 * b2 * a2 * xyz[6][1] + 0.5 * b2 * a1 * xyz[7][1];
375  dxdz(2, 0) =
376  -0.5 * b1 * a1 * xyz[0][2] - 0.5 * b1 * a2 * xyz[1][2] -
377  0.5 * b2 * a2 * xyz[2][2] - 0.5 * b2 * a1 * xyz[3][2] +
378  0.5 * b1 * a1 * xyz[5][2] + 0.5 * b1 * a2 * xyz[5][2] +
379  0.5 * b2 * a2 * xyz[6][2] + 0.5 * b2 * a1 * xyz[7][2];
380 
381  dxdz.Invert();
382  ret.push_back(dxdz);
383  }
384  }
385  }
386  }
387  else
388  {
389  ASSERTL0(false, "not coded");
390  }
391 
392  return ret;
393 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:216
std::shared_ptr< PointGeom > PointGeomSharedPtr
Definition: Geometry.h:59
NekMatrix< NekDouble, StandardMatrixTag > DNekMat
Definition: NekTypeDefs.hpp:51
double NekDouble

References ASSERTL0, Nektar::eFULL, Nektar::LibUtilities::eHexahedron, Nektar::LibUtilities::ePrism, Nektar::LibUtilities::eQuadrilateral, Nektar::LibUtilities::eTetrahedron, Nektar::LibUtilities::eTriangle, CG_Iterations::loc, and CellMLToNektar.cellml_metadata::p.

Referenced by Nektar::FieldUtils::ProcessQualityMetric::GetQ().

◆ MatSymEVals()

void Nektar::FieldUtils::MatSymEVals ( NekDouble  d1,
NekDouble  d2,
NekDouble  d3,
NekDouble  a,
NekDouble  b,
NekDouble  c,
NekDouble l1,
NekDouble l2,
NekDouble l3 
)

Calculates eigenvalues of a 3x3 Symmetric matrix.

Parameters
d1,d2,d3- matrix diagonal entries at [0,0], [1,1] and [2,2]
a- matrix value at [0,1] and [1,0]
b- matrix value at [0,2] and [2,0]
c- matrix value at [1,2] and [2,1]
l1,l2,l3the computed eigenvalues, ordered l3 >= l2 >= l1

Definition at line 72 of file ProcessL2Criterion.cpp.

75 {
76  NekDouble p = a * a + b * b + c * c;
77  if (p == 0)
78  {
79  l1 = d1;
80  l2 = d2;
81  l3 = d3;
82  if (l1 > l3)
83  {
84  swap(l1, l3);
85  }
86  if (l1 > l2)
87  {
88  swap(l1, l2);
89  }
90  if (l2 > l3)
91  {
92  swap(l2, l3);
93  }
94  }
95  else
96  {
97  NekDouble q = (d1 + d2 + d3) / 3.0;
98  p = (d1 - q) * (d1 - q) + (d2 - q) * (d2 - q) + (d3 - q) * (d3 - q) +
99  2.0 * p;
100  p = sqrt(p / 6.0);
101  NekDouble r =
102  -0.5 *
103  (a * a * d3 - a * a * q - 2.0 * a * b * c + b * b * d2 - b * b * q +
104  c * c * d1 - c * c * q - d1 * d2 * d3 + d1 * d2 * q + d1 * d3 * q -
105  d1 * q * q + d2 * d3 * q - d2 * q * q - d3 * q * q + q * q * q) /
106  (p * p * p);
107 
108  NekDouble phi = 0;
109  if (r <= -1)
110  {
111  phi = M_PI / 3.0;
112  }
113  else if (r >= 1)
114  {
115  phi = 0.0;
116  }
117  else
118  {
119  phi = acos(r) / 3.0;
120  }
121 
122  // the eigenvalues satisfy eig3 >= eig2 >= eig1
123  l3 = q + 2.0 * p * cos(phi);
124  l1 = q + 2.0 * p * cos(phi + (2.0 * M_PI / 3.0));
125  // since trace(A) = eig1 + eig2 + eig3
126  l2 = 3.0 * q - l1 - l3;
127  }
128 }
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:267

References CellMLToNektar.cellml_metadata::p, and tinysimd::sqrt().

Referenced by Nektar::FieldUtils::ProcessL2Criterion::Process().

◆ operator!=()

bool Nektar::FieldUtils::operator!= ( const IsoVertex x,
const IsoVertex y 
)

Definition at line 965 of file ProcessIsoContour.cpp.

966 {
967  return ((x.m_x-y.m_x)*(x.m_x-y.m_x) + (x.m_y-y.m_y)*(x.m_y-y.m_y) +
968  (x.m_z-y.m_z)*(x.m_z-y.m_z) < NekConstants::kNekZeroTol)? 0:1;
969 }
static const NekDouble kNekZeroTol

◆ operator<<()

FIELD_UTILS_EXPORT std::ostream & Nektar::FieldUtils::operator<< ( std::ostream &  os,
const ModuleKey rhs 
)

Prints a given module key to a stream.

Definition at line 58 of file Module.cpp.

59 {
60  return os << ModuleTypeMap[rhs.first] << ": " << rhs.second;
61 }
const std::string ModuleTypeMap[]
Definition: Module.h:74

References ModuleTypeMap.

◆ operator==() [1/2]

bool Nektar::FieldUtils::operator== ( const IsoVertex x,
const IsoVertex y 
)

Definition at line 958 of file ProcessIsoContour.cpp.

959 {
960  return ((x.m_x-y.m_x)*(x.m_x-y.m_x) + (x.m_y-y.m_y)*(x.m_y-y.m_y) +
961  (x.m_z-y.m_z)*(x.m_z-y.m_z) < NekConstants::kNekZeroTol)? true:false;
962 }

◆ operator==() [2/2]

bool Nektar::FieldUtils::operator== ( TriFaceIDs const &  p1,
TriFaceIDs const &  p2 
)

Definition at line 78 of file ProcessDisplacement.cpp.

79 {
80  std::vector<int> ids1(3), ids2(3);
81 
82  ids1[0] = p1.a;
83  ids1[1] = p1.b;
84  ids1[2] = p1.c;
85  ids2[0] = p2.a;
86  ids2[1] = p2.b;
87  ids2[2] = p2.c;
88 
89  std::sort(ids1.begin(), ids1.end());
90  std::sort(ids2.begin(), ids2.end());
91 
92  return ids1[0] == ids2[0] && ids1[1] == ids2[1] && ids1[2] == ids2[2];
93 }

References Nektar::FieldUtils::TriFaceIDs::a, Nektar::FieldUtils::TriFaceIDs::b, and Nektar::FieldUtils::TriFaceIDs::c.

◆ swap_endian() [1/2]

template<typename T >
void Nektar::FieldUtils::swap_endian ( std::vector< T > &  u)

Definition at line 121 of file Module.h.

122 {
123  size_t vecSize = u.size();
124  for (int i = 0; i < vecSize; ++i)
125  {
126  swap_endian(u[i]);
127  }
128 }
void swap_endian(std::vector< T > &u)
Definition: Module.h:121

References swap_endian().

◆ swap_endian() [2/2]

template<typename T >
void Nektar::FieldUtils::swap_endian ( T &  u)

Swap endian ordering of the input variable.

Definition at line 102 of file Module.h.

103 {
104  union
105  {
106  T u;
107  unsigned char u8[sizeof(T)];
108  } source, dest;
109 
110  source.u = u;
111 
112  for (size_t k = 0; k < sizeof(T); k++)
113  {
114  dest.u8[k] = source.u8[sizeof(T) - k - 1];
115  }
116 
117  u = dest.u;
118 }

Referenced by Nektar::FieldUtils::InputNek5000::Process(), Nektar::FieldUtils::InputSemtex::Process(), and swap_endian().

◆ ThreeSimilar()

void Nektar::FieldUtils::ThreeSimilar ( const int  i,
const int  j,
const int  k,
int &  pr,
int &  ps 
)

Definition at line 314 of file ProcessIsoContour.cpp.

318 {
319  switch (i + j + k)
320  {
321  case (3):
322  pr = 3;
323  ps = 4;
324  break;
325  case (4):
326  pr = 2;
327  ps = 4;
328  break;
329  case (5):
330  if (j == 1)
331  {
332  pr = 2;
333  ps = 3;
334  }
335  else
336  {
337  pr = 1;
338  ps = 4;
339  }
340  break;
341  case (6):
342  if (i == 0)
343  {
344  pr = 1;
345  ps = 3;
346  }
347  else
348  {
349  pr = 0;
350  ps = 4;
351  }
352  break;
353  case (7):
354  if (i == 0)
355  {
356  pr = 1;
357  ps = 2;
358  }
359  else
360  {
361  pr = 0;
362  ps = 3;
363  }
364  break;
365  case (8):
366  pr = 0;
367  ps = 2;
368  break;
369  case (9):
370  pr = 0;
371  ps = 1;
372  break;
373  default:
374  ASSERTL0(false,"Error in 5-point triangulation in ThreeSimilar");
375  break;
376  }
377 }

References ASSERTL0.

Referenced by Nektar::FieldUtils::ProcessIsoContour::ExtractContour().

◆ TwoPairs()

void Nektar::FieldUtils::TwoPairs ( Array< OneD, NekDouble > &  cx,
Array< OneD, NekDouble > &  cy,
Array< OneD, NekDouble > &  cz,
int &  pr 
)

Definition at line 282 of file ProcessIsoContour.cpp.

288 {
289  if (((cx[0]-cx[1])==0.0)&&
290  ((cy[0]-cy[1])==0.0)&&
291  ((cz[0]-cz[1])==0.0))
292  {
293  if (((cx[2]-cx[3])==0.0)&&
294  ((cy[2]-cy[3])==0.0)&&
295  ((cz[2]-cz[3])==0.0))
296  {
297  pr=4;
298  }
299  else
300  {
301  pr=3;
302  }
303  }
304  else
305  {
306  pr=1;
307  }
308 }

Referenced by Nektar::FieldUtils::ProcessIsoContour::ExtractContour().

◆ WriteStream() [1/4]

template<typename T >
void Nektar::FieldUtils::WriteStream ( std::ostream &  outfile,
Array< OneD, T >  data 
)

Specialisation of WriteStream to support writing Nektar::Array datatype.

Definition at line 138 of file OutputTecplot.cpp.

140 {
141  if (data.size())
142  {
143  outfile.write(reinterpret_cast<char *>(&data[0]),
144  data.size() * sizeof(T));
145  }
146 }

◆ WriteStream() [2/4]

template<>
void Nektar::FieldUtils::WriteStream ( std::ostream &  outfile,
std::string  data 
)

Specialisation of WriteStream to support writing std::string.

Tecplot binary formats represent all strings by writing out their characters as 32-bit integers, followed by a 32-bit integer null (0) character to denote the end of the string.

Definition at line 120 of file OutputTecplot.cpp.

121 {
122  // Convert string to array of int32_t
123  for (std::string::size_type i = 0; i < data.size(); ++i)
124  {
125  char strChar = data[i];
126  NekInt32 strCharInt = strChar;
127  WriteStream(outfile, strCharInt);
128  }
129 
130  // Now dump out zero character to terminate
131  WriteStream(outfile, 0);
132 }
void WriteStream(std::ostream &outfile, std::vector< T > data)
Specialisation of WriteStream to support writing std::vector datatype.
std::int32_t NekInt32

References WriteStream().

◆ WriteStream() [3/4]

template<typename T >
void Nektar::FieldUtils::WriteStream ( std::ostream &  outfile,
std::vector< T >  data 
)

Specialisation of WriteStream to support writing std::vector datatype.

Definition at line 151 of file OutputTecplot.cpp.

153 {
154  if (data.size())
155  {
156  outfile.write(reinterpret_cast<char *>(&data[0]),
157  data.size() * sizeof(T));
158  }
159 }

◆ WriteStream() [4/4]

template<typename T >
void Nektar::FieldUtils::WriteStream ( std::ostream &  outfile,
data 
)

Helper function to write binary data to stream.

Definition at line 107 of file OutputTecplot.cpp.

108 {
109  T tmp = data;
110  outfile.write(reinterpret_cast<char *>(&tmp), sizeof(T));
111 }

Referenced by Nektar::FieldUtils::OutputTecplotBinary::WriteDoubleOrFloat(), WriteStream(), Nektar::FieldUtils::OutputTecplotBinary::WriteTecplotConnectivity(), Nektar::FieldUtils::OutputTecplotBinary::WriteTecplotHeader(), and Nektar::FieldUtils::OutputTecplotBinary::WriteTecplotZone().

Variable Documentation

◆ ModulePriorityMap

const char* const Nektar::FieldUtils::ModulePriorityMap[]
Initial value:
= {
"CreateGraph", "CreateFieldData", "ModifyFieldData", "CreateExp",
"FillExp", "ModifyExp", "BndExtraction", "CreatePts",
"ConvertExpToPts", "ModifyPts", "Output"
}

Definition at line 92 of file Module.h.

◆ ModuleTypeMap

const std::string Nektar::FieldUtils::ModuleTypeMap[] = {"Input", "Process", "Output"}

Definition at line 74 of file Module.h.

Referenced by export_Module(), Module_Create(), Module_Register(), and operator<<().

◆ TecplotZoneTypeMap

std::string Nektar::FieldUtils::TecplotZoneTypeMap[]
Initial value:
= {
"ORDERED",
"LINESEG",
"TRIANGLE",
"QUADRILATERAL",
"TETRAHEDRON",
"BRICK",
"POLYGON",
"POLYHEDRON"
}

Definition at line 52 of file OutputTecplot.cpp.

Referenced by Nektar::FieldUtils::OutputTecplot::WriteTecplotZone().