Nektar++
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
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 59 of file Mapping.cpp.

61  : m_session(pSession), m_fields(pFields)
62 {
63  switch (m_fields[0]->GetExpType())
64  {
65  case MultiRegions::e1D:
66  {
68  }
69  break;
70 
71  case MultiRegions::e2D:
72  {
74  }
75  break;
76 
77  case MultiRegions::e3D:
80  {
82  }
83  break;
84 
85  default:
86  ASSERTL0(0, "Dimension not supported");
87  break;
88  }
89 
91 }
#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:197

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 543 of file Mapping.cpp.

546 {
547  if (inarray == outarray)
548  {
549  int physTot = m_fields[0]->GetTotPoints();
550  int nvel = m_nConvectiveFields;
551 
552  for (int i = 0; i < nvel; i++)
553  {
554  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
555  }
556  v_ContravarFromCartesian(m_tmp, outarray);
557  }
558  else
559  {
560  v_ContravarFromCartesian(inarray, outarray);
561  }
562 }
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 481 of file Mapping.cpp.

484 {
485  if (inarray == outarray)
486  {
487  int physTot = m_fields[0]->GetTotPoints();
488  int nvel = m_nConvectiveFields;
489 
490  for (int i = 0; i < nvel; i++)
491  {
492  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
493  }
494  v_ContravarToCartesian(m_tmp, outarray);
495  }
496  else
497  {
498  v_ContravarToCartesian(inarray, outarray);
499  }
500 }
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 574 of file Mapping.cpp.

577 {
578  if (inarray == outarray)
579  {
580  int physTot = m_fields[0]->GetTotPoints();
581  int nvel = m_nConvectiveFields;
582 
583  for (int i = 0; i < nvel; i++)
584  {
585  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
586  }
587  v_CovarFromCartesian(m_tmp, outarray);
588  }
589  else
590  {
591  v_CovarFromCartesian(inarray, outarray);
592  }
593 }
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 512 of file Mapping.cpp.

515 {
516  if (inarray == outarray)
517  {
518  int physTot = m_fields[0]->GetTotPoints();
519  int nvel = m_nConvectiveFields;
520 
521  for (int i = 0; i < nvel; i++)
522  {
523  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
524  }
525  v_CovarToCartesian(m_tmp, outarray);
526  }
527  else
528  {
529  v_CovarToCartesian(inarray, outarray);
530  }
531 }
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:989

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:752

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:681

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 402 of file Mapping.cpp.

407 {
408  ASSERTL0(pSession->DefinesFunction(pFunctionName),
409  "Function '" + pFunctionName + "' does not exist.");
410 
411  unsigned int nq = pFields[0]->GetNpoints();
412  if (pArray.size() != nq)
413  {
414  pArray = Array<OneD, NekDouble>(nq);
415  }
416 
418  vType = pSession->GetFunctionType(pFunctionName, pFieldName);
420  {
421  Array<OneD, NekDouble> x0(nq);
422  Array<OneD, NekDouble> x1(nq);
423  Array<OneD, NekDouble> x2(nq);
424 
425  pFields[0]->GetCoords(x0, x1, x2);
427  pSession->GetFunction(pFunctionName, pFieldName);
428 
429  ffunc->Evaluate(x0, x1, x2, pTime, pArray);
430  }
431  else if (vType == LibUtilities::eFunctionTypeFile)
432  {
433  std::string filename =
434  pSession->GetFunctionFilename(pFunctionName, pFieldName);
435 
436  std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
437  std::vector<std::vector<NekDouble>> FieldData;
438  Array<OneD, NekDouble> vCoeffs(pFields[0]->GetNcoeffs());
439  Vmath::Zero(vCoeffs.size(), vCoeffs, 1);
440 
442  LibUtilities::FieldIO::CreateForFile(pSession, filename);
443  fld->Import(filename, FieldDef, FieldData);
444 
445  int idx = -1;
446  for (int i = 0; i < FieldDef.size(); ++i)
447  {
448  for (int j = 0; j < FieldDef[i]->m_fields.size(); ++j)
449  {
450  if (FieldDef[i]->m_fields[j] == pFieldName)
451  {
452  idx = j;
453  }
454  }
455 
456  if (idx >= 0)
457  {
458  pFields[0]->ExtractDataToCoeffs(FieldDef[i], FieldData[i],
459  FieldDef[i]->m_fields[idx],
460  vCoeffs);
461  }
462  else
463  {
464  cout << "Field " + pFieldName + " not found." << endl;
465  }
466  }
467  pFields[0]->BwdTrans(vCoeffs, pArray);
468  }
469 }
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:226
std::shared_ptr< FieldIO > FieldIOSharedPtr
Definition: FieldIO.h:327
std::shared_ptr< Equation > EquationSharedPtr
Definition: Equation.h:129
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 384 of file Mapping.cpp.

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

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:655

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:670

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:881

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:101

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 270 of file Mapping.cpp.

273 {
274  if (!m_init)
275  {
276  TiXmlElement *vMapping = NULL;
277  string vType;
278  if (pSession->DefinesElement("Nektar/Mapping"))
279  {
280  vMapping = pSession->GetElement("Nektar/Mapping");
281  vType = vMapping->Attribute("TYPE");
282  m_isDefined = true;
283  }
284  else
285  {
286  vType = "Translation";
287  }
288 
289  m_mappingPtr = GetMappingFactory().CreateInstance(vType, pSession,
290  pFields, vMapping);
291 
292  m_init = true;
293  }
294 
295  return m_mappingPtr;
296 }
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:53

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 605 of file Mapping.cpp.

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

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 306 of file Mapping.cpp.

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

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

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 247 of file Mapping.cpp.

249 {
250  m_fields = pFields;
251 
252  TiXmlElement *vMapping = NULL;
253 
254  if (m_session->DefinesElement("Nektar/Mapping"))
255  {
256  vMapping = m_session->GetElement("Nektar/Mapping");
257  }
258  InitObject(pFields, vMapping);
259 }
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:1051

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:1209

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 989 of file Mapping.cpp.

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

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

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

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

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 670 of file Mapping.cpp.

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

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 881 of file Mapping.cpp.

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

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

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 710 of file Mapping.cpp.

712 {
713  int physTot = m_fields[0]->GetTotPoints();
714  int nvel = m_nConvectiveFields;
715 
716  Array<OneD, Array<OneD, NekDouble>> g(nvel * nvel);
717 
718  GetMetricTensor(g);
719 
720  for (int i = 0; i < nvel; i++)
721  {
722  outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
723  for (int j = 0; j < nvel; j++)
724  {
725  Vmath::Vvtvp(physTot, g[i * nvel + j], 1, inarray[j], 1,
726  outarray[i], 1, outarray[i], 1);
727  }
728  }
729 }
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 731 of file Mapping.cpp.

733 {
734  int physTot = m_fields[0]->GetTotPoints();
735  int nvel = m_nConvectiveFields;
736 
737  Array<OneD, Array<OneD, NekDouble>> g(nvel * nvel);
738 
740 
741  for (int i = 0; i < nvel; i++)
742  {
743  outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
744  for (int j = 0; j < nvel; j++)
745  {
746  Vmath::Vvtvp(physTot, g[i * nvel + j], 1, inarray[j], 1,
747  outarray[i], 1, outarray[i], 1);
748  }
749  }
750 }
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 1051 of file Mapping.cpp.

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

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

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 780 of file Mapping.cpp.

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

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