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 72 of file Mapping.h.

Constructor & Destructor Documentation

◆ ~Mapping()

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

Destructor.

Definition at line 76 of file Mapping.h.

77 {
78 }

◆ Mapping()

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

Constructor.

Definition at line 53 of file Mapping.cpp.

55 : m_session(pSession), m_fields(pFields)
56{
57 switch (m_fields[0]->GetExpType())
58 {
60 {
62 }
63 break;
64
66 {
68 }
69 break;
70
74 {
76 }
77 break;
78
79 default:
80 ASSERTL0(0, "Dimension not supported");
81 break;
82 }
83
85}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:418
LibUtilities::SessionReaderSharedPtr m_session
Session reader.
Definition: Mapping.h:406
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:410
LibUtilities::FieldIOSharedPtr m_fld
Definition: Mapping.h:408
static std::shared_ptr< FieldIO > CreateDefault(const LibUtilities::SessionReaderSharedPtr session)
Returns an object for the default FieldIO method.
Definition: FieldIO.cpp:194

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 215 of file Mapping.h.

218 {
219 v_ApplyChristoffelContravar(inarray, outarray);
220 }
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 233 of file Mapping.h.

236 {
237 v_ApplyChristoffelCovar(inarray, outarray);
238 }
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 539 of file Mapping.cpp.

542{
543 if (inarray == outarray)
544 {
545 int physTot = m_fields[0]->GetTotPoints();
546 int nvel = m_nConvectiveFields;
547
548 for (int i = 0; i < nvel; i++)
549 {
550 Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
551 }
553 }
554 else
555 {
556 v_ContravarFromCartesian(inarray, outarray);
557 }
558}
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:441
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

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

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

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

511{
512 if (inarray == outarray)
513 {
514 int physTot = m_fields[0]->GetTotPoints();
515 int nvel = m_nConvectiveFields;
516
517 for (int i = 0; i < nvel; i++)
518 {
519 Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
520 }
521 v_CovarToCartesian(m_tmp, outarray);
522 }
523 else
524 {
525 v_CovarToCartesian(inarray, outarray);
526 }
527}
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 328 of file Mapping.h.

331 {
332 v_CurlCurlField(inarray, outarray, generalized);
333 }
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:985

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

