Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Friends | List of all members
Nektar::GlobalMapping::MappingXYofXY Class Reference

#include <MappingXYofXY.h>

Inheritance diagram for Nektar::GlobalMapping::MappingXYofXY:
Inheritance graph
[legend]
Collaboration diagram for Nektar::GlobalMapping::MappingXYofXY:
Collaboration graph
[legend]

Static Public Member Functions

static GLOBAL_MAPPING_EXPORT
MappingSharedPtr 
create (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
 Creates an instance of this class. More...
 
- Static Public Member Functions inherited from Nektar::GlobalMapping::Mapping
static GLOBAL_MAPPING_EXPORT
MappingSharedPtr 
Load (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
 Return a pointer to the mapping, creating it on first call. More...
 

Static Public Attributes

static std::string className
 Name of the class. More...
 

Protected Member Functions

void CalculateMetricTensor ()
 
void CalculateChristoffel ()
 
 MappingXYofXY (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
 
virtual GLOBAL_MAPPING_EXPORT void v_InitObject (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
 
virtual GLOBAL_MAPPING_EXPORT void v_ContravarToCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_CovarToCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_ContravarFromCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_CovarFromCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_GetJacobian (Array< OneD, NekDouble > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_GetMetricTensor (Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_GetInvMetricTensor (Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelContravar (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelCovar (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_UpdateGeomInfo ()
 
- Protected Member Functions inherited from Nektar::GlobalMapping::Mapping
GLOBAL_MAPPING_EXPORT Mapping (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
 Constructor. More...
 
GLOBAL_MAPPING_EXPORT void EvaluateFunction (Array< OneD, MultiRegions::ExpListSharedPtr > pFields, LibUtilities::SessionReaderSharedPtr pSession, std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, NekDouble pTime=NekDouble(0))
 
GLOBAL_MAPPING_EXPORT void EvaluateTimeFunction (LibUtilities::SessionReaderSharedPtr pSession, std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, NekDouble pTime=NekDouble(0))
 
virtual GLOBAL_MAPPING_EXPORT void v_GetCartesianCoordinates (Array< OneD, NekDouble > &out0, Array< OneD, NekDouble > &out1, Array< OneD, NekDouble > &out2)
 
virtual GLOBAL_MAPPING_EXPORT void v_GetCoordVelocity (Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_DotGradJacobian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_LowerIndex (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_RaiseIndex (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_Divergence (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_VelocityLaplacian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble alpha)
 
virtual GLOBAL_MAPPING_EXPORT void v_gradgradU (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_CurlCurlField (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const bool generalized)
 
virtual GLOBAL_MAPPING_EXPORT void v_UpdateMapping (const NekDouble time, const Array< OneD, Array< OneD, NekDouble > > &coords=NullNekDoubleArrayofArray, const Array< OneD, Array< OneD, NekDouble > > &coordsVel=NullNekDoubleArrayofArray)
 
virtual GLOBAL_MAPPING_EXPORT void v_UpdateBCs (const NekDouble time)
 

Protected Attributes

Array< OneD, Array< OneD,
NekDouble > > 
m_metricTensor
 
Array< OneD, Array< OneD,
NekDouble > > 
m_Christoffel
 
- Protected Attributes inherited from Nektar::GlobalMapping::Mapping
LibUtilities::SessionReaderSharedPtr m_session
 Session reader. More...
 
LibUtilities::FieldIOSharedPtr m_fld
 
Array< OneD,
MultiRegions::ExpListSharedPtr
m_fields
 
Array< OneD, Array< OneD,
NekDouble > > 
m_coords
 Array with the Cartesian coordinates. More...
 
Array< OneD, Array< OneD,
NekDouble > > 
m_coordsVel
 Array with the velocity of the coordinates. More...
 
Array< OneD, Array< OneD,
NekDouble > > 
m_GeometricInfo
 Array with metric terms of the mapping. More...
 
int m_nConvectiveFields
 Number of velocity components. More...
 
string m_funcName
 Name of the function containing the coordinates. More...
 
string m_velFuncName
 Name of the function containing the velocity of the coordinates. More...
 
bool m_constantJacobian
 Flag defining if the Jacobian is constant. More...
 
bool m_timeDependent
 Flag defining if the Mapping is time-dependent. More...
 
bool m_fromFunction
 Flag defining if the Mapping is defined by a function. More...
 
Array< OneD, Array< OneD,
NekDouble > > 
m_wk1
 
Array< OneD, Array< OneD,
NekDouble > > 
m_wk2
 
Array< OneD, Array< OneD,
NekDouble > > 
m_tmp
 

Friends

class MemoryManager< MappingXYofXY >
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::GlobalMapping::Mapping
virtual GLOBAL_MAPPING_EXPORT ~Mapping ()
 Destructor. More...
 
GLOBAL_MAPPING_EXPORT void InitObject (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
 Initialise the mapping object. More...
 
GLOBAL_MAPPING_EXPORT void ReplaceField (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
 Replace the Expansion List used by the mapping. More...
 
GLOBAL_MAPPING_EXPORT void Output (LibUtilities::FieldMetaDataMap &fieldMetaDataMap, const std::string &outname)
 Output function called when a chk or fld file is written. More...
 
GLOBAL_MAPPING_EXPORT void ContravarToCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Convert a contravariant vector to the Cartesian system. More...
 
GLOBAL_MAPPING_EXPORT void CovarToCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Convert a covariant vector to the Cartesian system. More...
 
GLOBAL_MAPPING_EXPORT void ContravarFromCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Convert a contravariant vector to the transformed system. More...
 
GLOBAL_MAPPING_EXPORT void CovarFromCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Convert a covariant vector to the transformed system. More...
 
GLOBAL_MAPPING_EXPORT void GetCartesianCoordinates (Array< OneD, NekDouble > &out0, Array< OneD, NekDouble > &out1, Array< OneD, NekDouble > &out2)
 Get the Cartesian coordinates in the field. More...
 
GLOBAL_MAPPING_EXPORT void GetJacobian (Array< OneD, NekDouble > &outarray)
 Get the Jacobian of the transformation. More...
 
GLOBAL_MAPPING_EXPORT void DotGradJacobian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
 Calculate the dot product with the gradient of the Jacobian. More...
 
GLOBAL_MAPPING_EXPORT void GetMetricTensor (Array< OneD, Array< OneD, NekDouble > > &outarray)
 Get the metric tensor $g_{ij}$. More...
 
GLOBAL_MAPPING_EXPORT void GetInvMetricTensor (Array< OneD, Array< OneD, NekDouble > > &outarray)
 Get the inverse of metric tensor $g^{ij}$. More...
 
GLOBAL_MAPPING_EXPORT void LowerIndex (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Lower index of vector: $v_{i} = g_{ij}*v^{j}$. More...
 
GLOBAL_MAPPING_EXPORT void RaiseIndex (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Raise index of vector: $v^{i} = g^{ij}*v_{j}$. More...
 
GLOBAL_MAPPING_EXPORT void ApplyChristoffelContravar (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Apply the Christoffel symbols to a contravariant vector. More...
 
GLOBAL_MAPPING_EXPORT void ApplyChristoffelCovar (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Apply the Christoffel symbols to a covariant vector. More...
 
GLOBAL_MAPPING_EXPORT void GetCoordVelocity (Array< OneD, Array< OneD, NekDouble > > &outarray)
 Obtain the velocity of the coordinates. More...
 
GLOBAL_MAPPING_EXPORT void Divergence (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
 Calculate the generalised divergence operator. More...
 
GLOBAL_MAPPING_EXPORT void VelocityLaplacian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble alpha=0.0)
 Generalised (correction to the) velocity Laplacian operator. More...
 
GLOBAL_MAPPING_EXPORT void gradgradU (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Second order covariant derivatives of a contravariant vector. More...
 
GLOBAL_MAPPING_EXPORT void CurlCurlField (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const bool generalized)
 CurlCurl calculated on the whole field. More...
 
GLOBAL_MAPPING_EXPORT bool IsTimeDependent ()
 Get flag defining if mapping is time-dependent. More...
 
GLOBAL_MAPPING_EXPORT void SetTimeDependent (const bool value)
 Set flag defining if mapping is time-dependent. More...
 
GLOBAL_MAPPING_EXPORT bool IsFromFunction ()
 Get flag defining if mapping is defined by a function. More...
 
GLOBAL_MAPPING_EXPORT void SetFromFunction (const bool value)
 Set flag defining if mapping is defined by a function. More...
 
GLOBAL_MAPPING_EXPORT bool HasConstantJacobian ()
 Get flag defining if mapping has constant Jacobian. More...
 
GLOBAL_MAPPING_EXPORT bool IsDefined ()
 Get flag determining if the mapping was defined or is trivial. More...
 
GLOBAL_MAPPING_EXPORT void UpdateBCs (const NekDouble time)
 Update the Dirichlet Boundary Conditions when using Mappings. More...
 
GLOBAL_MAPPING_EXPORT void UpdateMapping (const NekDouble time, const Array< OneD, Array< OneD, NekDouble > > &coords=NullNekDoubleArrayofArray, const Array< OneD, Array< OneD, NekDouble > > &coordsVel=NullNekDoubleArrayofArray)
 Update the Mapping with new coordinates. More...
 
GLOBAL_MAPPING_EXPORT void UpdateGeomInfo ()
 Recompute the metric terms of the Mapping. More...
 
- Static Protected Attributes inherited from Nektar::GlobalMapping::Mapping
static MappingSharedPtr m_mappingPtr = MappingSharedPtr()
 
static bool m_init = false
 
static bool m_isDefined = false
 

Detailed Description

This class implements a mapping defined by the transformation

\[ \bar{x} = \bar{x}(x,y) \]

\[ \bar{y} = \bar{y}(x,y) \]

\[ \bar{z} = z \]

where $(\bar{x},\bar{y},\bar{z})$ are the Cartesian (physical) coordinates and $(x,y,z)$ are the transformed (computational) coordinates.

Definition at line 52 of file MappingXYofXY.h.

Constructor & Destructor Documentation

Nektar::GlobalMapping::MappingXYofXY::MappingXYofXY ( const LibUtilities::SessionReaderSharedPtr pSession,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields 
)
protected

Definition at line 58 of file MappingXYofXY.cpp.

61  : Mapping(pSession, pFields)
62 {
63 }
GLOBAL_MAPPING_EXPORT Mapping(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
Constructor.
Definition: Mapping.cpp:55

Member Function Documentation

void Nektar::GlobalMapping::MappingXYofXY::CalculateChristoffel ( )
protected

Definition at line 421 of file MappingXYofXY.cpp.

References Nektar::MultiRegions::DirCartesianMap, Nektar::GlobalMapping::Mapping::GetInvMetricTensor(), Nektar::GlobalMapping::Mapping::GetMetricTensor(), m_Christoffel, Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, Vmath::Smul(), Vmath::Vadd(), Vmath::Vsub(), and Vmath::Vvtvp().

Referenced by v_UpdateGeomInfo().

422 {
423  int physTot = m_fields[0]->GetTotPoints();
424  int nvel = m_nConvectiveFields;
425 
426  Array<OneD, Array<OneD, NekDouble> > G(nvel*nvel);
427  Array<OneD, Array<OneD, NekDouble> > G_inv(nvel*nvel);
428  Array<OneD, Array<OneD, NekDouble> > gradG(2*2*2);
431  // Allocate memory
432  for (int i = 0; i < gradG.num_elements(); i++)
433  {
434  gradG[i] = Array<OneD, NekDouble>(physTot, 0.0);
435  tmp[i] = Array<OneD, NekDouble>(physTot, 0.0);
436  }
437  for (int i = 0; i < G.num_elements(); i++)
438  {
439  G[i] = Array<OneD, NekDouble>(physTot, 0.0);
440  G_inv[i] = Array<OneD, NekDouble>(physTot, 0.0);
441  }
442 
443  // Get the metric tensor and its inverse
444  GetMetricTensor(G);
445  GetInvMetricTensor(G_inv);
446 
447  bool waveSpace = m_fields[0]->GetWaveSpace();
448  m_fields[0]->SetWaveSpace(false);
449  //Calculate gradients of g
450  // consider only 2 dimensions, since the 3rd is trivial
451  for (int i = 0; i <2; i++)
452  {
453  for(int j=0; j<2; j++)
454  {
455  for(int k=0; k<2; k++)
456  {
457  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[k],
458  G[i*nvel+j],gradG[i*2*2 + j*2 + k]);
459  }
460  }
461  }
462 
463  // Calculate tmp[p,j,k] = 1/2( gradG[pj,k]+ gradG[pk,j]-gradG[jk,p])
464  for (int p = 0; p <2; p++)
465  {
466  for (int j = 0; j < 2; j++)
467  {
468  for (int k = 0; k < 2; k++)
469  {
470  Vmath::Vadd(physTot, gradG[p*2*2 + j*2 + k], 1,
471  gradG[p*2*2 + k*2 + j], 1,
472  tmp[p*2*2 + j*2 + k], 1);
473  Vmath::Vsub(physTot, tmp[p*2*2 + j*2 + k], 1,
474  gradG[j*2*2 + k*2 + p], 1,
475  tmp[p*2*2 + j*2 + k], 1);
476  Vmath::Smul(physTot, 0.5, tmp[p*2*2 + j*2 + k], 1,
477  tmp[p*2*2 + j*2 + k], 1);
478  }
479  }
480  }
481 
482  // Calculate Christoffel symbols = g^ip tmp[p,j,k]
483  int n=0;
484  for (int i = 0; i <2; i++)
485  {
486  for (int j = 0; j < 2; j++)
487  {
488  for (int k = 0; k <= j; k++)
489  {
490  m_Christoffel[n] = Array<OneD, NekDouble>(physTot, 0.0);
491  for (int p = 0; p < 2; p++)
492  {
493  Vmath::Vvtvp(physTot, G_inv[i*nvel+p], 1,
494  tmp[p*2*2 + j*2 + k], 1,
495  m_Christoffel[n], 1,
496  m_Christoffel[n], 1);
497  }
498  n = n+1;
499  }
500  }
501  }
502 
503  m_fields[0]->SetWaveSpace(waveSpace);
504 }
GLOBAL_MAPPING_EXPORT void GetMetricTensor(Array< OneD, Array< OneD, NekDouble > > &outarray)
Get the metric tensor .
Definition: Mapping.h:178
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:428
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:199
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
Definition: Vmath.cpp:329
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
Array< OneD, Array< OneD, NekDouble > > m_Christoffel
Definition: MappingXYofXY.h:83
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:285
GLOBAL_MAPPING_EXPORT void GetInvMetricTensor(Array< OneD, Array< OneD, NekDouble > > &outarray)
Get the inverse of metric tensor .
Definition: Mapping.h:185
void Nektar::GlobalMapping::MappingXYofXY::CalculateMetricTensor ( )
protected

Definition at line 392 of file MappingXYofXY.cpp.

References Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_GeometricInfo, m_metricTensor, Vmath::Vmul(), and Vmath::Vvtvp().

Referenced by v_UpdateGeomInfo().

393 {
394  int physTot = m_fields[0]->GetTotPoints();
395  // Allocate memory
397  for (int i = 0; i < m_metricTensor.num_elements(); i++)
398  {
399  m_metricTensor[i] = Array<OneD, NekDouble>(physTot, 0.0);
400  }
401  // g_{1,1} = fx^2+gx^2
402  Vmath::Vmul(physTot, m_GeometricInfo[0], 1, m_GeometricInfo[0], 1,
403  m_metricTensor[0], 1);
404  Vmath::Vvtvp(physTot, m_GeometricInfo[2], 1, m_GeometricInfo[2], 1,
405  m_metricTensor[0], 1,
406  m_metricTensor[0], 1);
407  //g_{2,2} = fy^2+gy^2
408  Vmath::Vmul(physTot, m_GeometricInfo[1], 1, m_GeometricInfo[1], 1,
409  m_metricTensor[1], 1);
410  Vmath::Vvtvp(physTot, m_GeometricInfo[3], 1, m_GeometricInfo[3], 1,
411  m_metricTensor[1], 1,
412  m_metricTensor[1], 1);
413  //g_{1,2} = g_{2,1} = fy*fx+gx*gy
414  Vmath::Vmul(physTot, m_GeometricInfo[0], 1, m_GeometricInfo[1], 1,
415  m_metricTensor[2], 1);
416  Vmath::Vvtvp(physTot, m_GeometricInfo[2], 1, m_GeometricInfo[3], 1,
417  m_metricTensor[2], 1,
418  m_metricTensor[2], 1);
419 }
Array< OneD, Array< OneD, NekDouble > > m_GeometricInfo
Array with metric terms of the mapping.
Definition: Mapping.h:415
Array< OneD, Array< OneD, NekDouble > > m_metricTensor
Definition: MappingXYofXY.h:82
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:428
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169
static GLOBAL_MAPPING_EXPORT MappingSharedPtr Nektar::GlobalMapping::MappingXYofXY::create ( const LibUtilities::SessionReaderSharedPtr pSession,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const TiXmlElement *  pMapping 
)
inlinestatic

Creates an instance of this class.

Definition at line 60 of file MappingXYofXY.h.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and Nektar::GlobalMapping::MappingSharedPtr.

64  {
67  pFields);
68  p->InitObject(pFields, pMapping);
69  return p;
70  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
GLOBAL_MAPPING_EXPORT typedef boost::shared_ptr< Mapping > MappingSharedPtr
A shared pointer to a Mapping object.
Definition: Mapping.h:51
void Nektar::GlobalMapping::MappingXYofXY::v_ApplyChristoffelContravar ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 274 of file MappingXYofXY.cpp.

References m_Christoffel, Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, Vmath::Vmul(), and Vmath::Vvtvp().

277 {
278  int physTot = m_fields[0]->GetTotPoints();
279  int nvel = m_nConvectiveFields;
280 
281  for (int i = 0; i< nvel; i++)
282  {
283  for (int j = 0; j< nvel; j++)
284  {
285  outarray[i*nvel+j] = Array<OneD, NekDouble>(physTot,0.0);
286  }
287  }
288 
289  // Calculate non-zero terms
290 
291  // outarray(0,0) = U1 * m_Christoffel[0] + U2 * m_Christoffel[1]
292  Vmath::Vmul(physTot,m_Christoffel[0],1,inarray[0],1,
293  outarray[0*nvel+0],1);
294  Vmath::Vvtvp(physTot,m_Christoffel[1],1,inarray[1],1,
295  outarray[0*nvel+0],1,outarray[0*nvel+0],1);
296 
297  // outarray(0,1) = U1 * m_Christoffel[1] + U2 * m_Christoffel[2]
298  Vmath::Vmul(physTot,m_Christoffel[1],1,inarray[0],1,
299  outarray[0*nvel+1],1);
300  Vmath::Vvtvp(physTot,m_Christoffel[2],1,inarray[1],1,
301  outarray[0*nvel+1],1,outarray[0*nvel+1],1);
302 
303  // outarray(1,0) = U1 * m_Christoffel[3] + U2 * m_Christoffel[4]
304  Vmath::Vmul(physTot,m_Christoffel[3],1,inarray[0],1,
305  outarray[1*nvel+0],1);
306  Vmath::Vvtvp(physTot,m_Christoffel[4],1,inarray[1],1,
307  outarray[1*nvel+0],1,outarray[1*nvel+0],1);
308 
309  // outarray(1,1) = U1 * m_Christoffel[4] + U2 * m_Christoffel[5]
310  Vmath::Vmul(physTot,m_Christoffel[4],1,inarray[0],1,
311  outarray[1*nvel+1],1);
312  Vmath::Vvtvp(physTot,m_Christoffel[5],1,inarray[1],1,
313  outarray[1*nvel+1],1,outarray[1*nvel+1],1);
314 
315 }
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:428
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
Array< OneD, Array< OneD, NekDouble > > m_Christoffel
Definition: MappingXYofXY.h:83
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169
void Nektar::GlobalMapping::MappingXYofXY::v_ApplyChristoffelCovar ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 317 of file MappingXYofXY.cpp.

References m_Christoffel, Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, Vmath::Vmul(), and Vmath::Vvtvp().

320 {
321  int physTot = m_fields[0]->GetTotPoints();
322  int nvel = m_nConvectiveFields;
323 
324  for (int i = 0; i< nvel; i++)
325  {
326  for (int j = 0; j< nvel; j++)
327  {
328  outarray[i*nvel+j] = Array<OneD, NekDouble>(physTot,0.0);
329  }
330  }
331 
332  // Calculate non-zero terms
333 
334  // outarray(0,0) = U1 * m_Christoffel[0] + U2 * m_Christoffel[3]
335  Vmath::Vmul(physTot,m_Christoffel[0],1,inarray[0],1,
336  outarray[0*nvel+0],1);
337  Vmath::Vvtvp(physTot,m_Christoffel[3],1,inarray[1],1,
338  outarray[0*nvel+0],1,outarray[0*nvel+0],1);
339 
340  // outarray(0,1) = U1 * m_Christoffel[1] + U2 * m_Christoffel[4]
341  Vmath::Vmul(physTot,m_Christoffel[1],1,inarray[0],1,
342  outarray[0*nvel+1],1);
343  Vmath::Vvtvp(physTot,m_Christoffel[4],1,inarray[1],1,
344  outarray[0*nvel+1],1,outarray[0*nvel+1],1);
345 
346  // outarray(1,0) = U1 * m_Christoffel[1] + U2 * m_Christoffel[4]
347  Vmath::Vmul(physTot,m_Christoffel[1],1,inarray[0],1,
348  outarray[1*nvel+0],1);
349  Vmath::Vvtvp(physTot,m_Christoffel[4],1,inarray[1],1,
350  outarray[1*nvel+0],1,outarray[1*nvel+0],1);
351 
352  // outarray(1,1) = U1 * m_Christoffel[2] + U2 * m_Christoffel[5]
353  Vmath::Vmul(physTot,m_Christoffel[2],1,inarray[0],1,
354  outarray[1*nvel+1],1);
355  Vmath::Vvtvp(physTot,m_Christoffel[5],1,inarray[1],1,
356  outarray[1*nvel+1],1,outarray[1*nvel+1],1);
357 }
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:428
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
Array< OneD, Array< OneD, NekDouble > > m_Christoffel
Definition: MappingXYofXY.h:83
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169
void Nektar::GlobalMapping::MappingXYofXY::v_ContravarFromCartesian ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 138 of file MappingXYofXY.cpp.

References Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_GeometricInfo, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, Vmath::Vcopy(), Vmath::Vdiv(), Vmath::Vmul(), and Vmath::Vvtvm().

141 {
142  int physTot = m_fields[0]->GetTotPoints();
143  Array<OneD, NekDouble> wk(physTot, 0.0);
144 
145  // U1 = [gy*u1-fy*u2]/(fx*gy-gx*fy)
146  Vmath::Vmul(physTot, inarray[1], 1, m_GeometricInfo[1], 1,
147  outarray[0], 1);
148  Vmath::Vvtvm(physTot, inarray[0], 1, m_GeometricInfo[3], 1,
149  outarray[0], 1,
150  outarray[0], 1);
151  Vmath::Vdiv(physTot, outarray[0], 1, m_GeometricInfo[4], 1,
152  outarray[0], 1);
153 
154  // U2 = [fx*u2-gx*u1]/(fx*gy-gx*fy)
155  Vmath::Vmul(physTot, inarray[0], 1, m_GeometricInfo[2], 1,
156  outarray[1], 1);
157  Vmath::Vvtvm(physTot, inarray[1], 1, m_GeometricInfo[0], 1,
158  outarray[1], 1,
159  outarray[1], 1);
160  Vmath::Vdiv(physTot, outarray[1], 1, m_GeometricInfo[4], 1,
161  outarray[1], 1);
162 
163  // U3 = u3
164  if (m_nConvectiveFields ==3)
165  {
166  Vmath::Vcopy(physTot, inarray[2], 1, outarray[2], 1);
167  }
168 }
Array< OneD, Array< OneD, NekDouble > > m_GeometricInfo
Array with metric terms of the mapping.
Definition: Mapping.h:415
void Vdiv(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x/y.
Definition: Vmath.cpp:227
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Vvtvm(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvm (vector times vector plus vector): z = w*x - y
Definition: Vmath.cpp:451
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169
void Nektar::GlobalMapping::MappingXYofXY::v_ContravarToCartesian ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 81 of file MappingXYofXY.cpp.

References Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_GeometricInfo, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, Vmath::Vcopy(), Vmath::Vmul(), and Vmath::Vvtvp().

84 {
85  int physTot = m_fields[0]->GetTotPoints();
86 
87  // U1 = fx*u1 + fy*u2
88  Vmath::Vmul(physTot, m_GeometricInfo[0], 1, inarray[0], 1,
89  outarray[0], 1);
90  Vmath::Vvtvp(physTot, m_GeometricInfo[1], 1, inarray[1], 1,
91  outarray[0], 1, outarray[0],1);
92 
93  // U2 = gx*u1+gy*u2
94  Vmath::Vmul(physTot, m_GeometricInfo[2], 1, inarray[0], 1,
95  outarray[1], 1);
96  Vmath::Vvtvp(physTot, m_GeometricInfo[3], 1, inarray[1], 1,
97  outarray[1], 1, outarray[1],1);
98 
99  // U3 = u3
100  if (m_nConvectiveFields ==3)
101  {
102  Vmath::Vcopy(physTot, inarray[2], 1, outarray[2], 1);
103  }
104 }
Array< OneD, Array< OneD, NekDouble > > m_GeometricInfo
Array with metric terms of the mapping.
Definition: Mapping.h:415
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:428
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169
void Nektar::GlobalMapping::MappingXYofXY::v_CovarFromCartesian ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 170 of file MappingXYofXY.cpp.

References Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_GeometricInfo, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, Vmath::Vcopy(), Vmath::Vmul(), and Vmath::Vvtvp().

173 {
174  int physTot = m_fields[0]->GetTotPoints();
175 
176  // U1 = u1*fx +gx*u2
177  Vmath::Vmul(physTot, m_GeometricInfo[0], 1, inarray[0], 1,
178  outarray[0], 1);
179  Vmath::Vvtvp(physTot, m_GeometricInfo[2], 1, inarray[1], 1,
180  outarray[0], 1, outarray[0],1);
181 
182  // U2 = fy*u1 + gy*u2
183  Vmath::Vmul(physTot, m_GeometricInfo[1], 1, inarray[0], 1,
184  outarray[1], 1);
185  Vmath::Vvtvp(physTot, m_GeometricInfo[3], 1, inarray[1], 1,
186  outarray[1], 1, outarray[1],1);
187 
188  // U3 = u3
189  if (m_nConvectiveFields ==3)
190  {
191  Vmath::Vcopy(physTot, inarray[2], 1, outarray[2], 1);
192  }
193 }
Array< OneD, Array< OneD, NekDouble > > m_GeometricInfo
Array with metric terms of the mapping.
Definition: Mapping.h:415
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:428
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169
void Nektar::GlobalMapping::MappingXYofXY::v_CovarToCartesian ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 106 of file MappingXYofXY.cpp.

References Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_GeometricInfo, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, Vmath::Vcopy(), Vmath::Vdiv(), Vmath::Vmul(), and Vmath::Vvtvm().

109 {
110  int physTot = m_fields[0]->GetTotPoints();
111  Array<OneD, NekDouble> wk(physTot, 0.0);
112 
113  // U1 = [gy*u1-gx*u2]/(fx*gy-gx*fy)
114  Vmath::Vmul(physTot, inarray[1], 1, m_GeometricInfo[2], 1,
115  outarray[0], 1);
116  Vmath::Vvtvm(physTot, inarray[0], 1, m_GeometricInfo[3], 1,
117  outarray[0], 1,
118  outarray[0], 1);
119  Vmath::Vdiv(physTot, outarray[0], 1, m_GeometricInfo[4], 1,
120  outarray[0], 1);
121 
122  // U2 = [fx*u2 - fy*u1]/(fx*gy-gx*fy)
123  Vmath::Vmul(physTot, inarray[0], 1, m_GeometricInfo[1], 1,
124  outarray[1], 1);
125  Vmath::Vvtvm(physTot, inarray[1], 1, m_GeometricInfo[0], 1,
126  outarray[1], 1,
127  outarray[1], 1);
128  Vmath::Vdiv(physTot, outarray[1], 1, m_GeometricInfo[4], 1,
129  outarray[1], 1);
130 
131  // U3 = u3
132  if (m_nConvectiveFields ==3)
133  {
134  Vmath::Vcopy(physTot, inarray[2], 1, outarray[2], 1);
135  }
136 }
Array< OneD, Array< OneD, NekDouble > > m_GeometricInfo
Array with metric terms of the mapping.
Definition: Mapping.h:415
void Vdiv(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x/y.
Definition: Vmath.cpp:227
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Vvtvm(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvm (vector times vector plus vector): z = w*x - y
Definition: Vmath.cpp:451
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169
void Nektar::GlobalMapping::MappingXYofXY::v_GetInvMetricTensor ( Array< OneD, Array< OneD, NekDouble > > &  outarray)
protectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 231 of file MappingXYofXY.cpp.

References Nektar::GlobalMapping::Mapping::GetJacobian(), Nektar::GlobalMapping::Mapping::m_fields, m_metricTensor, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, Vmath::Neg(), Vmath::Sadd(), Vmath::Vcopy(), Vmath::Vdiv(), and Vmath::Vmul().

233 {
234  int physTot = m_fields[0]->GetTotPoints();
235  int nvel = m_nConvectiveFields;
236 
237  for (int i=0; i<nvel*nvel; i++)
238  {
239  outarray[i] = Array<OneD, NekDouble> (physTot, 0.0);
240  }
241 
242  // Get Jacobian
243  Array<OneD, NekDouble> Jac(physTot, 0.0);
244  GetJacobian(Jac);
245 
246  // Get Jacobian squared
247  Array<OneD, NekDouble> wk(physTot, 0.0);
248  Vmath::Vmul(physTot, Jac, 1, Jac, 1, wk, 1);
249  // G^{1,1} = m_metricTensor[1]/Jac^2
250  Vmath::Vcopy(physTot, m_metricTensor[1], 1, outarray[0*nvel+0], 1);
251  Vmath::Vdiv(physTot, outarray[0*nvel+0], 1, wk,1,
252  outarray[0*nvel+0], 1);
253 
254  // G^{2,2} = m_metricTensor[0]/Jac^2
255  Vmath::Vcopy(physTot, m_metricTensor[0], 1, outarray[1*nvel+1], 1);
256  Vmath::Vdiv(physTot, outarray[1*nvel+1], 1, wk,1,
257  outarray[1*nvel+1], 1);
258 
259  // G^{1,2} = G^{2,1} = -m_metricTensor[2]/Jac^2
260  Vmath::Vcopy(physTot, m_metricTensor[2], 1, outarray[0*nvel+1], 1);
261  Vmath::Neg(physTot, outarray[0*nvel+1], 1);
262  Vmath::Vdiv(physTot, outarray[0*nvel+1], 1, wk,1,
263  outarray[0*nvel+1], 1);
264  Vmath::Vcopy(physTot, outarray[0*nvel+1], 1, outarray[1*nvel+0], 1);
265 
266  // G^{3,3} = 1
267  if (m_nConvectiveFields ==3)
268  {
269  Vmath::Sadd(physTot, 1.0, outarray[2*nvel+2], 1,
270  outarray[2*nvel+2], 1);
271  }
272 }
Array< OneD, Array< OneD, NekDouble > > m_metricTensor
Definition: MappingXYofXY.h:82
void Vdiv(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x/y.
Definition: Vmath.cpp:227
GLOBAL_MAPPING_EXPORT void GetJacobian(Array< OneD, NekDouble > &outarray)
Get the Jacobian of the transformation.
Definition: Mapping.h:155
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:382
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
Definition: Vmath.cpp:301
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169
void Nektar::GlobalMapping::MappingXYofXY::v_GetJacobian ( Array< OneD, NekDouble > &  outarray)
protectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 195 of file MappingXYofXY.cpp.

References Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_GeometricInfo, and Vmath::Vabs().

197 {
198  int physTot = m_fields[0]->GetTotPoints();
199  Vmath::Vabs(physTot, m_GeometricInfo[4], 1, outarray, 1);
200 }
Array< OneD, Array< OneD, NekDouble > > m_GeometricInfo
Array with metric terms of the mapping.
Definition: Mapping.h:415
void Vabs(int n, const T *x, const int incx, T *y, const int incy)
vabs: y = |x|
Definition: Vmath.cpp:410
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Nektar::GlobalMapping::MappingXYofXY::v_GetMetricTensor ( Array< OneD, Array< OneD, NekDouble > > &  outarray)
protectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 202 of file MappingXYofXY.cpp.

References Nektar::GlobalMapping::Mapping::m_fields, m_metricTensor, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, Vmath::Sadd(), and Vmath::Vcopy().

204 {
205  int physTot = m_fields[0]->GetTotPoints();
206  int nvel = m_nConvectiveFields;
207 
208  for (int i=0; i<nvel*nvel; i++)
209  {
210  outarray[i] = Array<OneD, NekDouble> (physTot, 0.0);
211  }
212 
213  // g_{1,1} = m_metricTensor[0]
214  Vmath::Vcopy(physTot, m_metricTensor[0], 1, outarray[0*nvel+0], 1);
215 
216  // g_{2,2} = m_metricTensor[1]
217  Vmath::Vcopy(physTot, m_metricTensor[1], 1, outarray[1*nvel+1], 1);
218 
219  // g_{1,2}=g{2,1} = m_metricTensor[2]
220  Vmath::Vcopy(physTot, m_metricTensor[2], 1, outarray[0*nvel+1], 1);
221  Vmath::Vcopy(physTot, m_metricTensor[2], 1, outarray[1*nvel+0], 1);
222 
223  // g_{3,3} = 1
224  if (m_nConvectiveFields ==3)
225  {
226  Vmath::Sadd(physTot, 1.0, outarray[2*nvel+2], 1,
227  outarray[2*nvel+2], 1);
228  }
229 }
Array< OneD, Array< OneD, NekDouble > > m_metricTensor
Definition: MappingXYofXY.h:82
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
Definition: Vmath.cpp:301
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Nektar::GlobalMapping::MappingXYofXY::v_InitObject ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const TiXmlElement *  pMapping 
)
protectedvirtual

This function initialises the Mapping object. It computes the coordinates and velocity coordinates, initialises the workspace variables, and calls UpdateGeomInfo, which will perform the calculations specific for each type of Mapping.

Parameters
pFieldsExpList array used in the mapping
pMappingxml element describing the mapping

Reimplemented from Nektar::GlobalMapping::Mapping.

Definition at line 69 of file MappingXYofXY.cpp.

References ASSERTL0, Nektar::GlobalMapping::Mapping::m_constantJacobian, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, and Nektar::GlobalMapping::Mapping::v_InitObject().

72 {
73  Mapping::v_InitObject(pFields, pMapping);
74 
75  m_constantJacobian = false;
76 
78  "Mapping X = X(x,y), Y = Y(x,y) needs 2 velocity components.");
79 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
virtual GLOBAL_MAPPING_EXPORT void v_InitObject(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
Definition: Mapping.cpp:98
bool m_constantJacobian
Flag defining if the Jacobian is constant.
Definition: Mapping.h:427
void Nektar::GlobalMapping::MappingXYofXY::v_UpdateGeomInfo ( )
protectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 359 of file MappingXYofXY.cpp.

References CalculateChristoffel(), CalculateMetricTensor(), Nektar::MultiRegions::DirCartesianMap, Nektar::GlobalMapping::Mapping::m_coords, Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_GeometricInfo, Vmath::Vmul(), and Vmath::Vvtvm().

360 {
361  int phystot = m_fields[0]->GetTotPoints();
362  // Allocation of geometry memory
364  for (int i = 0; i < m_GeometricInfo.num_elements(); i++)
365  {
366  m_GeometricInfo[i] = Array<OneD, NekDouble>(phystot, 0.0);
367  }
368 
369  bool waveSpace = m_fields[0]->GetWaveSpace();
370  m_fields[0]->SetWaveSpace(false);
371 
372  // Calculate derivatives of x transformation --> m_GeometricInfo 0-1
375 
376  // Calculate derivatives of y transformation m_GeometricInfo 2-3
379 
380  // Calculate fx*gy-gx*fy --> m_GeometricInfo4
381  Vmath::Vmul(phystot, m_GeometricInfo[1], 1, m_GeometricInfo[2], 1, m_GeometricInfo[4], 1);
382  Vmath::Vvtvm(phystot, m_GeometricInfo[0], 1, m_GeometricInfo[3], 1,
383  m_GeometricInfo[4], 1,
384  m_GeometricInfo[4], 1);
385  //
388 
389  m_fields[0]->SetWaveSpace(waveSpace);
390 }
Array< OneD, Array< OneD, NekDouble > > m_GeometricInfo
Array with metric terms of the mapping.
Definition: Mapping.h:415
Array< OneD, Array< OneD, NekDouble > > m_coords
Array with the Cartesian coordinates.
Definition: Mapping.h:411
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
void Vvtvm(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvm (vector times vector plus vector): z = w*x - y
Definition: Vmath.cpp:451
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169

Friends And Related Function Documentation

friend class MemoryManager< MappingXYofXY >
friend

Definition at line 56 of file MappingXYofXY.h.

Member Data Documentation

std::string Nektar::GlobalMapping::MappingXYofXY::className
static
Initial value:
=
MappingXYofXY::create, "X = X(x,y), Y = Y(x,y)")

Name of the class.

Definition at line 73 of file MappingXYofXY.h.

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::MappingXYofXY::m_Christoffel
protected
Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::MappingXYofXY::m_metricTensor
protected

Definition at line 82 of file MappingXYofXY.h.

Referenced by CalculateMetricTensor(), v_GetInvMetricTensor(), and v_GetMetricTensor().