Nektar++
Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
Nektar::GlobalMapping::Mapping Class Referenceabstract

Base class for mapping to be applied to the coordinate system. More...

#include <Mapping.h>

Inheritance diagram for Nektar::GlobalMapping::Mapping:
[legend]

Public Member Functions

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 (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 Public Member Functions

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...
 

Protected Member Functions

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_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)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_CovarToCartesian (const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_ContravarFromCartesian (const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_CovarFromCartesian (const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=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_GetJacobian (Array< OneD, NekDouble > &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_DotGradJacobian (const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, NekDouble > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_GetMetricTensor (Array< OneD, Array< OneD, NekDouble >> &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_GetInvMetricTensor (Array< OneD, Array< OneD, NekDouble >> &outarray)=0
 
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_ApplyChristoffelContravar (const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelCovar (const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=0
 
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 (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_UpdateGeomInfo ()=0
 
virtual GLOBAL_MAPPING_EXPORT void v_UpdateBCs (const NekDouble time)
 

Protected Attributes

LibUtilities::SessionReaderSharedPtr m_session
 Session reader. More...
 
LibUtilities::FieldIOSharedPtr m_fld
 
Array< OneD, MultiRegions::ExpListSharedPtrm_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...
 
std::string m_funcName
 Name of the function containing the coordinates. More...
 
std::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
 

Static Protected Attributes

static MappingSharedPtr m_mappingPtr = MappingSharedPtr()
 
static bool m_init = false
 
static bool m_isDefined = false
 

Detailed Description

Base class for mapping to be applied to the coordinate system.

Definition at line 68 of file Mapping.h.

Constructor & Destructor Documentation

◆ ~Mapping()

virtual GLOBAL_MAPPING_EXPORT Nektar::GlobalMapping::Mapping::~Mapping ( )
inlinevirtual

Destructor.

Definition at line 72 of file Mapping.h.

73  {
74  }

◆ Mapping()

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

Constructor.

Definition at line 58 of file Mapping.cpp.

60  : m_session(pSession), m_fields(pFields)
61 {
62  switch (m_fields[0]->GetExpType())
63  {
64  case MultiRegions::e1D:
65  {
67  }
68  break;
69 
70  case MultiRegions::e2D:
71  {
73  }
74  break;
75 
76  case MultiRegions::e3D:
79  {
81  }
82  break;
83 
84  default:
85  ASSERTL0(0, "Dimension not supported");
86  break;
87  }
88 
90 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:215
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:414
LibUtilities::SessionReaderSharedPtr m_session
Session reader.
Definition: Mapping.h:402
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:406
LibUtilities::FieldIOSharedPtr m_fld
Definition: Mapping.h:404
static std::shared_ptr< FieldIO > CreateDefault(const LibUtilities::SessionReaderSharedPtr session)
Returns an object for the default FieldIO method.
Definition: FieldIO.cpp:198

References ASSERTL0, Nektar::LibUtilities::FieldIO::CreateDefault(), Nektar::MultiRegions::e1D, Nektar::MultiRegions::e2D, Nektar::MultiRegions::e3D, Nektar::MultiRegions::e3DH1D, Nektar::MultiRegions::e3DH2D, m_fields, m_fld, and m_nConvectiveFields.

Member Function Documentation

◆ ApplyChristoffelContravar()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::ApplyChristoffelContravar ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
inline

Apply the Christoffel symbols to a contravariant vector.

This function is used apply the Christoffel symbols \( \left( i,pk\right)\) to a contravariant vector \(u^p\), following the relation

\[ (out)^{i}_{k} = \left( i,pk\right)u^p \]

Parameters
inarrayContravariant vector \(u^p\)
outarrayResult of applying Christoffel symbols to \(u^p\)

Definition at line 211 of file Mapping.h.

214  {
215  v_ApplyChristoffelContravar(inarray, outarray);
216  }
virtual GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelContravar(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=0

References v_ApplyChristoffelContravar().

Referenced by v_gradgradU(), and v_VelocityLaplacian().

◆ ApplyChristoffelCovar()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::ApplyChristoffelCovar ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
inline

Apply the Christoffel symbols to a covariant vector.

This function is used apply the Christoffel symbols \( \left( p,ik\right)\) to a covariant vector \(u_p\), following the relation

\[ (out)_{ik} = \left( p,ik\right)u_p \]

Parameters
inarrayContravariant vector \(u_p\)
outarrayResult of applying Christoffel symbols to \(u_p\)

Definition at line 229 of file Mapping.h.

232  {
233  v_ApplyChristoffelCovar(inarray, outarray);
234  }
virtual GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelCovar(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=0

References v_ApplyChristoffelCovar().

Referenced by v_gradgradU().

◆ ContravarFromCartesian()

void Nektar::GlobalMapping::Mapping::ContravarFromCartesian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)

Convert a contravariant vector to the transformed system.

This function converts a contravariant vector in Cartesian space \(\bar{v}^{i}\) to the corresponding \(v^{j}\) in transformed space using the relation

\[v^{j} = \frac{\partial x^j}{\partial \bar{x}^i}\bar{v}^{i}\]

Parameters
inarrayComponents of the vector in Cartesian space ( \(\bar{v}^{i}\))
outarrayComponents of the vector in transformed space ( \(v^{j}\))

Definition at line 542 of file Mapping.cpp.

545 {
546  if (inarray == outarray)
547  {
548  int physTot = m_fields[0]->GetTotPoints();
549  int nvel = m_nConvectiveFields;
550 
551  for (int i = 0; i < nvel; i++)
552  {
553  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
554  }
555  v_ContravarFromCartesian(m_tmp, outarray);
556  }
557  else
558  {
559  v_ContravarFromCartesian(inarray, outarray);
560  }
561 }
virtual GLOBAL_MAPPING_EXPORT void v_ContravarFromCartesian(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=0
Array< OneD, Array< OneD, NekDouble > > m_tmp
Definition: Mapping.h:437
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1255

References m_fields, m_nConvectiveFields, m_tmp, v_ContravarFromCartesian(), and Vmath::Vcopy().

Referenced by v_UpdateBCs().

◆ ContravarToCartesian()

void Nektar::GlobalMapping::Mapping::ContravarToCartesian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)

Convert a contravariant vector to the Cartesian system.

This function converts a contravariant vector in transformed space \(v^{j}\) to the corresponding \(\bar{v}^{i}\) in Cartesian (physical) space using the relation

\[\bar{v}^{i} = \frac{\partial \bar{x}^i}{\partial x^j}v^{j}\]

Parameters
inarrayComponents of the vector in transformed space ( \(v^{j}\))
outarrayComponents of the vector in Cartesian space ( \(\bar{v}^{i}\))

Definition at line 480 of file Mapping.cpp.

483 {
484  if (inarray == outarray)
485  {
486  int physTot = m_fields[0]->GetTotPoints();
487  int nvel = m_nConvectiveFields;
488 
489  for (int i = 0; i < nvel; i++)
490  {
491  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
492  }
493  v_ContravarToCartesian(m_tmp, outarray);
494  }
495  else
496  {
497  v_ContravarToCartesian(inarray, outarray);
498  }
499 }
virtual GLOBAL_MAPPING_EXPORT void v_ContravarToCartesian(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=0

References m_fields, m_nConvectiveFields, m_tmp, v_ContravarToCartesian(), and Vmath::Vcopy().

◆ CovarFromCartesian()

void Nektar::GlobalMapping::Mapping::CovarFromCartesian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)

Convert a covariant vector to the transformed system.

This function converts a covariant vector in Cartesian space \(\bar{v}_{i}\) to the corresponding \(v_{j}\) in transformed space using the relation

\[\bar{v}_{j} = \frac{\partial \bar{x}^i}{\partial x^j}v_{i}\]

Parameters
inarrayComponents of the vector in Cartesian space ( \(\bar{v}_{i}\))
outarrayComponents of the vector in transformed space ( \(v_{j}\))

Definition at line 573 of file Mapping.cpp.

576 {
577  if (inarray == outarray)
578  {
579  int physTot = m_fields[0]->GetTotPoints();
580  int nvel = m_nConvectiveFields;
581 
582  for (int i = 0; i < nvel; i++)
583  {
584  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
585  }
586  v_CovarFromCartesian(m_tmp, outarray);
587  }
588  else
589  {
590  v_CovarFromCartesian(inarray, outarray);
591  }
592 }
virtual GLOBAL_MAPPING_EXPORT void v_CovarFromCartesian(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=0

References m_fields, m_nConvectiveFields, m_tmp, v_CovarFromCartesian(), and Vmath::Vcopy().

◆ CovarToCartesian()

void Nektar::GlobalMapping::Mapping::CovarToCartesian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)

Convert a covariant vector to the Cartesian system.

This function converts a covariant vector in transformed space \(v_{j}\) to the corresponding \(\bar{v}_{i}\) in Cartesian (physical) space using the relation

\[\bar{v}_{i} = \frac{\partial x^j}{\partial \bar{x}^i}v_{j}\]

Parameters
inarrayComponents of the vector in transformed space ( \(v_{j}\))
outarrayComponents of the vector in Cartesian space ( \(\bar{v}_{i}\))

Definition at line 511 of file Mapping.cpp.

514 {
515  if (inarray == outarray)
516  {
517  int physTot = m_fields[0]->GetTotPoints();
518  int nvel = m_nConvectiveFields;
519 
520  for (int i = 0; i < nvel; i++)
521  {
522  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
523  }
524  v_CovarToCartesian(m_tmp, outarray);
525  }
526  else
527  {
528  v_CovarToCartesian(inarray, outarray);
529  }
530 }
virtual GLOBAL_MAPPING_EXPORT void v_CovarToCartesian(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)=0

References m_fields, m_nConvectiveFields, m_tmp, v_CovarToCartesian(), and Vmath::Vcopy().

◆ CurlCurlField()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::CurlCurlField ( Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray,
const bool  generalized 
)
inline

CurlCurl calculated on the whole field.

This function can be used to compute both the generalised CurlCurl or the typical (Cartesian) CurlCurl, depending on the flag generalized

Parameters
inarrayContravariant vector \(u^i\)
outarrayCurlCurl of \(u\)
generalizedFlag defining if generalised or typical CurlCurl

Definition at line 324 of file Mapping.h.

327  {
328  v_CurlCurlField(inarray, outarray, generalized);
329  }
virtual GLOBAL_MAPPING_EXPORT void v_CurlCurlField(Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, const bool generalized)
Definition: Mapping.cpp:988

References v_CurlCurlField().

◆ Divergence()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::Divergence ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, NekDouble > &  outarray 
)
inline

Calculate the generalised divergence operator.

This function is used to calculate the generalised divergence of a contravariant vector, defined as

\[ D = u^i_{,i} = \frac{1}{J}\frac{\partial(Ju^i)}{\partial x^i}\]

Parameters
inarrayContravariant vector \(u^i\)
outarrayDivergence of \(u^i\)

Definition at line 266 of file Mapping.h.

269  {
270  v_Divergence(inarray, outarray);
271  }
virtual GLOBAL_MAPPING_EXPORT void v_Divergence(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, NekDouble > &outarray)
Definition: Mapping.cpp:751

References v_Divergence().

◆ DotGradJacobian()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::DotGradJacobian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, NekDouble > &  outarray 
)
inline

Calculate the dot product with the gradient of the Jacobian.

This function calculates the dot product of an array against the gradient of the Jacobian of the Mapping.

Parameters
inarrayInput array
outarrayOutput array

Definition at line 169 of file Mapping.h.

172  {
173  v_DotGradJacobian(inarray, outarray);
174  }
virtual GLOBAL_MAPPING_EXPORT void v_DotGradJacobian(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, NekDouble > &outarray)
Definition: Mapping.cpp:680

References v_DotGradJacobian().

◆ EvaluateFunction()

void Nektar::GlobalMapping::Mapping::EvaluateFunction ( Array< OneD, MultiRegions::ExpListSharedPtr pFields,
LibUtilities::SessionReaderSharedPtr  pSession,
std::string  pFieldName,
Array< OneD, NekDouble > &  pArray,
const std::string &  pFunctionName,
NekDouble  pTime = NekDouble(0) 
)
protected

Definition at line 401 of file Mapping.cpp.

406 {
407  ASSERTL0(pSession->DefinesFunction(pFunctionName),
408  "Function '" + pFunctionName + "' does not exist.");
409 
410  unsigned int nq = pFields[0]->GetNpoints();
411  if (pArray.size() != nq)
412  {
413  pArray = Array<OneD, NekDouble>(nq);
414  }
415 
417  vType = pSession->GetFunctionType(pFunctionName, pFieldName);
419  {
420  Array<OneD, NekDouble> x0(nq);
421  Array<OneD, NekDouble> x1(nq);
422  Array<OneD, NekDouble> x2(nq);
423 
424  pFields[0]->GetCoords(x0, x1, x2);
426  pSession->GetFunction(pFunctionName, pFieldName);
427 
428  ffunc->Evaluate(x0, x1, x2, pTime, pArray);
429  }
430  else if (vType == LibUtilities::eFunctionTypeFile)
431  {
432  std::string filename =
433  pSession->GetFunctionFilename(pFunctionName, pFieldName);
434 
435  std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
436  std::vector<std::vector<NekDouble>> FieldData;
437  Array<OneD, NekDouble> vCoeffs(pFields[0]->GetNcoeffs());
438  Vmath::Zero(vCoeffs.size(), vCoeffs, 1);
439 
441  LibUtilities::FieldIO::CreateForFile(pSession, filename);
442  fld->Import(filename, FieldDef, FieldData);
443 
444  int idx = -1;
445  for (int i = 0; i < FieldDef.size(); ++i)
446  {
447  for (int j = 0; j < FieldDef[i]->m_fields.size(); ++j)
448  {
449  if (FieldDef[i]->m_fields[j] == pFieldName)
450  {
451  idx = j;
452  }
453  }
454 
455  if (idx >= 0)
456  {
457  pFields[0]->ExtractDataToCoeffs(FieldDef[i], FieldData[i],
458  FieldDef[i]->m_fields[idx],
459  vCoeffs);
460  }
461  else
462  {
463  cout << "Field " + pFieldName + " not found." << endl;
464  }
465  }
466  pFields[0]->BwdTrans(vCoeffs, pArray);
467  }
468 }
static std::shared_ptr< FieldIO > CreateForFile(const LibUtilities::SessionReaderSharedPtr session, const std::string &filename)
Construct a FieldIO object for a given input filename.
Definition: FieldIO.cpp:227
std::shared_ptr< FieldIO > FieldIOSharedPtr
Definition: FieldIO.h:301
std::shared_ptr< Equation > EquationSharedPtr
Definition: Equation.h:130
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:492

References ASSERTL0, Nektar::LibUtilities::FieldIO::CreateForFile(), Nektar::LibUtilities::eFunctionTypeExpression, Nektar::LibUtilities::eFunctionTypeFile, m_fields, and Vmath::Zero().

Referenced by v_InitObject(), and v_UpdateMapping().

◆ EvaluateTimeFunction()

void Nektar::GlobalMapping::Mapping::EvaluateTimeFunction ( LibUtilities::SessionReaderSharedPtr  pSession,
std::string  pFieldName,
Array< OneD, NekDouble > &  pArray,
const std::string &  pFunctionName,
NekDouble  pTime = NekDouble(0) 
)
protected

Definition at line 383 of file Mapping.cpp.

387 {
388  ASSERTL0(pSession->DefinesFunction(pFunctionName),
389  "Function '" + pFunctionName + "' does not exist.");
390 
392  pSession->GetFunction(pFunctionName, pFieldName);
393 
394  Array<OneD, NekDouble> x0(1, 0.0);
395  Array<OneD, NekDouble> x1(1, 0.0);
396  Array<OneD, NekDouble> x2(1, 0.0);
397 
398  ffunc->Evaluate(x0, x1, x2, pTime, pArray);
399 }

References ASSERTL0.

◆ GetCartesianCoordinates()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::GetCartesianCoordinates ( Array< OneD, NekDouble > &  out0,
Array< OneD, NekDouble > &  out1,
Array< OneD, NekDouble > &  out2 
)
inline

Get the Cartesian coordinates in the field.

This function is used to obtain the Cartesian coordinates associated withthe Mapping

Parameters
out0Coordinates in the x-direction
out1Coordinates in the y-direction
out2Coordinates in the z-direction

Definition at line 135 of file Mapping.h.

138  {
139  v_GetCartesianCoordinates(out0, out1, out2);
140  }
virtual GLOBAL_MAPPING_EXPORT void v_GetCartesianCoordinates(Array< OneD, NekDouble > &out0, Array< OneD, NekDouble > &out1, Array< OneD, NekDouble > &out2)
Definition: Mapping.cpp:654

References v_GetCartesianCoordinates().

Referenced by v_UpdateBCs().

◆ GetCoordVelocity()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::GetCoordVelocity ( Array< OneD, Array< OneD, NekDouble >> &  outarray)
inline

Obtain the velocity of the coordinates.

This function is used to obtain the velocity of the coordinates associated with the Mapping

Parameters
outarrayVelocity of the coordinates

Definition at line 244 of file Mapping.h.

246  {
247  v_GetCoordVelocity(outarray);
248  }
virtual GLOBAL_MAPPING_EXPORT void v_GetCoordVelocity(Array< OneD, Array< OneD, NekDouble >> &outarray)
Definition: Mapping.cpp:669

References v_GetCoordVelocity().

Referenced by v_UpdateBCs().

◆ GetInvMetricTensor()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::GetInvMetricTensor ( Array< OneD, Array< OneD, NekDouble >> &  outarray)
inline

Get the inverse of metric tensor \(g^{ij}\).

Definition at line 184 of file Mapping.h.

186  {
187  v_GetInvMetricTensor(outarray);
188  }
virtual GLOBAL_MAPPING_EXPORT void v_GetInvMetricTensor(Array< OneD, Array< OneD, NekDouble >> &outarray)=0

References v_GetInvMetricTensor().

Referenced by Nektar::GlobalMapping::MappingXYofXY::CalculateChristoffel(), and v_RaiseIndex().

◆ GetJacobian()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::GetJacobian ( Array< OneD, NekDouble > &  outarray)
inline

Get the Jacobian of the transformation.

This function is used to obtain the Jacobian of the Mapping.

Parameters
outarrayArray containing the Jacobian

Definition at line 155 of file Mapping.h.

156  {
157  v_GetJacobian(outarray);
158  }
virtual GLOBAL_MAPPING_EXPORT void v_GetJacobian(Array< OneD, NekDouble > &outarray)=0

References v_GetJacobian().

Referenced by v_Divergence(), v_DotGradJacobian(), and Nektar::GlobalMapping::MappingXYofXY::v_GetInvMetricTensor().

◆ GetMetricTensor()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::GetMetricTensor ( Array< OneD, Array< OneD, NekDouble >> &  outarray)
inline

Get the metric tensor \(g_{ij}\).

Definition at line 177 of file Mapping.h.

179  {
180  v_GetMetricTensor(outarray);
181  }
virtual GLOBAL_MAPPING_EXPORT void v_GetMetricTensor(Array< OneD, Array< OneD, NekDouble >> &outarray)=0

References v_GetMetricTensor().

Referenced by Nektar::GlobalMapping::MappingXYofXY::CalculateChristoffel(), and v_LowerIndex().

◆ gradgradU()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::gradgradU ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
inline

Second order covariant derivatives of a contravariant vector.

This function computes the second order covariant derivatives of a contravariant vector, resulting in \(u^{i}_{,jk}\)

Parameters
inarrayContravariant vector \(u^i\)
outarraySecond order derivatives \(u^{i}_{,jk}\)

Definition at line 306 of file Mapping.h.

309  {
310  v_gradgradU(inarray, outarray);
311  }
virtual GLOBAL_MAPPING_EXPORT void v_gradgradU(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
Definition: Mapping.cpp:880

References v_gradgradU().

Referenced by v_CurlCurlField().

◆ HasConstantJacobian()

GLOBAL_MAPPING_EXPORT bool Nektar::GlobalMapping::Mapping::HasConstantJacobian ( )
inline

Get flag defining if mapping has constant Jacobian.

Definition at line 362 of file Mapping.h.

363  {
364  return m_constantJacobian;
365  }
bool m_constantJacobian
Flag defining if the Jacobian is constant.
Definition: Mapping.h:423

References m_constantJacobian.

Referenced by v_DotGradJacobian().

◆ InitObject()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::InitObject ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const TiXmlElement *  pMapping 
)
inline

Initialise the mapping object.

Definition at line 77 of file Mapping.h.

80  {
81  v_InitObject(pFields, pMapping);
82  }
virtual GLOBAL_MAPPING_EXPORT void v_InitObject(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
Definition: Mapping.cpp:100

References v_InitObject().

Referenced by ReplaceField().

◆ IsDefined()

GLOBAL_MAPPING_EXPORT bool Nektar::GlobalMapping::Mapping::IsDefined ( )
inline

Get flag determining if the mapping was defined or is trivial.

Definition at line 368 of file Mapping.h.

369  {
370  return m_isDefined;
371  }

References m_isDefined.

◆ IsFromFunction()

GLOBAL_MAPPING_EXPORT bool Nektar::GlobalMapping::Mapping::IsFromFunction ( )
inline

Get flag defining if mapping is defined by a function.

Definition at line 350 of file Mapping.h.

351  {
352  return m_fromFunction;
353  }
bool m_fromFunction
Flag defining if the Mapping is defined by a function.
Definition: Mapping.h:427

References m_fromFunction.

◆ IsTimeDependent()

GLOBAL_MAPPING_EXPORT bool Nektar::GlobalMapping::Mapping::IsTimeDependent ( void  )
inline

Get flag defining if mapping is time-dependent.

Definition at line 338 of file Mapping.h.

339  {
340  return m_timeDependent;
341  }
bool m_timeDependent
Flag defining if the Mapping is time-dependent.
Definition: Mapping.h:425

References m_timeDependent.

Referenced by v_UpdateBCs().

◆ Load()

MappingSharedPtr Nektar::GlobalMapping::Mapping::Load ( const LibUtilities::SessionReaderSharedPtr pSession,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields 
)
static

Return a pointer to the mapping, creating it on first call.

This function is responsible for loading the Mapping, guaranteeing that a single instance of this class exists. When it is first called, it creates a Mapping and returns a pointer to it. On subsequent calls, it just returns the pointer.

Parameters
pSessionSession reader object
pFieldsFields which will be used by the Mapping
Returns
Pointer to the Mapping

Definition at line 269 of file Mapping.cpp.

272 {
273  if (!m_init)
274  {
275  TiXmlElement *vMapping = NULL;
276  string vType;
277  if (pSession->DefinesElement("Nektar/Mapping"))
278  {
279  vMapping = pSession->GetElement("Nektar/Mapping");
280  vType = vMapping->Attribute("TYPE");
281  m_isDefined = true;
282  }
283  else
284  {
285  vType = "Translation";
286  }
287 
288  m_mappingPtr = GetMappingFactory().CreateInstance(vType, pSession,
289  pFields, vMapping);
290 
291  m_init = true;
292  }
293 
294  return m_mappingPtr;
295 }
static MappingSharedPtr m_mappingPtr
Definition: Mapping.h:430
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:144
MappingFactory & GetMappingFactory()
Declaration of the mapping factory singleton.
Definition: Mapping.cpp:52

References Nektar::LibUtilities::NekFactory< tKey, tBase, tParam >::CreateInstance(), Nektar::GlobalMapping::GetMappingFactory(), m_init, m_isDefined, and m_mappingPtr.

Referenced by Nektar::FieldUtils::ProcessMapping::GetMapping(), Nektar::MappingExtrapolate::MappingExtrapolate(), Nektar::SolverUtils::DriverAdaptive::v_Execute(), Nektar::SolverUtils::FilterAeroForces::v_Initialise(), Nektar::VCSMapping::v_InitObject(), Nektar::ForcingMovingBody::v_InitObject(), and Nektar::SolverUtils::EquationSystem::WriteFld().

◆ LowerIndex()

void Nektar::GlobalMapping::Mapping::LowerIndex ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)

Lower index of vector: \(v_{i} = g_{ij}*v^{j}\).

This function lowers the index of the contravariant vector \(v^{i}\), transforming it into its associated covariant vector \(v_{j}\) according to the relation

\[v_{j} = g_{ij}v^{i}\]

where \(g_{ij}\) is the metric tensor.

Parameters
inarrayComponents of the contravariant vector \(v^{i}\)
outarrayComponents of the covariant vector \(v_{j}\)

Definition at line 604 of file Mapping.cpp.

606 {
607  if (inarray == outarray)
608  {
609  int physTot = m_fields[0]->GetTotPoints();
610  int nvel = m_nConvectiveFields;
611 
612  for (int i = 0; i < nvel; i++)
613  {
614  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
615  }
616  v_LowerIndex(m_tmp, outarray);
617  }
618  else
619  {
620  v_LowerIndex(inarray, outarray);
621  }
622 }
virtual GLOBAL_MAPPING_EXPORT void v_LowerIndex(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
Definition: Mapping.cpp:709

References m_fields, m_nConvectiveFields, m_tmp, v_LowerIndex(), and Vmath::Vcopy().

◆ Output()

void Nektar::GlobalMapping::Mapping::Output ( LibUtilities::FieldMetaDataMap fieldMetaDataMap,
const std::string &  outname 
)

Output function called when a chk or fld file is written.

This function should be called before writing a fld or chk file. It updates the metadata with information from the Mapping. Also, if the mapping is not defined by a function, it writes a .map file containing the coordinates and velocity of the coordinates.

Parameters
fieldMetaDataMapMetadata of the output file
outnameName of the output file

Definition at line 305 of file Mapping.cpp.

307 {
308  // Only do anything if mapping exists
309  if (m_isDefined)
310  {
311  fieldMetaDataMap["MappingCartesianVel"] = std::string("False");
312  if (m_fromFunction)
313  {
314  // Add metadata
315  fieldMetaDataMap["MappingType"] = std::string("Expression");
316  fieldMetaDataMap["MappingExpression"] = m_funcName;
317  if (m_timeDependent)
318  {
319  fieldMetaDataMap["MappingVelExpression"] = m_velFuncName;
320  }
321  }
322  else
323  {
324  int expdim = m_fields[0]->GetGraph()->GetMeshDimension();
325  string fieldNames[3] = {"x", "y", "z"};
326  string velFieldNames[3] = {"vx", "vy", "vz"};
327 
328  std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef =
329  m_fields[0]->GetFieldDefinitions();
330  std::vector<std::vector<NekDouble>> FieldData(FieldDef.size());
331 
332  int ncoeffs = m_fields[0]->GetNcoeffs();
333  Array<OneD, NekDouble> fieldcoeffs(ncoeffs);
334 
335  bool wavespace = m_fields[0]->GetWaveSpace();
336  m_fields[0]->SetWaveSpace(false);
337  // copy coordinates Data into FieldData and set variable
338  for (int j = 0; j < expdim; ++j)
339  {
340  m_fields[0]->FwdTransLocalElmt(m_coords[j], fieldcoeffs);
341 
342  for (int i = 0; i < FieldDef.size(); ++i)
343  {
344  // Could do a search here to find correct variable
345  FieldDef[i]->m_fields.push_back(fieldNames[j]);
346  m_fields[0]->AppendFieldData(FieldDef[i], FieldData[i],
347  fieldcoeffs);
348  }
349  }
350  if (m_timeDependent)
351  {
352  // copy coordinates velocity Data into FieldData and set
353  // variable
354  for (int j = 0; j < expdim; ++j)
355  {
356  m_fields[0]->FwdTrans(m_coordsVel[j], fieldcoeffs);
357 
358  for (int i = 0; i < FieldDef.size(); ++i)
359  {
360  // Could do a search here to find correct variable
361  FieldDef[i]->m_fields.push_back(velFieldNames[j]);
362  m_fields[0]->AppendFieldData(FieldDef[i], FieldData[i],
363  fieldcoeffs);
364  }
365  }
366  }
367 
368  std::string outfile = outname;
369  outfile.erase(outfile.end() - 4, outfile.end());
370  outfile += ".map";
371 
372  m_fld->Write(outfile, FieldDef, FieldData, fieldMetaDataMap);
373 
374  // Write metadata to orginal output
375  fieldMetaDataMap["MappingType"] = std::string("File");
376  fieldMetaDataMap["FileName"] = outfile;
377 
378  m_fields[0]->SetWaveSpace(wavespace);
379  }
380  }
381 }
Array< OneD, Array< OneD, NekDouble > > m_coords
Array with the Cartesian coordinates.
Definition: Mapping.h:408
Array< OneD, Array< OneD, NekDouble > > m_coordsVel
Array with the velocity of the coordinates.
Definition: Mapping.h:410
std::string m_velFuncName
Name of the function containing the velocity of the coordinates.
Definition: Mapping.h:419
std::string m_funcName
Name of the function containing the coordinates.
Definition: Mapping.h:417

References m_coords, m_coordsVel, m_fields, m_fld, m_fromFunction, m_funcName, m_isDefined, m_timeDependent, and m_velFuncName.

◆ RaiseIndex()

void Nektar::GlobalMapping::Mapping::RaiseIndex ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)

Raise index of vector: \(v^{i} = g^{ij}*v_{j}\).

This function raises the index of the covariant vector \(v_{j}\), transforming it into its associated contravariant vector \(v^{i}\) according to the relation

\[v^{i} = g^{ij}v_{j}\]

where \(g^{ij}\) is the inverse of the metric tensor.

Parameters
inarrayComponents of the contravariant vector \(v^{i}\)
outarrayComponents of the covariant vector \(v_{j}\)

Definition at line 634 of file Mapping.cpp.

636 {
637  if (inarray == outarray)
638  {
639  int physTot = m_fields[0]->GetTotPoints();
640  int nvel = m_nConvectiveFields;
641 
642  for (int i = 0; i < nvel; i++)
643  {
644  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
645  }
646  v_RaiseIndex(m_tmp, outarray);
647  }
648  else
649  {
650  v_RaiseIndex(inarray, outarray);
651  }
652 }
virtual GLOBAL_MAPPING_EXPORT void v_RaiseIndex(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
Definition: Mapping.cpp:730

References m_fields, m_nConvectiveFields, m_tmp, v_RaiseIndex(), and Vmath::Vcopy().

Referenced by v_CurlCurlField(), and v_VelocityLaplacian().

◆ ReplaceField()

void Nektar::GlobalMapping::Mapping::ReplaceField ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields)

Replace the Expansion List used by the mapping.

This function replaces the expansion list contained in m_fields, and then proceeds reinitialising the Mapping with this new field.

Parameters
pFieldsNew field to be used by the Mapping

Definition at line 246 of file Mapping.cpp.

248 {
249  m_fields = pFields;
250 
251  TiXmlElement *vMapping = NULL;
252 
253  if (m_session->DefinesElement("Nektar/Mapping"))
254  {
255  vMapping = m_session->GetElement("Nektar/Mapping");
256  }
257  InitObject(pFields, vMapping);
258 }
GLOBAL_MAPPING_EXPORT void InitObject(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
Initialise the mapping object.
Definition: Mapping.h:77

References InitObject(), m_fields, and m_session.

◆ SetFromFunction()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::SetFromFunction ( const bool  value)
inline

Set flag defining if mapping is defined by a function.

Definition at line 356 of file Mapping.h.

357  {
358  m_fromFunction = value;
359  }

References m_fromFunction.

◆ SetTimeDependent()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::SetTimeDependent ( const bool  value)
inline

Set flag defining if mapping is time-dependent.

Definition at line 344 of file Mapping.h.

345  {
346  m_timeDependent = value;
347  }

References m_timeDependent.

◆ UpdateBCs()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::UpdateBCs ( const NekDouble  time)
inline

Update the Dirichlet Boundary Conditions when using Mappings.

Definition at line 378 of file Mapping.h.

379  {
380  v_UpdateBCs(time);
381  }
virtual GLOBAL_MAPPING_EXPORT void v_UpdateBCs(const NekDouble time)
Definition: Mapping.cpp:1050

References v_UpdateBCs().

◆ UpdateGeomInfo()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::UpdateGeomInfo ( )
inline

Recompute the metric terms of the Mapping.

Definition at line 395 of file Mapping.h.

396  {
398  }
virtual GLOBAL_MAPPING_EXPORT void v_UpdateGeomInfo()=0

References v_UpdateGeomInfo().

Referenced by v_InitObject(), and v_UpdateMapping().

◆ UpdateMapping()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::UpdateMapping ( const NekDouble  time,
const Array< OneD, Array< OneD, NekDouble >> &  coords = NullNekDoubleArrayOfArray,
const Array< OneD, Array< OneD, NekDouble >> &  coordsVel = NullNekDoubleArrayOfArray 
)
inline

Update the Mapping with new coordinates.

Definition at line 384 of file Mapping.h.

390  {
391  v_UpdateMapping(time, coords, coordsVel);
392  }
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)
Definition: Mapping.cpp:1206

References v_UpdateMapping().

◆ v_ApplyChristoffelContravar()

virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_ApplyChristoffelContravar ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
protectedpure virtual

◆ v_ApplyChristoffelCovar()

virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_ApplyChristoffelCovar ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
protectedpure virtual

◆ v_ContravarFromCartesian()

virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_ContravarFromCartesian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
protectedpure virtual

◆ v_ContravarToCartesian()

virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_ContravarToCartesian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
protectedpure virtual

◆ v_CovarFromCartesian()

virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_CovarFromCartesian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
protectedpure virtual

◆ v_CovarToCartesian()

virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_CovarToCartesian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
protectedpure virtual

◆ v_CurlCurlField()

void Nektar::GlobalMapping::Mapping::v_CurlCurlField ( Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray,
const bool  generalized 
)
protectedvirtual

Definition at line 988 of file Mapping.cpp.

991 {
992  int physTot = m_fields[0]->GetTotPoints();
993  int nvel = m_nConvectiveFields;
994 
995  // Set wavespace to false and store current value
996  bool wavespace = m_fields[0]->GetWaveSpace();
997  m_fields[0]->SetWaveSpace(false);
998 
999  // For implicit treatment of viscous terms, we want the generalized
1000  // curlcurl and for explicit treatment, we want the cartesian one.
1001  if (generalized)
1002  {
1003  // Get the second derivatives u^i_{,jk}
1004  Array<OneD, Array<OneD, NekDouble>> tmp(nvel);
1005  Array<OneD, Array<OneD, NekDouble>> ddU(nvel * nvel * nvel);
1006  gradgradU(inarray, ddU);
1007 
1008  // Raise index to obtain A^{ip}_{k} = g^pj u^i_{,jk}
1009  for (int i = 0; i < nvel; ++i)
1010  {
1011  for (int k = 0; k < nvel; ++k)
1012  {
1013  // Copy to wk
1014  for (int j = 0; j < nvel; ++j)
1015  {
1016  tmp[j] = ddU[i * nvel * nvel + j * nvel + k];
1017  }
1018  RaiseIndex(tmp, m_tmp);
1019  for (int p = 0; p < nvel; ++p)
1020  {
1021  Vmath::Vcopy(physTot, m_tmp[p], 1,
1022  ddU[i * nvel * nvel + p * nvel + k], 1);
1023  }
1024  }
1025  }
1026  // The curlcurl is g^ji u^k_{kj} - g^jk u^i_kj = A^{ki}_k - A^{ik}_k
1027  for (int i = 0; i < nvel; ++i)
1028  {
1029  outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
1030  for (int k = 0; k < nvel; ++k)
1031  {
1032  Vmath::Vadd(physTot, outarray[i], 1,
1033  ddU[k * nvel * nvel + i * nvel + k], 1, outarray[i],
1034  1);
1035  Vmath::Vsub(physTot, outarray[i], 1,
1036  ddU[i * nvel * nvel + k * nvel + k], 1, outarray[i],
1037  1);
1038  }
1039  }
1040  }
1041  else
1042  {
1043  m_fields[0]->CurlCurl(inarray, outarray);
1044  }
1045 
1046  // Restore value of wavespace
1047  m_fields[0]->SetWaveSpace(wavespace);
1048 }
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.
Definition: Mapping.h:306
GLOBAL_MAPPING_EXPORT void RaiseIndex(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
Raise index of vector: .
Definition: Mapping.cpp:634
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:359
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:419

References gradgradU(), m_fields, m_nConvectiveFields, m_tmp, CellMLToNektar.cellml_metadata::p, RaiseIndex(), Vmath::Vadd(), Vmath::Vcopy(), and Vmath::Vsub().

Referenced by CurlCurlField().

◆ v_Divergence()

void Nektar::GlobalMapping::Mapping::v_Divergence ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Definition at line 751 of file Mapping.cpp.

753 {
754  int physTot = m_fields[0]->GetTotPoints();
755  Array<OneD, NekDouble> wk(physTot, 0.0);
756 
757  Vmath::Zero(physTot, outarray, 1);
758 
759  // Set wavespace to false and store current value
760  bool wavespace = m_fields[0]->GetWaveSpace();
761  m_fields[0]->SetWaveSpace(false);
762 
763  // Get Mapping Jacobian
764  Array<OneD, NekDouble> Jac(physTot, 0.0);
765  GetJacobian(Jac);
766 
767  for (int i = 0; i < m_nConvectiveFields; ++i)
768  {
769  Vmath::Vmul(physTot, Jac, 1, inarray[i], 1, wk, 1); // J*Ui
770  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[i], wk,
771  wk); // (J*Ui)_i
772  Vmath::Vadd(physTot, wk, 1, outarray, 1, outarray, 1);
773  }
774  Vmath::Vdiv(physTot, outarray, 1, Jac, 1, outarray, 1); // 1/J*(J*Ui)_i
775 
776  m_fields[0]->SetWaveSpace(wavespace);
777 }
GLOBAL_MAPPING_EXPORT void GetJacobian(Array< OneD, NekDouble > &outarray)
Get the Jacobian of the transformation.
Definition: Mapping.h:155
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:89
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:209
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:284

References Nektar::MultiRegions::DirCartesianMap, GetJacobian(), m_fields, m_nConvectiveFields, Vmath::Vadd(), Vmath::Vdiv(), Vmath::Vmul(), and Vmath::Zero().

Referenced by Divergence().

◆ v_DotGradJacobian()

void Nektar::GlobalMapping::Mapping::v_DotGradJacobian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Reimplemented in Nektar::GlobalMapping::MappingXYofZ, Nektar::GlobalMapping::MappingXofZ, Nektar::GlobalMapping::MappingXofXZ, and Nektar::GlobalMapping::MappingTranslation.

Definition at line 680 of file Mapping.cpp.

683 {
684  int physTot = m_fields[0]->GetTotPoints();
685 
686  outarray = Array<OneD, NekDouble>(physTot, 0.0);
687  if (!HasConstantJacobian())
688  {
689  // Set wavespace to false and store current value
690  bool wavespace = m_fields[0]->GetWaveSpace();
691  m_fields[0]->SetWaveSpace(false);
692 
693  // Get Mapping Jacobian
694  Array<OneD, NekDouble> Jac(physTot, 0.0);
695  GetJacobian(Jac);
696 
697  // Calculate inarray . grad(Jac)
698  Array<OneD, NekDouble> wk(physTot, 0.0);
699  for (int i = 0; i < m_nConvectiveFields; ++i)
700  {
701  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[i], Jac, wk);
702  Vmath::Vvtvp(physTot, inarray[i], 1, wk, 1, outarray, 1, outarray,
703  1);
704  }
705  m_fields[0]->SetWaveSpace(wavespace);
706  }
707 }
GLOBAL_MAPPING_EXPORT bool HasConstantJacobian()
Get flag defining if mapping has constant Jacobian.
Definition: Mapping.h:362
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:574

References Nektar::MultiRegions::DirCartesianMap, GetJacobian(), HasConstantJacobian(), m_fields, m_nConvectiveFields, and Vmath::Vvtvp().

Referenced by DotGradJacobian().

◆ v_GetCartesianCoordinates()

void Nektar::GlobalMapping::Mapping::v_GetCartesianCoordinates ( Array< OneD, NekDouble > &  out0,
Array< OneD, NekDouble > &  out1,
Array< OneD, NekDouble > &  out2 
)
protectedvirtual

Definition at line 654 of file Mapping.cpp.

657 {
658  int physTot = m_fields[0]->GetTotPoints();
659 
660  out0 = Array<OneD, NekDouble>(physTot, 0.0);
661  out1 = Array<OneD, NekDouble>(physTot, 0.0);
662  out2 = Array<OneD, NekDouble>(physTot, 0.0);
663 
664  Vmath::Vcopy(physTot, m_coords[0], 1, out0, 1);
665  Vmath::Vcopy(physTot, m_coords[1], 1, out1, 1);
666  Vmath::Vcopy(physTot, m_coords[2], 1, out2, 1);
667 }

References m_coords, m_fields, and Vmath::Vcopy().

Referenced by GetCartesianCoordinates().

◆ v_GetCoordVelocity()

void Nektar::GlobalMapping::Mapping::v_GetCoordVelocity ( Array< OneD, Array< OneD, NekDouble >> &  outarray)
protectedvirtual

Definition at line 669 of file Mapping.cpp.

670 {
671  int physTot = m_fields[0]->GetTotPoints();
672 
673  for (int i = 0; i < m_nConvectiveFields; ++i)
674  {
675  outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
676  Vmath::Vcopy(physTot, m_coordsVel[i], 1, outarray[i], 1);
677  }
678 }

References m_coordsVel, m_fields, m_nConvectiveFields, and Vmath::Vcopy().

Referenced by GetCoordVelocity().

◆ v_GetInvMetricTensor()

virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_GetInvMetricTensor ( Array< OneD, Array< OneD, NekDouble >> &  outarray)
protectedpure virtual

◆ v_GetJacobian()

virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_GetJacobian ( Array< OneD, NekDouble > &  outarray)
protectedpure virtual

◆ v_GetMetricTensor()

virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_GetMetricTensor ( Array< OneD, Array< OneD, NekDouble >> &  outarray)
protectedpure virtual

◆ v_gradgradU()

void Nektar::GlobalMapping::Mapping::v_gradgradU ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
protectedvirtual

Definition at line 880 of file Mapping.cpp.

882 {
883  int physTot = m_fields[0]->GetTotPoints();
884  int nvel = m_nConvectiveFields;
885 
886  // Declare variables
887  outarray = Array<OneD, Array<OneD, NekDouble>>(nvel * nvel * nvel);
888  Array<OneD, Array<OneD, NekDouble>> tmp(nvel);
889  for (int i = 0; i < nvel * nvel * nvel; i++)
890  {
891  outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
892  }
893 
894  // Set wavespace to false and store current value
895  bool wavespace = m_fields[0]->GetWaveSpace();
896  m_fields[0]->SetWaveSpace(false);
897 
898  // Calculate vector gradient u^i_(,j) = du^i/dx^j + {i,pj}*u^p
900  for (int i = 0; i < nvel; i++)
901  {
902  if (nvel == 2)
903  {
904  m_fields[0]->PhysDeriv(inarray[i], m_wk2[i * nvel + 0],
905  m_wk2[i * nvel + 1]);
906  }
907  else
908  {
909  m_fields[0]->PhysDeriv(inarray[i], m_wk2[i * nvel + 0],
910  m_wk2[i * nvel + 1], m_wk2[i * nvel + 2]);
911  }
912  for (int j = 0; j < nvel; j++)
913  {
914  Vmath::Vadd(physTot, m_wk1[i * nvel + j], 1, m_wk2[i * nvel + j], 1,
915  m_wk1[i * nvel + j], 1);
916  }
917  }
918 
919  //
920  // Calculate (u^i_,j),k
921  //
922 
923  // Step 1 : d(u^i_,j))/d(x^k)
924  for (int i = 0; i < nvel; i++)
925  {
926  for (int j = 0; j < nvel; j++)
927  {
928  if (nvel == 2)
929  {
930  m_fields[0]->PhysDeriv(
931  m_wk1[i * nvel + j],
932  outarray[i * nvel * nvel + j * nvel + 0],
933  outarray[i * nvel * nvel + j * nvel + 1]);
934  }
935  else
936  {
937  m_fields[0]->PhysDeriv(
938  m_wk1[i * nvel + j],
939  outarray[i * nvel * nvel + j * nvel + 0],
940  outarray[i * nvel * nvel + j * nvel + 1],
941  outarray[i * nvel * nvel + j * nvel + 2]);
942  }
943  }
944  }
945 
946  // Step 2: d(u^i_,j)/d(x^k) - {p,jk}*u^i_,p
947  for (int i = 0; i < nvel; i++)
948  {
949  for (int p = 0; p < nvel; p++)
950  {
951  tmp[p] = m_wk1[i * nvel + p];
952  }
954  for (int j = 0; j < nvel; j++)
955  {
956  for (int k = 0; k < nvel; k++)
957  {
958  Vmath::Vsub(physTot, outarray[i * nvel * nvel + j * nvel + k],
959  1, m_wk2[j * nvel + k], 1,
960  outarray[i * nvel * nvel + j * nvel + k], 1);
961  }
962  }
963  }
964 
965  // Step 3: d(u^i_,j)/d(x^k) - {p,jk}*u^i_,p + {i,pk} u^p_,j
966  for (int j = 0; j < nvel; j++)
967  {
968  for (int p = 0; p < nvel; p++)
969  {
970  tmp[p] = m_wk1[p * nvel + j];
971  }
973  for (int i = 0; i < nvel; i++)
974  {
975  for (int k = 0; k < nvel; k++)
976  {
977  Vmath::Vadd(physTot, outarray[i * nvel * nvel + j * nvel + k],
978  1, m_wk2[i * nvel + k], 1,
979  outarray[i * nvel * nvel + j * nvel + k], 1);
980  }
981  }
982  }
983 
984  // Restore value of wavespace
985  m_fields[0]->SetWaveSpace(wavespace);
986 }
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.
Definition: Mapping.h:229
Array< OneD, Array< OneD, NekDouble > > m_wk2
Definition: Mapping.h:436
Array< OneD, Array< OneD, NekDouble > > m_wk1
Definition: Mapping.h:435
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.
Definition: Mapping.h:211

References ApplyChristoffelContravar(), ApplyChristoffelCovar(), m_fields, m_nConvectiveFields, m_wk1, m_wk2, CellMLToNektar.cellml_metadata::p, Vmath::Vadd(), and Vmath::Vsub().

Referenced by gradgradU().

◆ v_InitObject()

void Nektar::GlobalMapping::Mapping::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 in Nektar::GlobalMapping::MappingXYofZ, Nektar::GlobalMapping::MappingXYofXY, Nektar::GlobalMapping::MappingXofZ, Nektar::GlobalMapping::MappingXofXZ, Nektar::GlobalMapping::MappingTranslation, and Nektar::GlobalMapping::MappingGeneral.

Definition at line 100 of file Mapping.cpp.

103 {
104  boost::ignore_unused(pFields);
105 
106  int phystot = m_fields[0]->GetTotPoints();
107  m_fromFunction = true;
108  // Initialise variables
109  m_coords = Array<OneD, Array<OneD, NekDouble>>(3);
110  m_coordsVel = Array<OneD, Array<OneD, NekDouble>>(3);
111  Array<OneD, Array<OneD, NekDouble>> coords(3);
112  for (int i = 0; i < 3; i++)
113  {
114  m_coords[i] = Array<OneD, NekDouble>(phystot);
115  m_coordsVel[i] = Array<OneD, NekDouble>(phystot);
116  coords[i] = Array<OneD, NekDouble>(phystot);
117  }
118 
119  // Check if mapping is defined as time-dependent
120  const TiXmlElement *timeDep = pMapping->FirstChildElement("TIMEDEPENDENT");
121  if (timeDep)
122  {
123  string sTimeDep = timeDep->GetText();
124  m_timeDependent = (boost::iequals(sTimeDep, "true")) ||
125  (boost::iequals(sTimeDep, "yes"));
126  }
127  else
128  {
129  m_timeDependent = false;
130  }
131 
132  // Load coordinates
133  string fieldNames[3] = {"x", "y", "z"};
134  const TiXmlElement *funcNameElmt = pMapping->FirstChildElement("COORDS");
135  if (funcNameElmt)
136  {
137  m_funcName = funcNameElmt->GetText();
138  ASSERTL0(m_session->DefinesFunction(m_funcName),
139  "Function '" + m_funcName + "' not defined.");
140 
141  // Get coordinates in the domain
142  m_fields[0]->GetCoords(coords[0], coords[1], coords[2]);
143 
144  std::string s_FieldStr;
145  // Check if function from session file defines each component
146  // and evaluate them, otherwise use trivial transformation
147  for (int i = 0; i < 3; i++)
148  {
149  s_FieldStr = fieldNames[i];
150  if (m_session->DefinesFunction(m_funcName, s_FieldStr))
151  {
152  EvaluateFunction(m_fields, m_session, s_FieldStr, m_coords[i],
153  m_funcName);
154  if (i == 2 && m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
155  {
156  ASSERTL0(
157  false,
158  "3DH1D does not support mapping in the z-direction.");
159  }
160  }
161  else
162  {
163  // This coordinate is not defined, so use (x^i)' = x^i
164  Vmath::Vcopy(phystot, coords[i], 1, m_coords[i], 1);
165  }
166  }
167  }
168  else
169  {
170  m_fields[0]->GetCoords(coords[0], coords[1], coords[2]);
171  for (int i = 0; i < 3; i++)
172  {
173  // Use (x^i)' = x^i as default. This can be useful if we
174  // have a time-dependent mapping, and then only the
175  // initial mapping will be trivial
176  Vmath::Vcopy(phystot, coords[i], 1, m_coords[i], 1);
177  }
178  }
179 
180  // Load coordinate velocity if they are defined,
181  // otherwise use zero to make it general
182  string velFieldNames[3] = {"vx", "vy", "vz"};
183  const TiXmlElement *velFuncNameElmt = pMapping->FirstChildElement("VEL");
184  if (velFuncNameElmt)
185  {
186  m_velFuncName = velFuncNameElmt->GetText();
187  ASSERTL0(m_session->DefinesFunction(m_velFuncName),
188  "Function '" + m_velFuncName + "' not defined.");
189 
190  std::string s_FieldStr;
191  // Check if function from session file defines each component
192  // and evaluate them, otherwise use 0
193  for (int i = 0; i < 3; i++)
194  {
195  s_FieldStr = velFieldNames[i];
196  if (m_session->DefinesFunction(m_velFuncName, s_FieldStr))
197  {
198  EvaluateFunction(m_fields, m_session, s_FieldStr,
200  if (i == 2 && m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
201  {
202  ASSERTL0(
203  false,
204  "3DH1D does not support mapping in the z-direction.");
205  }
206  }
207  else
208  {
209  // This coordinate velocity is not defined, so use 0
210  Vmath::Zero(phystot, m_coordsVel[i], 1);
211  }
212  }
213  }
214  else
215  {
216  for (int i = 0; i < 3; i++)
217  {
218  Vmath::Zero(phystot, m_coordsVel[i], 1);
219  }
220  }
221 
222  // Initialise workspace variables
223  int nvel = m_nConvectiveFields;
224  m_wk1 = Array<OneD, Array<OneD, NekDouble>>(nvel * nvel);
225  m_wk2 = Array<OneD, Array<OneD, NekDouble>>(nvel * nvel);
226  m_tmp = Array<OneD, Array<OneD, NekDouble>>(nvel);
227  for (int i = 0; i < nvel; i++)
228  {
229  m_tmp[i] = Array<OneD, NekDouble>(phystot, 0.0);
230  for (int j = 0; j < nvel; j++)
231  {
232  m_wk1[i * nvel + j] = Array<OneD, NekDouble>(phystot, 0.0);
233  m_wk2[i * nvel + j] = Array<OneD, NekDouble>(phystot, 0.0);
234  }
235  }
236 
237  // Calculate information required by the particular mapping
238  UpdateGeomInfo();
239 }
GLOBAL_MAPPING_EXPORT void UpdateGeomInfo()
Recompute the metric terms of the Mapping.
Definition: Mapping.h:395
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))
Definition: Mapping.cpp:401

References ASSERTL0, Nektar::MultiRegions::e3DH1D, EvaluateFunction(), m_coords, m_coordsVel, m_fields, m_fromFunction, m_funcName, m_nConvectiveFields, m_session, m_timeDependent, m_tmp, m_velFuncName, m_wk1, m_wk2, UpdateGeomInfo(), Vmath::Vcopy(), and Vmath::Zero().

Referenced by InitObject(), Nektar::GlobalMapping::MappingGeneral::v_InitObject(), Nektar::GlobalMapping::MappingTranslation::v_InitObject(), Nektar::GlobalMapping::MappingXofXZ::v_InitObject(), Nektar::GlobalMapping::MappingXofZ::v_InitObject(), Nektar::GlobalMapping::MappingXYofXY::v_InitObject(), and Nektar::GlobalMapping::MappingXYofZ::v_InitObject().

◆ v_LowerIndex()

void Nektar::GlobalMapping::Mapping::v_LowerIndex ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
protectedvirtual

Reimplemented in Nektar::GlobalMapping::MappingXofZ, Nektar::GlobalMapping::MappingXofXZ, and Nektar::GlobalMapping::MappingTranslation.

Definition at line 709 of file Mapping.cpp.

711 {
712  int physTot = m_fields[0]->GetTotPoints();
713  int nvel = m_nConvectiveFields;
714 
715  Array<OneD, Array<OneD, NekDouble>> g(nvel * nvel);
716 
717  GetMetricTensor(g);
718 
719  for (int i = 0; i < nvel; i++)
720  {
721  outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
722  for (int j = 0; j < nvel; j++)
723  {
724  Vmath::Vvtvp(physTot, g[i * nvel + j], 1, inarray[j], 1,
725  outarray[i], 1, outarray[i], 1);
726  }
727  }
728 }
GLOBAL_MAPPING_EXPORT void GetMetricTensor(Array< OneD, Array< OneD, NekDouble >> &outarray)
Get the metric tensor .
Definition: Mapping.h:177

References GetMetricTensor(), m_fields, m_nConvectiveFields, and Vmath::Vvtvp().

Referenced by LowerIndex().

◆ v_RaiseIndex()

void Nektar::GlobalMapping::Mapping::v_RaiseIndex ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray 
)
protectedvirtual

Reimplemented in Nektar::GlobalMapping::MappingXofZ, Nektar::GlobalMapping::MappingXofXZ, and Nektar::GlobalMapping::MappingTranslation.

Definition at line 730 of file Mapping.cpp.

732 {
733  int physTot = m_fields[0]->GetTotPoints();
734  int nvel = m_nConvectiveFields;
735 
736  Array<OneD, Array<OneD, NekDouble>> g(nvel * nvel);
737 
739 
740  for (int i = 0; i < nvel; i++)
741  {
742  outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
743  for (int j = 0; j < nvel; j++)
744  {
745  Vmath::Vvtvp(physTot, g[i * nvel + j], 1, inarray[j], 1,
746  outarray[i], 1, outarray[i], 1);
747  }
748  }
749 }
GLOBAL_MAPPING_EXPORT void GetInvMetricTensor(Array< OneD, Array< OneD, NekDouble >> &outarray)
Get the inverse of metric tensor .
Definition: Mapping.h:184

References GetInvMetricTensor(), m_fields, m_nConvectiveFields, and Vmath::Vvtvp().

Referenced by RaiseIndex().

◆ v_UpdateBCs()

void Nektar::GlobalMapping::Mapping::v_UpdateBCs ( const NekDouble  time)
protectedvirtual

Definition at line 1050 of file Mapping.cpp.

1051 {
1052  int physTot = m_fields[0]->GetTotPoints();
1053  int nvel = m_nConvectiveFields;
1054  int nfields = m_fields.size();
1055  int nbnds = m_fields[0]->GetBndConditions().size();
1056 
1057  // Declare variables
1058  Array<OneD, const SpatialDomains::BoundaryConditionShPtr> BndConds;
1059  Array<OneD, MultiRegions::ExpListSharedPtr> BndExp;
1060 
1061  Array<OneD, bool> isDirichlet(nfields);
1062 
1063  Array<OneD, Array<OneD, NekDouble>> values(nfields);
1064  for (int i = 0; i < nfields; i++)
1065  {
1066  values[i] = Array<OneD, NekDouble>(physTot, 0.0);
1067  }
1068 
1069  Array<OneD, Array<OneD, NekDouble>> tmp(nvel);
1070  Array<OneD, Array<OneD, NekDouble>> tmp2(nvel);
1071  Array<OneD, Array<OneD, NekDouble>> coordVel(nvel);
1072  for (int i = 0; i < nvel; i++)
1073  {
1074  tmp[i] = Array<OneD, NekDouble>(physTot, 0.0);
1075  tmp2[i] = Array<OneD, NekDouble>(physTot, 0.0);
1076  coordVel[i] = Array<OneD, NekDouble>(physTot, 0.0);
1077  }
1078 
1079  // Get coordinates velocity in transformed system (for MovingBody regions)
1080  GetCoordVelocity(tmp);
1081  ContravarFromCartesian(tmp, coordVel);
1082 
1083  // Get Cartesian coordinates for evaluating boundary conditions
1084  Array<OneD, Array<OneD, NekDouble>> coords(3);
1085  for (int dir = 0; dir < 3; dir++)
1086  {
1087  coords[dir] = Array<OneD, NekDouble>(physTot, 0.0);
1088  }
1089  GetCartesianCoordinates(coords[0], coords[1], coords[2]);
1090 
1091  // Loop boundary conditions looking for Dirichlet bc's
1092  for (int n = 0; n < nbnds; ++n)
1093  {
1094  // Evaluate original Dirichlet boundary conditions in whole domain
1095  for (int i = 0; i < nfields; ++i)
1096  {
1097  BndConds = m_fields[i]->GetBndConditions();
1098  BndExp = m_fields[i]->GetBndCondExpansions();
1099  if (BndConds[n]->GetBoundaryConditionType() ==
1101  {
1102  isDirichlet[i] = true;
1103  // If we have the a velocity component
1104  // check if all vel bc's are also Dirichlet
1105  if (i < nvel)
1106  {
1107  for (int j = 0; j < nvel; ++j)
1108  {
1109  ASSERTL0(m_fields[j]
1110  ->GetBndConditions()[n]
1111  ->GetBoundaryConditionType() ==
1113  "Mapping only supported when all velocity "
1114  "components have the same type of boundary "
1115  "conditions");
1116  }
1117  }
1118  // Check if bc is time-dependent
1119  ASSERTL0(!BndConds[n]->IsTimeDependent(),
1120  "Time-dependent Dirichlet boundary conditions not "
1121  "supported with mapping yet.");
1122 
1123  // Get boundary condition
1124  LibUtilities::Equation condition =
1125  std::static_pointer_cast<
1126  SpatialDomains::DirichletBoundaryCondition>(BndConds[n])
1127  ->m_dirichletCondition;
1128  // Evaluate
1129  condition.Evaluate(coords[0], coords[1], coords[2], time,
1130  values[i]);
1131  }
1132  else
1133  {
1134  isDirichlet[i] = false;
1135  }
1136  }
1137  // Convert velocity vector to transformed system
1138  if (isDirichlet[0])
1139  {
1140  for (int i = 0; i < nvel; ++i)
1141  {
1142  Vmath::Vcopy(physTot, values[i], 1, tmp[i], 1);
1143  }
1144  ContravarFromCartesian(tmp, tmp2);
1145  for (int i = 0; i < nvel; ++i)
1146  {
1147  Vmath::Vcopy(physTot, tmp2[i], 1, values[i], 1);
1148  }
1149  }
1150 
1151  // Now, project result to boundary
1152  for (int i = 0; i < nfields; ++i)
1153  {
1154  BndConds = m_fields[i]->GetBndConditions();
1155  BndExp = m_fields[i]->GetBndCondExpansions();
1156  if (BndConds[n]->GetUserDefined() == "" ||
1157  BndConds[n]->GetUserDefined() == "MovingBody")
1158  {
1159  m_fields[i]->ExtractPhysToBnd(n, values[i],
1160  BndExp[n]->UpdatePhys());
1161 
1162  // Apply MovingBody correction
1163  if ((i < nvel) && BndConds[n]->GetUserDefined() == "MovingBody")
1164  {
1165  // Get coordinate velocity on boundary
1166  Array<OneD, NekDouble> coordVelBnd(
1167  BndExp[n]->GetTotPoints());
1168  m_fields[i]->ExtractPhysToBnd(n, coordVel[i], coordVelBnd);
1169 
1170  // Apply correction
1171  Vmath::Vadd(BndExp[n]->GetTotPoints(), coordVelBnd, 1,
1172  BndExp[n]->UpdatePhys(), 1,
1173  BndExp[n]->UpdatePhys(), 1);
1174  }
1175  }
1176  }
1177  }
1178 
1179  // Finally, perform FwdTrans in all fields
1180  for (int i = 0; i < m_fields.size(); ++i)
1181  {
1182  // Get boundary condition information
1183  BndConds = m_fields[i]->GetBndConditions();
1184  BndExp = m_fields[i]->GetBndCondExpansions();
1185  for (int n = 0; n < BndConds.size(); ++n)
1186  {
1187  if (BndConds[n]->GetBoundaryConditionType() ==
1189  {
1190  if (BndConds[n]->GetUserDefined() == "" ||
1191  BndConds[n]->GetUserDefined() == "MovingBody")
1192  {
1193  if (m_fields[i]->GetExpType() == MultiRegions::e3DH1D)
1194  {
1195  BndExp[n]->SetWaveSpace(false);
1196  }
1197 
1198  BndExp[n]->FwdTransBndConstrained(
1199  BndExp[n]->GetPhys(), BndExp[n]->UpdateCoeffs());
1200  }
1201  }
1202  }
1203  }
1204 }
GLOBAL_MAPPING_EXPORT void GetCartesianCoordinates(Array< OneD, NekDouble > &out0, Array< OneD, NekDouble > &out1, Array< OneD, NekDouble > &out2)
Get the Cartesian coordinates in the field.
Definition: Mapping.h:135
GLOBAL_MAPPING_EXPORT bool IsTimeDependent()
Get flag defining if mapping is time-dependent.
Definition: Mapping.h:338
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.
Definition: Mapping.cpp:542
GLOBAL_MAPPING_EXPORT void GetCoordVelocity(Array< OneD, Array< OneD, NekDouble >> &outarray)
Obtain the velocity of the coordinates.
Definition: Mapping.h:244

References ASSERTL0, ContravarFromCartesian(), Nektar::MultiRegions::e3DH1D, Nektar::SpatialDomains::eDirichlet, Nektar::LibUtilities::Equation::Evaluate(), GetCartesianCoordinates(), GetCoordVelocity(), IsTimeDependent(), m_fields, m_nConvectiveFields, Vmath::Vadd(), and Vmath::Vcopy().

Referenced by UpdateBCs().

◆ v_UpdateGeomInfo()

virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_UpdateGeomInfo ( )
protectedpure virtual

◆ v_UpdateMapping()

void Nektar::GlobalMapping::Mapping::v_UpdateMapping ( const NekDouble  time,
const Array< OneD, Array< OneD, NekDouble >> &  coords = NullNekDoubleArrayOfArray,
const Array< OneD, Array< OneD, NekDouble >> &  coordsVel = NullNekDoubleArrayOfArray 
)
protectedvirtual

Definition at line 1206 of file Mapping.cpp.

1209 {
1210  if (m_fromFunction)
1211  {
1212  std::string s_FieldStr;
1213  string fieldNames[3] = {"x", "y", "z"};
1214  string velFieldNames[3] = {"vx", "vy", "vz"};
1215  // Check if function from session file defines each component
1216  // and evaluate them, otherwise there is no need to update
1217  // coords
1218  for (int i = 0; i < 3; i++)
1219  {
1220  s_FieldStr = fieldNames[i];
1221  if (m_session->DefinesFunction(m_funcName, s_FieldStr))
1222  {
1223  EvaluateFunction(m_fields, m_session, s_FieldStr, m_coords[i],
1224  m_funcName, time);
1225  if (i == 2 && m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
1226  {
1227  ASSERTL0(
1228  false,
1229  "3DH1D does not support mapping in the z-direction.");
1230  }
1231  }
1232  s_FieldStr = velFieldNames[i];
1233  if (m_session->DefinesFunction(m_velFuncName, s_FieldStr))
1234  {
1235  EvaluateFunction(m_fields, m_session, s_FieldStr,
1236  m_coordsVel[i], m_velFuncName, time);
1237  if (i == 2 && m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
1238  {
1239  ASSERTL0(
1240  false,
1241  "3DH1D does not support mapping in the z-direction.");
1242  }
1243  }
1244  }
1245  }
1246  else
1247  {
1248  int physTot = m_fields[0]->GetTotPoints();
1249  int nvel = m_nConvectiveFields;
1250  // Copy coordinates
1251  for (int i = 0; i < 3; i++)
1252  {
1253  Vmath::Vcopy(physTot, coords[i], 1, m_coords[i], 1);
1254  }
1255 
1256  for (int i = 0; i < nvel; i++)
1257  {
1258  Vmath::Vcopy(physTot, coordsVel[i], 1, m_coordsVel[i], 1);
1259  }
1260  }
1261 
1262  // Update the information required by the specific mapping
1263  UpdateGeomInfo();
1264 }

References ASSERTL0, Nektar::MultiRegions::e3DH1D, EvaluateFunction(), m_coords, m_coordsVel, m_fields, m_fromFunction, m_funcName, m_nConvectiveFields, m_session, m_velFuncName, UpdateGeomInfo(), and Vmath::Vcopy().

Referenced by UpdateMapping().

◆ v_VelocityLaplacian()

void Nektar::GlobalMapping::Mapping::v_VelocityLaplacian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray,
const NekDouble  alpha 
)
protectedvirtual

Definition at line 779 of file Mapping.cpp.

782 {
783  int physTot = m_fields[0]->GetTotPoints();
784  int nvel = m_nConvectiveFields;
785 
786  Array<OneD, Array<OneD, NekDouble>> tmp(nvel);
787 
788  // Set wavespace to false and store current value
789  bool wavespace = m_fields[0]->GetWaveSpace();
790  m_fields[0]->SetWaveSpace(false);
791 
792  // Calculate vector gradient wk2 = u^i_(,k) = du^i/dx^k + {i,jk}*u^j
794  for (int i = 0; i < nvel; i++)
795  {
796  if (nvel == 2)
797  {
798  m_fields[0]->PhysDeriv(inarray[i], m_wk2[i * nvel + 0],
799  m_wk2[i * nvel + 1]);
800  }
801  else
802  {
803  m_fields[0]->PhysDeriv(inarray[i], m_wk2[i * nvel + 0],
804  m_wk2[i * nvel + 1], m_wk2[i * nvel + 2]);
805  }
806  for (int k = 0; k < nvel; k++)
807  {
808  Vmath::Vadd(physTot, m_wk1[i * nvel + k], 1, m_wk2[i * nvel + k], 1,
809  m_wk1[i * nvel + k], 1);
810  }
811  }
812  // Calculate wk1 = A^(ij) = g^(jk)*u^i_(,k)
813  for (int i = 0; i < nvel; i++)
814  {
815  for (int k = 0; k < nvel; k++)
816  {
817  tmp[k] = m_wk1[i * nvel + k];
818  }
819  RaiseIndex(tmp, m_tmp);
820  for (int j = 0; j < nvel; j++)
821  {
822  Vmath::Vcopy(physTot, m_tmp[j], 1, m_wk1[i * nvel + j], 1);
823  }
824  }
825  //
826  // Calculate L(U)^i = (A^(ij))_(,j) - alpha*d^2(u^i)/dx^jdx^j
827  //
828 
829  // Step 1 :
830  // d(A^(ij) - alpha*du^i/dx^j)/d(x^j)
831  for (int i = 0; i < nvel; i++)
832  {
833  outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
834  for (int j = 0; j < nvel; j++)
835  {
836  Vmath::Smul(physTot, alpha, m_wk2[i * nvel + j], 1, m_tmp[0], 1);
837  Vmath::Vsub(physTot, m_wk1[i * nvel + j], 1, m_tmp[0], 1, m_tmp[0],
838  1);
839 
840  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[j], m_tmp[0],
841  m_tmp[1]);
842  Vmath::Vadd(physTot, outarray[i], 1, m_tmp[1], 1, outarray[i], 1);
843  }
844  }
845 
846  // Step 2: d(A^(ij))/d(x^j) + {j,pj}*A^(ip)
847  for (int i = 0; i < nvel; i++)
848  {
849  for (int p = 0; p < nvel; p++)
850  {
851  tmp[p] = m_wk1[i * nvel + p];
852  }
854  for (int j = 0; j < nvel; j++)
855  {
856  Vmath::Vadd(physTot, outarray[i], 1, m_wk2[j * nvel + j], 1,
857  outarray[i], 1);
858  }
859  }
860 
861  // Step 3: d(A^(ij))/d(x^j) + {j,pj}*A^(ip) + {i,pj} A^(pj)
862  for (int j = 0; j < nvel; j++)
863  {
864  for (int p = 0; p < nvel; p++)
865  {
866  tmp[p] = m_wk1[p * nvel + j];
867  }
869  for (int i = 0; i < nvel; i++)
870  {
871  Vmath::Vadd(physTot, outarray[i], 1, m_wk2[i * nvel + j], 1,
872  outarray[i], 1);
873  }
874  }
875 
876  // Restore value of wavespace
877  m_fields[0]->SetWaveSpace(wavespace);
878 }
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.cpp:248

References ApplyChristoffelContravar(), Nektar::MultiRegions::DirCartesianMap, m_fields, m_nConvectiveFields, m_tmp, m_wk1, m_wk2, CellMLToNektar.cellml_metadata::p, RaiseIndex(), Vmath::Smul(), Vmath::Vadd(), Vmath::Vcopy(), and Vmath::Vsub().

Referenced by VelocityLaplacian().

◆ VelocityLaplacian()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::VelocityLaplacian ( const Array< OneD, Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray,
const NekDouble  alpha = 0.0 
)
inline

Generalised (correction to the) velocity Laplacian operator.

This function is used to calculate a correction defined as the difference between the generalised Laplacian and the original Laplacian multiplied by a constant \(\alpha\), resulting in

\[ L^i = g^{jk}u^{i}_{,jk} - \alpha \frac{\partial^2 x^i} {\partial x^j \partial x^j}\]

By default, \(\alpha\) is zero, resulting in the generalised Laplacian.

Parameters
inarrayContravariant vector \(u^i\)
outarrayResult of the operation
alphaThe constant \(\alpha\)

Definition at line 289 of file Mapping.h.

293  {
294  v_VelocityLaplacian(inarray, outarray, alpha);
295  }
virtual GLOBAL_MAPPING_EXPORT void v_VelocityLaplacian(const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, const NekDouble alpha)
Definition: Mapping.cpp:779

References v_VelocityLaplacian().

Member Data Documentation

◆ m_constantJacobian

bool Nektar::GlobalMapping::Mapping::m_constantJacobian
protected

◆ m_coords

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::Mapping::m_coords
protected

◆ m_coordsVel

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::Mapping::m_coordsVel
protected

Array with the velocity of the coordinates.

Definition at line 410 of file Mapping.h.

Referenced by Output(), v_GetCoordVelocity(), v_InitObject(), Nektar::GlobalMapping::MappingTranslation::v_InitObject(), and v_UpdateMapping().

◆ m_fields

Array<OneD, MultiRegions::ExpListSharedPtr> Nektar::GlobalMapping::Mapping::m_fields
protected

Definition at line 406 of file Mapping.h.

Referenced by Nektar::GlobalMapping::MappingGeneral::CalculateChristoffel(), Nektar::GlobalMapping::MappingXYofXY::CalculateChristoffel(), Nektar::GlobalMapping::MappingXYofXY::CalculateMetricTensor(), Nektar::GlobalMapping::MappingGeneral::CalculateMetricTerms(), ContravarFromCartesian(), ContravarToCartesian(), CovarFromCartesian(), CovarToCartesian(), EvaluateFunction(), LowerIndex(), Mapping(), Output(), RaiseIndex(), ReplaceField(), Nektar::GlobalMapping::MappingGeneral::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingTranslation::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXYofXY::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXYofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingGeneral::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingTranslation::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXYofXY::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXYofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingGeneral::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingTranslation::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXofZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingGeneral::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingTranslation::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXofZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingGeneral::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingTranslation::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXofZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingGeneral::v_CovarToCartesian(), Nektar::GlobalMapping::MappingTranslation::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXofZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_CovarToCartesian(), v_CurlCurlField(), v_Divergence(), v_DotGradJacobian(), Nektar::GlobalMapping::MappingTranslation::v_DotGradJacobian(), Nektar::GlobalMapping::MappingXofXZ::v_DotGradJacobian(), Nektar::GlobalMapping::MappingXofZ::v_DotGradJacobian(), Nektar::GlobalMapping::MappingXYofZ::v_DotGradJacobian(), v_GetCartesianCoordinates(), v_GetCoordVelocity(), Nektar::GlobalMapping::MappingGeneral::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingTranslation::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofXY::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingGeneral::v_GetJacobian(), Nektar::GlobalMapping::MappingTranslation::v_GetJacobian(), Nektar::GlobalMapping::MappingXofXZ::v_GetJacobian(), Nektar::GlobalMapping::MappingXofZ::v_GetJacobian(), Nektar::GlobalMapping::MappingXYofXY::v_GetJacobian(), Nektar::GlobalMapping::MappingXYofZ::v_GetJacobian(), Nektar::GlobalMapping::MappingGeneral::v_GetMetricTensor(), Nektar::GlobalMapping::MappingTranslation::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXYofXY::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetMetricTensor(), v_gradgradU(), v_InitObject(), Nektar::GlobalMapping::MappingTranslation::v_InitObject(), v_LowerIndex(), Nektar::GlobalMapping::MappingTranslation::v_LowerIndex(), Nektar::GlobalMapping::MappingXofXZ::v_LowerIndex(), Nektar::GlobalMapping::MappingXofZ::v_LowerIndex(), v_RaiseIndex(), Nektar::GlobalMapping::MappingTranslation::v_RaiseIndex(), Nektar::GlobalMapping::MappingXofXZ::v_RaiseIndex(), Nektar::GlobalMapping::MappingXofZ::v_RaiseIndex(), v_UpdateBCs(), Nektar::GlobalMapping::MappingXofXZ::v_UpdateGeomInfo(), Nektar::GlobalMapping::MappingXofZ::v_UpdateGeomInfo(), Nektar::GlobalMapping::MappingXYofXY::v_UpdateGeomInfo(), Nektar::GlobalMapping::MappingXYofZ::v_UpdateGeomInfo(), v_UpdateMapping(), and v_VelocityLaplacian().

◆ m_fld

LibUtilities::FieldIOSharedPtr Nektar::GlobalMapping::Mapping::m_fld
protected

Definition at line 404 of file Mapping.h.

Referenced by Mapping(), and Output().

◆ m_fromFunction

bool Nektar::GlobalMapping::Mapping::m_fromFunction
protected

Flag defining if the Mapping is defined by a function.

Definition at line 427 of file Mapping.h.

Referenced by IsFromFunction(), Output(), SetFromFunction(), v_InitObject(), and v_UpdateMapping().

◆ m_funcName

std::string Nektar::GlobalMapping::Mapping::m_funcName
protected

Name of the function containing the coordinates.

Definition at line 417 of file Mapping.h.

Referenced by Output(), v_InitObject(), and v_UpdateMapping().

◆ m_GeometricInfo

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::Mapping::m_GeometricInfo
protected

Array with metric terms of the mapping.

Definition at line 412 of file Mapping.h.

Referenced by Nektar::GlobalMapping::MappingXYofXY::CalculateMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXYofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXYofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofXZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXofZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXofZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXofZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXofZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_DotGradJacobian(), Nektar::GlobalMapping::MappingXofXZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_GetJacobian(), Nektar::GlobalMapping::MappingXYofXY::v_GetJacobian(), Nektar::GlobalMapping::MappingXofXZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_LowerIndex(), Nektar::GlobalMapping::MappingXofZ::v_LowerIndex(), Nektar::GlobalMapping::MappingXofXZ::v_RaiseIndex(), Nektar::GlobalMapping::MappingXofZ::v_RaiseIndex(), Nektar::GlobalMapping::MappingXofXZ::v_UpdateGeomInfo(), Nektar::GlobalMapping::MappingXofZ::v_UpdateGeomInfo(), Nektar::GlobalMapping::MappingXYofXY::v_UpdateGeomInfo(), and Nektar::GlobalMapping::MappingXYofZ::v_UpdateGeomInfo().

◆ m_init

bool Nektar::GlobalMapping::Mapping::m_init = false
staticprotected

Definition at line 431 of file Mapping.h.

Referenced by Load().

◆ m_isDefined

bool Nektar::GlobalMapping::Mapping::m_isDefined = false
staticprotected

Definition at line 432 of file Mapping.h.

Referenced by IsDefined(), Load(), and Output().

◆ m_mappingPtr

MappingSharedPtr Nektar::GlobalMapping::Mapping::m_mappingPtr = MappingSharedPtr()
staticprotected

Definition at line 430 of file Mapping.h.

Referenced by Load().

◆ m_nConvectiveFields

int Nektar::GlobalMapping::Mapping::m_nConvectiveFields
protected

Number of velocity components.

Definition at line 414 of file Mapping.h.

Referenced by Nektar::GlobalMapping::MappingGeneral::CalculateChristoffel(), Nektar::GlobalMapping::MappingXYofXY::CalculateChristoffel(), Nektar::GlobalMapping::MappingGeneral::CalculateMetricTerms(), ContravarFromCartesian(), ContravarToCartesian(), CovarFromCartesian(), CovarToCartesian(), LowerIndex(), Mapping(), RaiseIndex(), Nektar::GlobalMapping::MappingGeneral::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingTranslation::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXYofXY::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXYofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingGeneral::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingTranslation::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXYofXY::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXYofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingGeneral::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingTranslation::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingGeneral::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingTranslation::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingGeneral::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingTranslation::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingGeneral::v_CovarToCartesian(), Nektar::GlobalMapping::MappingTranslation::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarToCartesian(), v_CurlCurlField(), v_Divergence(), v_DotGradJacobian(), v_GetCoordVelocity(), Nektar::GlobalMapping::MappingGeneral::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingTranslation::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofXY::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingGeneral::v_GetMetricTensor(), Nektar::GlobalMapping::MappingTranslation::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXYofXY::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetMetricTensor(), v_gradgradU(), v_InitObject(), Nektar::GlobalMapping::MappingGeneral::v_InitObject(), Nektar::GlobalMapping::MappingTranslation::v_InitObject(), Nektar::GlobalMapping::MappingXofXZ::v_InitObject(), Nektar::GlobalMapping::MappingXofZ::v_InitObject(), Nektar::GlobalMapping::MappingXYofXY::v_InitObject(), Nektar::GlobalMapping::MappingXYofZ::v_InitObject(), v_LowerIndex(), Nektar::GlobalMapping::MappingTranslation::v_LowerIndex(), v_RaiseIndex(), Nektar::GlobalMapping::MappingTranslation::v_RaiseIndex(), v_UpdateBCs(), v_UpdateMapping(), and v_VelocityLaplacian().

◆ m_session

LibUtilities::SessionReaderSharedPtr Nektar::GlobalMapping::Mapping::m_session
protected

◆ m_timeDependent

bool Nektar::GlobalMapping::Mapping::m_timeDependent
protected

Flag defining if the Mapping is time-dependent.

Definition at line 425 of file Mapping.h.

Referenced by IsTimeDependent(), Output(), SetTimeDependent(), v_InitObject(), and Nektar::GlobalMapping::MappingTranslation::v_InitObject().

◆ m_tmp

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::Mapping::m_tmp
protected

◆ m_velFuncName

std::string Nektar::GlobalMapping::Mapping::m_velFuncName
protected

Name of the function containing the velocity of the coordinates.

Definition at line 419 of file Mapping.h.

Referenced by Output(), v_InitObject(), and v_UpdateMapping().

◆ m_wk1

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::Mapping::m_wk1
protected

Definition at line 435 of file Mapping.h.

Referenced by v_gradgradU(), v_InitObject(), and v_VelocityLaplacian().

◆ m_wk2

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::Mapping::m_wk2
protected

Definition at line 436 of file Mapping.h.

Referenced by v_gradgradU(), v_InitObject(), and v_VelocityLaplacian().