273 {
274 v_Divergence(inarray, outarray);
275 }
virtual GLOBAL_MAPPING_EXPORT void v_Divergence(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Mapping.cpp:748

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 173 of file Mapping.h.

176 {
177 v_DotGradJacobian(inarray, outarray);
178 }
virtual GLOBAL_MAPPING_EXPORT void v_DotGradJacobian(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Mapping.cpp:677

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

403{
404 ASSERTL0(pSession->DefinesFunction(pFunctionName),
405 "Function '" + pFunctionName + "' does not exist.");
406
407 unsigned int nq = pFields[0]->GetNpoints();
408 if (pArray.size() != nq)
409 {
410 pArray = Array<OneD, NekDouble>(nq);
411 }
412
414 vType = pSession->GetFunctionType(pFunctionName, pFieldName);
416 {
417 Array<OneD, NekDouble> x0(nq);
418 Array<OneD, NekDouble> x1(nq);
419 Array<OneD, NekDouble> x2(nq);
420
421 pFields[0]->GetCoords(x0, x1, x2);
423 pSession->GetFunction(pFunctionName, pFieldName);
424
425 ffunc->Evaluate(x0, x1, x2, pTime, pArray);
426 }
427 else if (vType == LibUtilities::eFunctionTypeFile)
428 {
429 std::string filename =
430 pSession->GetFunctionFilename(pFunctionName, pFieldName);
431
432 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
433 std::vector<std::vector<NekDouble>> FieldData;
434 Array<OneD, NekDouble> vCoeffs(pFields[0]->GetNcoeffs());
435 Vmath::Zero(vCoeffs.size(), vCoeffs, 1);
436
438 LibUtilities::FieldIO::CreateForFile(pSession, filename);
439 fld->Import(filename, FieldDef, FieldData);
440
441 int idx = -1;
442 for (int i = 0; i < FieldDef.size(); ++i)
443 {
444 for (int j = 0; j < FieldDef[i]->m_fields.size(); ++j)
445 {
446 if (FieldDef[i]->m_fields[j] == pFieldName)
447 {
448 idx = j;
449 }
450 }
451
452 if (idx >= 0)
453 {
454 pFields[0]->ExtractDataToCoeffs(FieldDef[i], FieldData[i],
455 FieldDef[i]->m_fields[idx],
456 vCoeffs);
457 }
458 else
459 {
460 std::cout << "Field " + pFieldName + " not found." << std::endl;
461 }
462 }
463 pFields[0]->BwdTrans(vCoeffs, pArray);
464 }
465}
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:223
std::shared_ptr< FieldIO > FieldIOSharedPtr
Definition: FieldIO.h:322
std::shared_ptr< Equation > EquationSharedPtr
Definition: Equation.h:125
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273

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

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

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 139 of file Mapping.h.

142 {
143 v_GetCartesianCoordinates(out0, out1, out2);
144 }
virtual GLOBAL_MAPPING_EXPORT void v_GetCartesianCoordinates(Array< OneD, NekDouble > &out0, Array< OneD, NekDouble > &out1, Array< OneD, NekDouble > &out2)
Definition: Mapping.cpp:651

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 248 of file Mapping.h.

250 {
251 v_GetCoordVelocity(outarray);
252 }
virtual GLOBAL_MAPPING_EXPORT void v_GetCoordVelocity(Array< OneD, Array< OneD, NekDouble > > &outarray)
Definition: Mapping.cpp:666

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 188 of file Mapping.h.

190 {
191 v_GetInvMetricTensor(outarray);
192 }
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 159 of file Mapping.h.

160 {
161 v_GetJacobian(outarray);
162 }
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 181 of file Mapping.h.

183 {
184 v_GetMetricTensor(outarray);
185 }
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 310 of file Mapping.h.

313 {
314 v_gradgradU(inarray, outarray);
315 }
virtual GLOBAL_MAPPING_EXPORT void v_gradgradU(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
Definition: Mapping.cpp:877

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 366 of file Mapping.h.

367 {
368 return m_constantJacobian;
369 }
bool m_constantJacobian
Flag defining if the Jacobian is constant.
Definition: Mapping.h:427

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 81 of file Mapping.h.

84 {
85 v_InitObject(pFields, pMapping);
86 }
virtual GLOBAL_MAPPING_EXPORT void v_InitObject(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
Definition: Mapping.cpp:95

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 372 of file Mapping.h.

373 {
374 return m_isDefined;
375 }

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 354 of file Mapping.h.

355 {
356 return m_fromFunction;
357 }
bool m_fromFunction
Flag defining if the Mapping is defined by a function.
Definition: Mapping.h:431

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 342 of file Mapping.h.

343 {
344 return m_timeDependent;
345 }
bool m_timeDependent
Flag defining if the Mapping is time-dependent.
Definition: Mapping.h:429

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

267{
268 if (!m_init)
269 {
270 TiXmlElement *vMapping = nullptr;
271 std::string vType;
272 if (pSession->DefinesElement("Nektar/Mapping"))
273 {
274 vMapping = pSession->GetElement("Nektar/Mapping");
275 vType = vMapping->Attribute("TYPE");
276 m_isDefined = true;
278 vMapping, pSession->GetTimeLevel());
279 }
280 else
281 {
282 vType = "Translation";
283 }
284
286 pFields, vMapping);
287
288 m_init = true;
289 }
290
291 return m_mappingPtr;
292}
static MappingSharedPtr m_mappingPtr
Definition: Mapping.h:434
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
static void GetXMLElementTimeLevel(TiXmlElement *&element, const size_t timeLevel, const bool enableCheck=true)
Get XML elment time level (Parallel-in-Time)
MappingFactory & GetMappingFactory()
Declaration of the mapping factory singleton.
Definition: Mapping.cpp:47

References Nektar::LibUtilities::NekFactory< tKey, tBase, tParam >::CreateInstance(), Nektar::GlobalMapping::GetMappingFactory(), Nektar::LibUtilities::SessionReader::GetXMLElementTimeLevel(), 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 601 of file Mapping.cpp.

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

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

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

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

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

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

241{
242 m_fields = pFields;
243
244 TiXmlElement *vMapping = nullptr;
245
246 if (m_session->DefinesElement("Nektar/Mapping"))
247 {
248 vMapping = m_session->GetElement("Nektar/Mapping");
250 vMapping, m_session->GetTimeLevel());
251 }
252 InitObject(pFields, vMapping);
253}
GLOBAL_MAPPING_EXPORT void InitObject(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
Initialise the mapping object.
Definition: Mapping.h:81

References Nektar::LibUtilities::SessionReader::GetXMLElementTimeLevel(), 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 360 of file Mapping.h.

361 {
362 m_fromFunction = value;
363 }

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 348 of file Mapping.h.

349 {
350 m_timeDependent = value;
351 }

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 382 of file Mapping.h.

383 {
384 v_UpdateBCs(time);
385 }
virtual GLOBAL_MAPPING_EXPORT void v_UpdateBCs(const NekDouble time)
Definition: Mapping.cpp:1047

References v_UpdateBCs().

◆ UpdateGeomInfo()

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

Recompute the metric terms of the Mapping.

Definition at line 399 of file Mapping.h.

400 {
402 }
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 388 of file Mapping.h.

394 {
395 v_UpdateMapping(time, coords, coordsVel);
396 }
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:1205

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

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

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

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

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::MappingTranslation, Nektar::GlobalMapping::MappingXofXZ, Nektar::GlobalMapping::MappingXofZ, and Nektar::GlobalMapping::MappingXYofZ.

Definition at line 677 of file Mapping.cpp.

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

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

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

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

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

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

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

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

Definition at line 95 of file Mapping.cpp.

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

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::MappingTranslation, Nektar::GlobalMapping::MappingXofXZ, and Nektar::GlobalMapping::MappingXofZ.

Definition at line 706 of file Mapping.cpp.

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

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::MappingTranslation, Nektar::GlobalMapping::MappingXofXZ, and Nektar::GlobalMapping::MappingXofZ.

Definition at line 727 of file Mapping.cpp.

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

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

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

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

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

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

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

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 293 of file Mapping.h.

297 {
298 v_VelocityLaplacian(inarray, outarray, alpha);
299 }
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:776

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 414 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 410 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 408 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 431 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 421 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 416 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 435 of file Mapping.h.

Referenced by Load().

◆ m_isDefined

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

Definition at line 436 of file Mapping.h.

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

◆ m_mappingPtr

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

Definition at line 434 of file Mapping.h.

Referenced by Load().

◆ m_nConvectiveFields

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

Number of velocity components.

Definition at line 418 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 429 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 423 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 439 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 440 of file Mapping.h.

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