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:
Inheritance graph
[legend]
Collaboration diagram for Nektar::GlobalMapping::Mapping:
Collaboration graph
[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::ExpListSharedPtr
m_fields
 
Array< OneD, Array< OneD,
NekDouble > > 
m_coords
 Array with the Cartesian coordinates. More...
 
Array< OneD, Array< OneD,
NekDouble > > 
m_coordsVel
 Array with the velocity of the coordinates. More...
 
Array< OneD, Array< OneD,
NekDouble > > 
m_GeometricInfo
 Array with metric terms of the mapping. More...
 
int m_nConvectiveFields
 Number of velocity components. More...
 
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 69 of file Mapping.h.

Constructor & Destructor Documentation

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

Destructor.

Definition at line 73 of file Mapping.h.

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

Constructor.

Definition at line 57 of file Mapping.cpp.

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.

59  : m_session(pSession), m_fields(pFields)
60 {
61  switch (m_fields[0]->GetExpType())
62  {
63  case MultiRegions::e1D:
64  {
66  }
67  break;
68 
69  case MultiRegions::e2D:
70  {
72  }
73  break;
74 
75  case MultiRegions::e3D:
78  {
80  }
81  break;
82 
83  default:
84  ASSERTL0(0,"Dimension not supported");
85  break;
86  }
87 
89 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
LibUtilities::SessionReaderSharedPtr m_session
Session reader.
Definition: Mapping.h:405
LibUtilities::FieldIOSharedPtr m_fld
Definition: Mapping.h:407
static boost::shared_ptr< FieldIO > CreateDefault(const LibUtilities::SessionReaderSharedPtr session)
Returns an object for the default FieldIO method.
Definition: FieldIO.cpp:181

Member Function Documentation

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

References v_ApplyChristoffelContravar().

Referenced by v_gradgradU(), and v_VelocityLaplacian().

215  {
216  v_ApplyChristoffelContravar( inarray, outarray);
217  }
virtual GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelContravar(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
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 230 of file Mapping.h.

References v_ApplyChristoffelCovar().

Referenced by v_gradgradU().

233  {
234  v_ApplyChristoffelCovar( inarray, outarray);
235  }
virtual GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelCovar(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
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 546 of file Mapping.cpp.

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

Referenced by v_UpdateBCs().

549 {
550  if(inarray == outarray)
551  {
552  int physTot = m_fields[0]->GetTotPoints();
553  int nvel = m_nConvectiveFields;
554 
555  for (int i=0; i< nvel; i++)
556  {
557  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
558  }
559  v_ContravarFromCartesian( m_tmp, outarray);
560  }
561  else
562  {
563  v_ContravarFromCartesian( inarray, outarray);
564  }
565 }
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
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.cpp:1061
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 486 of file Mapping.cpp.

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

489 {
490  if(inarray == outarray)
491  {
492  int physTot = m_fields[0]->GetTotPoints();
493  int nvel = m_nConvectiveFields;
494 
495  for (int i=0; i< nvel; i++)
496  {
497  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
498  }
499  v_ContravarToCartesian( m_tmp, outarray);
500  }
501  else
502  {
503  v_ContravarToCartesian( inarray, outarray);
504  }
505 }
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
Array< OneD, Array< OneD, NekDouble > > m_tmp
Definition: Mapping.h:441
virtual GLOBAL_MAPPING_EXPORT void v_ContravarToCartesian(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
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 576 of file Mapping.cpp.

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

579 {
580  if(inarray == outarray)
581  {
582  int physTot = m_fields[0]->GetTotPoints();
583  int nvel = m_nConvectiveFields;
584 
585  for (int i=0; i< nvel; i++)
586  {
587  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
588  }
589  v_CovarFromCartesian( m_tmp, outarray);
590  }
591  else
592  {
593  v_CovarFromCartesian( inarray, outarray);
594  }
595 }
virtual GLOBAL_MAPPING_EXPORT void v_CovarFromCartesian(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
Array< OneD, Array< OneD, NekDouble > > m_tmp
Definition: Mapping.h:441
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
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 516 of file Mapping.cpp.

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

519 {
520  if(inarray == outarray)
521  {
522  int physTot = m_fields[0]->GetTotPoints();
523  int nvel = m_nConvectiveFields;
524 
525  for (int i=0; i< nvel; i++)
526  {
527  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
528  }
529  v_CovarToCartesian( m_tmp, outarray);
530  }
531  else
532  {
533  v_CovarToCartesian( inarray, outarray);
534  }
535 }
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
virtual GLOBAL_MAPPING_EXPORT void v_CovarToCartesian(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
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.cpp:1061
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 325 of file Mapping.h.

References v_CurlCurlField().

329  {
330  v_CurlCurlField( inarray, outarray, generalized);
331  }
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:1009
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 267 of file Mapping.h.

References v_Divergence().

270  {
271  v_Divergence( inarray, outarray);
272  }
virtual GLOBAL_MAPPING_EXPORT void v_Divergence(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Mapping.cpp:763
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 170 of file Mapping.h.

References v_DotGradJacobian().

173  {
174  v_DotGradJacobian( inarray, outarray);
175  }
virtual GLOBAL_MAPPING_EXPORT void v_DotGradJacobian(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Mapping.cpp:687
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 403 of file Mapping.cpp.

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

Referenced by v_InitObject(), and v_UpdateMapping().

410 {
411  ASSERTL0(pSession->DefinesFunction(pFunctionName),
412  "Function '" + pFunctionName + "' does not exist.");
413 
414  unsigned int nq = pFields[0]->GetNpoints();
415  if (pArray.num_elements() != nq)
416  {
417  pArray = Array<OneD, NekDouble> (nq);
418  }
419 
421  vType = pSession->GetFunctionType(pFunctionName, pFieldName);
423  {
424  Array<OneD, NekDouble> x0(nq);
425  Array<OneD, NekDouble> x1(nq);
426  Array<OneD, NekDouble> x2(nq);
427 
428  pFields[0]->GetCoords(x0, x1, x2);
430  pSession->GetFunction(pFunctionName, pFieldName);
431 
432  ffunc->Evaluate(x0, x1, x2, pTime, pArray);
433  }
434  else if (vType == LibUtilities::eFunctionTypeFile)
435  {
436  std::string filename = pSession->GetFunctionFilename(
437  pFunctionName,
438  pFieldName);
439 
440  std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
441  std::vector<std::vector<NekDouble> > FieldData;
442  Array<OneD, NekDouble> vCoeffs(pFields[0]->GetNcoeffs());
443  Vmath::Zero(vCoeffs.num_elements(), vCoeffs, 1);
444 
446  LibUtilities::FieldIO::CreateForFile(pSession, filename);
447  fld->Import(filename, FieldDef, FieldData);
448 
449  int idx = -1;
450  for (int i = 0; i < FieldDef.size(); ++i)
451  {
452  for (int j = 0; j < FieldDef[i]->m_fields.size(); ++j)
453  {
454  if (FieldDef[i]->m_fields[j] == pFieldName)
455  {
456  idx = j;
457  }
458  }
459 
460  if (idx >= 0)
461  {
462  pFields[0]->ExtractDataToCoeffs(
463  FieldDef[i],
464  FieldData[i],
465  FieldDef[i]->m_fields[idx],
466  vCoeffs);
467  }
468  else
469  {
470  cout << "Field " + pFieldName + " not found." << endl;
471  }
472  }
473  pFields[0]->BwdTrans_IterPerExp(vCoeffs, pArray);
474  }
475 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
boost::shared_ptr< FieldIO > FieldIOSharedPtr
Definition: FieldIO.h:309
boost::shared_ptr< Equation > EquationSharedPtr
static boost::shared_ptr< FieldIO > CreateForFile(const LibUtilities::SessionReaderSharedPtr session, const std::string &filename)
Construct a FieldIO object for a given input filename.
Definition: FieldIO.cpp:212
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:373
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 382 of file Mapping.cpp.

References ASSERTL0.

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 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
boost::shared_ptr< Equation > EquationSharedPtr
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 134 of file Mapping.h.

References v_GetCartesianCoordinates().

Referenced by v_UpdateBCs().

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:659
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 245 of file Mapping.h.

References v_GetCoordVelocity().

Referenced by v_UpdateBCs().

247  {
248  v_GetCoordVelocity( outarray);
249  }
virtual GLOBAL_MAPPING_EXPORT void v_GetCoordVelocity(Array< OneD, Array< OneD, NekDouble > > &outarray)
Definition: Mapping.cpp:675
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 185 of file Mapping.h.

References v_GetInvMetricTensor().

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

187  {
188  v_GetInvMetricTensor( outarray);
189  }
virtual GLOBAL_MAPPING_EXPORT void v_GetInvMetricTensor(Array< OneD, Array< OneD, NekDouble > > &outarray)=0
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.

References v_GetJacobian().

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

157  {
158  v_GetJacobian( outarray);
159  }
virtual GLOBAL_MAPPING_EXPORT void v_GetJacobian(Array< OneD, NekDouble > &outarray)=0
GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::GetMetricTensor ( Array< OneD, Array< OneD, NekDouble > > &  outarray)
inline

Get the metric tensor $g_{ij}$.

Definition at line 178 of file Mapping.h.

References v_GetMetricTensor().

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

180  {
181  v_GetMetricTensor( outarray);
182  }
virtual GLOBAL_MAPPING_EXPORT void v_GetMetricTensor(Array< OneD, Array< OneD, NekDouble > > &outarray)=0
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 307 of file Mapping.h.

References v_gradgradU().

Referenced by v_CurlCurlField().

310  {
311  v_gradgradU( inarray, outarray);
312  }
virtual GLOBAL_MAPPING_EXPORT void v_gradgradU(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
Definition: Mapping.cpp:899
GLOBAL_MAPPING_EXPORT bool Nektar::GlobalMapping::Mapping::HasConstantJacobian ( )
inline

Get flag defining if mapping has constant Jacobian.

Definition at line 365 of file Mapping.h.

References m_constantJacobian.

Referenced by v_DotGradJacobian().

366  {
367  return m_constantJacobian;
368  }
bool m_constantJacobian
Flag defining if the Jacobian is constant.
Definition: Mapping.h:427
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 76 of file Mapping.h.

References v_InitObject().

Referenced by ReplaceField().

79  {
80  v_InitObject( pFields, pMapping);
81  }
virtual GLOBAL_MAPPING_EXPORT void v_InitObject(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
Definition: Mapping.cpp:99
GLOBAL_MAPPING_EXPORT bool Nektar::GlobalMapping::Mapping::IsDefined ( )
inline

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

Definition at line 371 of file Mapping.h.

References m_isDefined.

372  {
373  return m_isDefined;
374  }
GLOBAL_MAPPING_EXPORT bool Nektar::GlobalMapping::Mapping::IsFromFunction ( )
inline

Get flag defining if mapping is defined by a function.

Definition at line 353 of file Mapping.h.

References m_fromFunction.

354  {
355  return m_fromFunction;
356  }
bool m_fromFunction
Flag defining if the Mapping is defined by a function.
Definition: Mapping.h:431
GLOBAL_MAPPING_EXPORT bool Nektar::GlobalMapping::Mapping::IsTimeDependent ( void  )
inline

Get flag defining if mapping is time-dependent.

Definition at line 341 of file Mapping.h.

References m_timeDependent.

Referenced by v_UpdateBCs().

342  {
343  return m_timeDependent;
344  }
bool m_timeDependent
Flag defining if the Mapping is time-dependent.
Definition: Mapping.h:429
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 265 of file Mapping.cpp.

References Nektar::LibUtilities::NekFactory< tKey, tBase, >::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().

268 {
269  if (!m_init)
270  {
271  TiXmlElement* vMapping = NULL;
272  string vType;
273  if (pSession->DefinesElement("Nektar/Mapping"))
274  {
275  vMapping = pSession->GetElement("Nektar/Mapping");
276  vType = vMapping->Attribute("TYPE");
277  m_isDefined = true;
278  }
279  else
280  {
281  vType = "Translation";
282  }
283 
285  vType, pSession, pFields,
286  vMapping);
287 
288  m_init = true;
289  }
290 
291  return m_mappingPtr;
292 }
static MappingSharedPtr m_mappingPtr
Definition: Mapping.h:434
MappingFactory & GetMappingFactory()
Declaration of the mapping factory singleton.
Definition: Mapping.cpp:49
tBaseSharedPtr CreateInstance(tKey idKey BOOST_PP_COMMA_IF(MAX_PARAM) BOOST_PP_ENUM_BINARY_PARAMS(MAX_PARAM, tParam, x))
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:162
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 607 of file Mapping.cpp.

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

610 {
611  if(inarray == outarray)
612  {
613  int physTot = m_fields[0]->GetTotPoints();
614  int nvel = m_nConvectiveFields;
615 
616  for (int i=0; i< nvel; i++)
617  {
618  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
619  }
620  v_LowerIndex( m_tmp, outarray);
621  }
622  else
623  {
624  v_LowerIndex( inarray, outarray);
625  }
626 }
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
Array< OneD, Array< OneD, NekDouble > > m_tmp
Definition: Mapping.h:441
virtual GLOBAL_MAPPING_EXPORT void v_LowerIndex(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
Definition: Mapping.cpp:717
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
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.

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

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

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

Referenced by v_CurlCurlField(), and v_VelocityLaplacian().

641 {
642  if(inarray == outarray)
643  {
644  int physTot = m_fields[0]->GetTotPoints();
645  int nvel = m_nConvectiveFields;
646 
647  for (int i=0; i< nvel; i++)
648  {
649  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
650  }
651  v_RaiseIndex( m_tmp, outarray);
652  }
653  else
654  {
655  v_RaiseIndex( inarray, outarray);
656  }
657 }
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
Array< OneD, Array< OneD, NekDouble > > m_tmp
Definition: Mapping.h:441
virtual GLOBAL_MAPPING_EXPORT void v_RaiseIndex(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
Definition: Mapping.cpp:740
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
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 242 of file Mapping.cpp.

References InitObject(), m_fields, and m_session.

244 {
245  m_fields = pFields;
246 
247  TiXmlElement* vMapping = NULL;
248 
249  if (m_session->DefinesElement("Nektar/Mapping"))
250  {
251  vMapping = m_session->GetElement("Nektar/Mapping");
252  }
253  InitObject(pFields, vMapping);
254 }
GLOBAL_MAPPING_EXPORT void InitObject(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
Initialise the mapping object.
Definition: Mapping.h:76
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
LibUtilities::SessionReaderSharedPtr m_session
Session reader.
Definition: Mapping.h:405
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 359 of file Mapping.h.

References m_fromFunction.

360  {
361  m_fromFunction = value;
362  }
bool m_fromFunction
Flag defining if the Mapping is defined by a function.
Definition: Mapping.h:431
GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::SetTimeDependent ( const bool  value)
inline

Set flag defining if mapping is time-dependent.

Definition at line 347 of file Mapping.h.

References m_timeDependent.

348  {
349  m_timeDependent = value;
350  }
bool m_timeDependent
Flag defining if the Mapping is time-dependent.
Definition: Mapping.h:429
GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::UpdateBCs ( const NekDouble  time)
inline

Update the Dirichlet Boundary Conditions when using Mappings.

Definition at line 381 of file Mapping.h.

References v_UpdateBCs().

382  {
383  v_UpdateBCs(time);
384  }
virtual GLOBAL_MAPPING_EXPORT void v_UpdateBCs(const NekDouble time)
Definition: Mapping.cpp:1072
GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::UpdateGeomInfo ( )
inline

Recompute the metric terms of the Mapping.

Definition at line 397 of file Mapping.h.

References v_UpdateGeomInfo().

Referenced by v_InitObject(), and v_UpdateMapping().

398  {
400  }
virtual GLOBAL_MAPPING_EXPORT void v_UpdateGeomInfo()=0
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 387 of file Mapping.h.

References v_UpdateMapping().

392  {
393  v_UpdateMapping( time, coords, coordsVel);
394  }
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:1273
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
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
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
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
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
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
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 1009 of file Mapping.cpp.

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

Referenced by CurlCurlField().

1013 {
1014  int physTot = m_fields[0]->GetTotPoints();
1015  int nvel = m_nConvectiveFields;
1016 
1017  // Set wavespace to false and store current value
1018  bool wavespace = m_fields[0]->GetWaveSpace();
1019  m_fields[0]->SetWaveSpace(false);
1020 
1021  // For implicit treatment of viscous terms, we want the generalized
1022  // curlcurl and for explicit treatment, we want the cartesian one.
1023  if (generalized)
1024  {
1025  // Get the second derivatives u^i_{,jk}
1026  Array<OneD, Array<OneD, NekDouble> > tmp(nvel);
1027  Array<OneD, Array<OneD, NekDouble> > ddU(nvel*nvel*nvel);
1028  gradgradU(inarray, ddU);
1029 
1030  // Raise index to obtain A^{ip}_{k} = g^pj u^i_{,jk}
1031  for (int i = 0; i < nvel; ++i)
1032  {
1033  for (int k = 0; k < nvel; ++k)
1034  {
1035  // Copy to wk
1036  for (int j = 0; j < nvel; ++j)
1037  {
1038  tmp[j] = ddU[i*nvel*nvel+j*nvel+k];
1039  }
1040  RaiseIndex(tmp, m_tmp);
1041  for (int p=0; p<nvel; ++p)
1042  {
1043  Vmath::Vcopy(physTot, m_tmp[p], 1,
1044  ddU[i*nvel*nvel+p*nvel+k], 1);
1045  }
1046  }
1047  }
1048  // The curlcurl is g^ji u^k_{kj} - g^jk u^i_kj = A^{ki}_k - A^{ik}_k
1049  for (int i = 0; i < nvel; ++i)
1050  {
1051  outarray[i] = Array<OneD, NekDouble> (physTot, 0.0);
1052  for (int k = 0; k < nvel; ++k)
1053  {
1054  Vmath::Vadd(physTot, outarray[i], 1,
1055  ddU[k*nvel*nvel+i*nvel+k], 1,
1056  outarray[i], 1);
1057  Vmath::Vsub(physTot, outarray[i], 1,
1058  ddU[i*nvel*nvel+k*nvel+k], 1,
1059  outarray[i], 1);
1060  }
1061  }
1062  }
1063  else
1064  {
1065  m_fields[0]->CurlCurl(inarray, outarray);
1066  }
1067 
1068  // Restore value of wavespace
1069  m_fields[0]->SetWaveSpace(wavespace);
1070 }
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
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:638
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
Array< OneD, Array< OneD, NekDouble > > m_tmp
Definition: Mapping.h:441
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:307
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:343
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
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:299
void Nektar::GlobalMapping::Mapping::v_Divergence ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Definition at line 763 of file Mapping.cpp.

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

Referenced by Divergence().

766 {
767  int physTot = m_fields[0]->GetTotPoints();
768  Array<OneD, NekDouble> wk(physTot, 0.0);
769 
770  Vmath::Zero(physTot, outarray, 1);
771 
772  // Set wavespace to false and store current value
773  bool wavespace = m_fields[0]->GetWaveSpace();
774  m_fields[0]->SetWaveSpace(false);
775 
776  // Get Mapping Jacobian
777  Array<OneD, NekDouble> Jac(physTot, 0.0);
778  GetJacobian(Jac);
779 
780  for(int i = 0; i < m_nConvectiveFields; ++i)
781  {
782  Vmath::Vmul(physTot,Jac, 1, inarray[i], 1, wk, 1); // J*Ui
783  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[i],
784  wk, wk); // (J*Ui)_i
785  Vmath::Vadd(physTot, wk, 1, outarray, 1, outarray, 1);
786  }
787  Vmath::Vdiv(physTot,outarray,1,Jac,1,outarray,1); //1/J*(J*Ui)_i
788 
789  m_fields[0]->SetWaveSpace(wavespace);
790 }
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:241
GLOBAL_MAPPING_EXPORT void GetJacobian(Array< OneD, NekDouble > &outarray)
Get the Jacobian of the transformation.
Definition: Mapping.h:155
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:373
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:299
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:183
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 687 of file Mapping.cpp.

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

Referenced by DotGradJacobian().

690 {
691  int physTot = m_fields[0]->GetTotPoints();
692 
693  outarray = Array<OneD, NekDouble>(physTot, 0.0);
694  if ( !HasConstantJacobian() )
695  {
696  // Set wavespace to false and store current value
697  bool wavespace = m_fields[0]->GetWaveSpace();
698  m_fields[0]->SetWaveSpace(false);
699 
700  // Get Mapping Jacobian
701  Array<OneD, NekDouble> Jac(physTot, 0.0);
702  GetJacobian(Jac);
703 
704  // Calculate inarray . grad(Jac)
705  Array<OneD, NekDouble> wk(physTot, 0.0);
706  for(int i = 0; i < m_nConvectiveFields; ++i)
707  {
708  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[i],
709  Jac, wk);
710  Vmath::Vvtvp(physTot, inarray[i], 1, wk, 1,
711  outarray, 1, outarray, 1);
712  }
713  m_fields[0]->SetWaveSpace(wavespace);
714  }
715 }
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:442
GLOBAL_MAPPING_EXPORT void GetJacobian(Array< OneD, NekDouble > &outarray)
Get the Jacobian of the transformation.
Definition: Mapping.h:155
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
GLOBAL_MAPPING_EXPORT bool HasConstantJacobian()
Get flag defining if mapping has constant Jacobian.
Definition: Mapping.h:365
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
void Nektar::GlobalMapping::Mapping::v_GetCartesianCoordinates ( Array< OneD, NekDouble > &  out0,
Array< OneD, NekDouble > &  out1,
Array< OneD, NekDouble > &  out2 
)
protectedvirtual

Definition at line 659 of file Mapping.cpp.

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

Referenced by GetCartesianCoordinates().

663 {
664  int physTot = m_fields[0]->GetTotPoints();
665 
666  out0 = Array<OneD, NekDouble>(physTot, 0.0);
667  out1 = Array<OneD, NekDouble>(physTot, 0.0);
668  out2 = Array<OneD, NekDouble>(physTot, 0.0);
669 
670  Vmath::Vcopy(physTot, m_coords[0], 1, out0, 1);
671  Vmath::Vcopy(physTot, m_coords[1], 1, out1, 1);
672  Vmath::Vcopy(physTot, m_coords[2], 1, out2, 1);
673 }
Array< OneD, Array< OneD, NekDouble > > m_coords
Array with the Cartesian coordinates.
Definition: Mapping.h:411
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
void Nektar::GlobalMapping::Mapping::v_GetCoordVelocity ( Array< OneD, Array< OneD, NekDouble > > &  outarray)
protectedvirtual

Definition at line 675 of file Mapping.cpp.

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

Referenced by GetCoordVelocity().

677 {
678  int physTot = m_fields[0]->GetTotPoints();
679 
680  for(int i = 0; i < m_nConvectiveFields; ++i)
681  {
682  outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
683  Vmath::Vcopy(physTot, m_coordsVel[i], 1, outarray[i], 1);
684  }
685 }
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
Array< OneD, Array< OneD, NekDouble > > m_coordsVel
Array with the velocity of the coordinates.
Definition: Mapping.h:413
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_GetInvMetricTensor ( Array< OneD, Array< OneD, NekDouble > > &  outarray)
protectedpure virtual
virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_GetJacobian ( Array< OneD, NekDouble > &  outarray)
protectedpure virtual
virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_GetMetricTensor ( Array< OneD, Array< OneD, NekDouble > > &  outarray)
protectedpure virtual
void Nektar::GlobalMapping::Mapping::v_gradgradU ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protectedvirtual

Definition at line 899 of file Mapping.cpp.

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

Referenced by gradgradU().

902 {
903  int physTot = m_fields[0]->GetTotPoints();
904  int nvel = m_nConvectiveFields;
905 
906  // Declare variables
907  outarray = Array<OneD, Array<OneD, NekDouble> > (nvel*nvel*nvel);
908  Array<OneD, Array<OneD, NekDouble> > tmp(nvel);
909  for (int i=0; i< nvel*nvel*nvel; i++)
910  {
911  outarray[i] = Array<OneD, NekDouble>(physTot,0.0);
912  }
913 
914  // Set wavespace to false and store current value
915  bool wavespace = m_fields[0]->GetWaveSpace();
916  m_fields[0]->SetWaveSpace(false);
917 
918  // Calculate vector gradient u^i_(,j) = du^i/dx^j + {i,pj}*u^p
919  ApplyChristoffelContravar(inarray, m_wk1);
920  for (int i=0; i< nvel; i++)
921  {
922  if (nvel == 2)
923  {
924  m_fields[0]->PhysDeriv(inarray[i],
925  m_wk2[i*nvel+0],
926  m_wk2[i*nvel+1]);
927  }
928  else
929  {
930  m_fields[0]->PhysDeriv(inarray[i],
931  m_wk2[i*nvel+0],
932  m_wk2[i*nvel+1],
933  m_wk2[i*nvel+2]);
934  }
935  for (int j=0; j< nvel; j++)
936  {
937  Vmath::Vadd(physTot,m_wk1[i*nvel+j],1,m_wk2[i*nvel+j],1,
938  m_wk1[i*nvel+j], 1);
939  }
940  }
941 
942  //
943  // Calculate (u^i_,j),k
944  //
945 
946  // Step 1 : d(u^i_,j))/d(x^k)
947  for (int i=0; i< nvel; i++)
948  {
949  for (int j=0; j< nvel; j++)
950  {
951  if (nvel == 2)
952  {
953  m_fields[0]->PhysDeriv(m_wk1[i*nvel+j],
954  outarray[i*nvel*nvel+j*nvel+0],
955  outarray[i*nvel*nvel+j*nvel+1]);
956  }
957  else
958  {
959  m_fields[0]->PhysDeriv(m_wk1[i*nvel+j],
960  outarray[i*nvel*nvel+j*nvel+0],
961  outarray[i*nvel*nvel+j*nvel+1],
962  outarray[i*nvel*nvel+j*nvel+2]);
963  }
964  }
965  }
966 
967  // Step 2: d(u^i_,j)/d(x^k) - {p,jk}*u^i_,p
968  for (int i=0; i< nvel; i++)
969  {
970  for (int p=0; p< nvel; p++)
971  {
972  tmp[p] = m_wk1[i*nvel+p];
973  }
975  for (int j=0; j< nvel; j++)
976  {
977  for (int k=0; k< nvel; k++)
978  {
979  Vmath::Vsub(physTot,outarray[i*nvel*nvel+j*nvel+k],1,
980  m_wk2[j*nvel+k],1,
981  outarray[i*nvel*nvel+j*nvel+k], 1);
982  }
983  }
984  }
985 
986  // Step 3: d(u^i_,j)/d(x^k) - {p,jk}*u^i_,p + {i,pk} u^p_,j
987  for (int j=0; j< nvel; j++)
988  {
989  for (int p=0; p< nvel; p++)
990  {
991  tmp[p] = m_wk1[p*nvel+j];
992  }
994  for (int i=0; i< nvel; i++)
995  {
996  for (int k=0; k< nvel; k++)
997  {
998  Vmath::Vadd(physTot,outarray[i*nvel*nvel+j*nvel+k],1,
999  m_wk2[i*nvel+k],1,
1000  outarray[i*nvel*nvel+j*nvel+k], 1);
1001  }
1002  }
1003  }
1004 
1005  // Restore value of wavespace
1006  m_fields[0]->SetWaveSpace(wavespace);
1007 }
Array< OneD, Array< OneD, NekDouble > > m_wk2
Definition: Mapping.h:440
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
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:212
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:343
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:230
Array< OneD, Array< OneD, NekDouble > > m_wk1
Definition: Mapping.h:439
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:299
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::MappingXYofXY, Nektar::GlobalMapping::MappingTranslation, Nektar::GlobalMapping::MappingXofXZ, Nektar::GlobalMapping::MappingXofZ, and Nektar::GlobalMapping::MappingXYofZ.

Definition at line 99 of file Mapping.cpp.

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

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

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

Referenced by LowerIndex().

720 {
721  int physTot = m_fields[0]->GetTotPoints();
722  int nvel = m_nConvectiveFields;
723 
724  Array<OneD, Array<OneD, NekDouble> > g(nvel*nvel);
725 
726  GetMetricTensor(g);
727 
728  for (int i=0; i< nvel; i++)
729  {
730  outarray[i] = Array<OneD, NekDouble> (physTot, 0.0);
731  for (int j=0; j< nvel; j++)
732  {
733  Vmath::Vvtvp(physTot, g[i*nvel+j], 1, inarray[j], 1,
734  outarray[i], 1,
735  outarray[i], 1);
736  }
737  }
738 }
GLOBAL_MAPPING_EXPORT void GetMetricTensor(Array< OneD, Array< OneD, NekDouble > > &outarray)
Get the metric tensor .
Definition: Mapping.h:178
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:442
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
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 740 of file Mapping.cpp.

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

Referenced by RaiseIndex().

743 {
744  int physTot = m_fields[0]->GetTotPoints();
745  int nvel = m_nConvectiveFields;
746 
747  Array<OneD, Array<OneD, NekDouble> > g(nvel*nvel);
748 
750 
751  for (int i=0; i< nvel; i++)
752  {
753  outarray[i] = Array<OneD, NekDouble> (physTot, 0.0);
754  for (int j=0; j< nvel; j++)
755  {
756  Vmath::Vvtvp(physTot, g[i*nvel+j], 1, inarray[j], 1,
757  outarray[i], 1,
758  outarray[i], 1);
759  }
760  }
761 }
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:442
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
GLOBAL_MAPPING_EXPORT void GetInvMetricTensor(Array< OneD, Array< OneD, NekDouble > > &outarray)
Get the inverse of metric tensor .
Definition: Mapping.h:185
void Nektar::GlobalMapping::Mapping::v_UpdateBCs ( const NekDouble  time)
protectedvirtual

Casting the bnd exp to the specific case

Definition at line 1072 of file Mapping.cpp.

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

1073 {
1074  int physTot = m_fields[0]->GetTotPoints();
1075  int nvel = m_nConvectiveFields;
1076  int nfields = m_fields.num_elements();
1077  int nbnds = m_fields[0]->GetBndConditions().num_elements();
1078 
1079  // Declare variables
1080  Array<OneD, int> BCtoElmtID;
1081  Array<OneD, int> BCtoTraceID;
1082  Array<OneD, const SpatialDomains::BoundaryConditionShPtr> BndConds;
1083  Array<OneD, MultiRegions::ExpListSharedPtr> BndExp;
1086 
1087  Array<OneD, NekDouble> ElmtVal(physTot, 0.0);
1088  Array<OneD, NekDouble> BndVal(physTot, 0.0);
1089  Array<OneD, NekDouble> coordVelElmt(physTot, 0.0);
1090  Array<OneD, NekDouble> coordVelBnd(physTot, 0.0);
1091  Array<OneD, NekDouble> Vals(physTot, 0.0);
1092 
1093  Array<OneD, bool> isDirichlet(nfields);
1094 
1095  Array<OneD, Array<OneD, NekDouble> > values(nfields);
1096  for (int i=0; i < nfields; i++)
1097  {
1098  values[i] = Array<OneD, NekDouble> (physTot, 0.0);
1099  }
1100 
1101  Array<OneD, Array<OneD, NekDouble> > tmp(nvel);
1102  Array<OneD, Array<OneD, NekDouble> > tmp2(nvel);
1103  Array<OneD, Array<OneD, NekDouble> > coordVel(nvel);
1104  for (int i = 0; i< nvel; i++)
1105  {
1106  tmp[i] = Array<OneD, NekDouble> (physTot, 0.0);
1107  tmp2[i] = Array<OneD, NekDouble> (physTot, 0.0);
1108  coordVel[i] = Array<OneD, NekDouble> (physTot, 0.0);
1109  }
1110 
1111  // Get coordinates velocity in transformed system (for MovingBody regions)
1112  GetCoordVelocity(tmp);
1113  ContravarFromCartesian(tmp, coordVel);
1114 
1115  // Get Cartesian coordinates for evaluating boundary conditions
1116  Array<OneD, Array<OneD, NekDouble> > coords(3);
1117  for (int dir=0; dir < 3; dir++)
1118  {
1119  coords[dir] = Array<OneD, NekDouble> (physTot, 0.0);
1120  }
1121  GetCartesianCoordinates(coords[0],coords[1],coords[2]);
1122 
1123  // Loop boundary conditions looking for Dirichlet bc's
1124  for(int n = 0 ; n < nbnds ; ++n)
1125  {
1126  // Evaluate original Dirichlet boundary conditions in whole domain
1127  for (int i = 0; i < nfields; ++i)
1128  {
1129  BndConds = m_fields[i]->GetBndConditions();
1130  BndExp = m_fields[i]->GetBndCondExpansions();
1131  if ( BndConds[n]->GetBoundaryConditionType() ==
1133  {
1134  isDirichlet[i] = true;
1135  // If we have the a velocity component
1136  // check if all vel bc's are also Dirichlet
1137  if ( i<nvel )
1138  {
1139  for (int j = 0; j < nvel; ++j)
1140  {
1141  ASSERTL0(m_fields[j]->GetBndConditions()[n]->
1142  GetBoundaryConditionType() ==
1144  "Mapping only supported when all velocity components have the same type of boundary conditions");
1145  }
1146  }
1147  // Check if bc is time-dependent
1148  ASSERTL0( !BndConds[n]->IsTimeDependent(),
1149  "Time-dependent Dirichlet boundary conditions not supported with mapping yet.");
1150 
1151  // Get boundary condition
1152  LibUtilities::Equation condition =
1153  boost::static_pointer_cast<
1154  SpatialDomains::DirichletBoundaryCondition>
1155  (BndConds[n])->
1156  m_dirichletCondition;
1157  // Evaluate
1158  condition.Evaluate(coords[0], coords[1], coords[2],
1159  time, values[i]);
1160  }
1161  else
1162  {
1163  isDirichlet[i] = false;
1164  }
1165  }
1166  // Convert velocity vector to transformed system
1167  if ( isDirichlet[0])
1168  {
1169  for (int i = 0; i < nvel; ++i)
1170  {
1171  Vmath::Vcopy(physTot, values[i], 1, tmp[i], 1);
1172  }
1173  ContravarFromCartesian(tmp, tmp2);
1174  for (int i = 0; i < nvel; ++i)
1175  {
1176  Vmath::Vcopy(physTot, tmp2[i], 1, values[i], 1);
1177  }
1178  }
1179 
1180  // Now, project result to boundary
1181  for (int i = 0; i < nfields; ++i)
1182  {
1183  BndConds = m_fields[i]->GetBndConditions();
1184  BndExp = m_fields[i]->GetBndCondExpansions();
1185 
1186  // Loop boundary conditions again to get correct
1187  // values for cnt
1188  int cnt = 0;
1189  for(int m = 0 ; m < nbnds; ++m)
1190  {
1191  int exp_size = BndExp[m]->GetExpSize();
1192  if (m==n && isDirichlet[i])
1193  {
1194  for (int j = 0; j < exp_size; ++j, cnt++)
1195  {
1196  m_fields[i]->GetBoundaryToElmtMap(BCtoElmtID,
1197  BCtoTraceID);
1198  /// Casting the bnd exp to the specific case
1199  Bc = boost::dynamic_pointer_cast<
1200  StdRegions::StdExpansion>
1201  (BndExp[n]->GetExp(j));
1202  // Get element expansion
1203  elmt = m_fields[i]->GetExp(BCtoElmtID[cnt]);
1204  // Get values on the element
1205  ElmtVal = values[i] +
1206  m_fields[i]->GetPhys_Offset(
1207  BCtoElmtID[cnt]);
1208  // Get values on boundary
1209  elmt->GetTracePhysVals(BCtoTraceID[cnt],
1210  Bc, ElmtVal, BndVal);
1211 
1212  // Pointer to value that should be updated
1213  Vals = BndExp[n]->UpdatePhys()
1214  + BndExp[n]->GetPhys_Offset(j);
1215 
1216  // Copy result
1217  Vmath::Vcopy(Bc->GetTotPoints(),
1218  BndVal, 1, Vals, 1);
1219 
1220  // Apply MovingBody correction
1221  if ( (i<nvel) &&
1222  BndConds[n]->GetUserDefined() ==
1223  "MovingBody" )
1224  {
1225  // get coordVel in the element
1226  coordVelElmt = coordVel[i] +
1227  m_fields[i]->GetPhys_Offset(
1228  BCtoElmtID[cnt]);
1229 
1230  // Get values on boundary
1231  elmt->GetTracePhysVals(
1232  BCtoTraceID[cnt], Bc,
1233  coordVelElmt, coordVelBnd);
1234 
1235  // Apply correction
1236  Vmath::Vadd(Bc->GetTotPoints(),
1237  coordVelBnd, 1,
1238  Vals, 1, Vals, 1);
1239  }
1240  }
1241  }
1242  else // setting if m!=n
1243  {
1244  cnt += exp_size;
1245  }
1246  }
1247  }
1248  }
1249 
1250  // Finally, perform FwdTrans in all fields
1251  for (int i = 0; i < m_fields.num_elements(); ++i)
1252  {
1253  // Get boundary condition information
1254  BndConds = m_fields[i]->GetBndConditions();
1255  BndExp = m_fields[i]->GetBndCondExpansions();
1256  for(int n = 0 ; n < BndConds.num_elements(); ++n)
1257  {
1258  if ( BndConds[n]->GetBoundaryConditionType() ==
1260  {
1261  BndExp[n]->FwdTrans_BndConstrained(BndExp[n]->GetPhys(),
1262  BndExp[n]->UpdateCoeffs());
1263  if (m_fields[i]->GetExpType() == MultiRegions::e3DH1D)
1264  {
1265  BndExp[n]->HomogeneousFwdTrans(BndExp[n]->GetCoeffs(),
1266  BndExp[n]->UpdateCoeffs());
1267  }
1268  }
1269  }
1270  }
1271 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
GLOBAL_MAPPING_EXPORT void GetCoordVelocity(Array< OneD, Array< OneD, NekDouble > > &outarray)
Obtain the velocity of the coordinates.
Definition: Mapping.h:245
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:134
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
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:546
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
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:299
GLOBAL_MAPPING_EXPORT bool IsTimeDependent()
Get flag defining if mapping is time-dependent.
Definition: Mapping.h:341
virtual GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::v_UpdateGeomInfo ( )
protectedpure virtual
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 1273 of file Mapping.cpp.

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

1277 {
1278  if (m_fromFunction)
1279  {
1280  std::string s_FieldStr;
1281  string fieldNames[3] = {"x", "y", "z"};
1282  string velFieldNames[3] = {"vx", "vy", "vz"};
1283  // Check if function from session file defines each component
1284  // and evaluate them, otherwise there is no need to update
1285  // coords
1286  for(int i = 0; i < 3; i++)
1287  {
1288  s_FieldStr = fieldNames[i];
1289  if ( m_session->DefinesFunction(m_funcName, s_FieldStr))
1290  {
1291  EvaluateFunction(m_fields, m_session, s_FieldStr, m_coords[i],
1292  m_funcName, time);
1293  if ( i==2 && m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
1294  {
1295  ASSERTL0 (false,
1296  "3DH1D does not support mapping in the z-direction.");
1297  }
1298  }
1299  s_FieldStr = velFieldNames[i];
1300  if ( m_session->DefinesFunction(m_velFuncName, s_FieldStr))
1301  {
1302  EvaluateFunction(m_fields, m_session, s_FieldStr,
1303  m_coordsVel[i], m_velFuncName, time);
1304  if ( i==2 && m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
1305  {
1306  ASSERTL0 (false,
1307  "3DH1D does not support mapping in the z-direction.");
1308  }
1309  }
1310  }
1311  }
1312  else
1313  {
1314  int physTot = m_fields[0]->GetTotPoints();
1315  int nvel = m_nConvectiveFields;
1316  // Copy coordinates
1317  for(int i = 0; i < 3; i++)
1318  {
1319  Vmath::Vcopy(physTot, coords[i], 1, m_coords[i], 1);
1320  }
1321 
1322  for(int i = 0; i < nvel; i++)
1323  {
1324  Vmath::Vcopy(physTot, coordsVel[i], 1, m_coordsVel[i], 1);
1325  }
1326  }
1327 
1328  // Update the information required by the specific mapping
1329  UpdateGeomInfo();
1330 }
std::string m_funcName
Name of the function containing the coordinates.
Definition: Mapping.h:420
std::string m_velFuncName
Name of the function containing the velocity of the coordinates.
Definition: Mapping.h:422
Array< OneD, Array< OneD, NekDouble > > m_coords
Array with the Cartesian coordinates.
Definition: Mapping.h:411
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
LibUtilities::SessionReaderSharedPtr m_session
Session reader.
Definition: Mapping.h:405
Array< OneD, Array< OneD, NekDouble > > m_coordsVel
Array with the velocity of the coordinates.
Definition: Mapping.h:413
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:403
bool m_fromFunction
Flag defining if the Mapping is defined by a function.
Definition: Mapping.h:431
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
GLOBAL_MAPPING_EXPORT void UpdateGeomInfo()
Recompute the metric terms of the Mapping.
Definition: Mapping.h:397
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 792 of file Mapping.cpp.

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

796 {
797  int physTot = m_fields[0]->GetTotPoints();
798  int nvel = m_nConvectiveFields;
799 
800  Array<OneD, Array<OneD, NekDouble> > tmp(nvel);
801 
802  // Set wavespace to false and store current value
803  bool wavespace = m_fields[0]->GetWaveSpace();
804  m_fields[0]->SetWaveSpace(false);
805 
806  // Calculate vector gradient wk2 = u^i_(,k) = du^i/dx^k + {i,jk}*u^j
807  ApplyChristoffelContravar(inarray, m_wk1);
808  for (int i=0; i< nvel; i++)
809  {
810  if(nvel == 2)
811  {
812  m_fields[0]->PhysDeriv(inarray[i],
813  m_wk2[i*nvel+0],
814  m_wk2[i*nvel+1]);
815  }
816  else
817  {
818  m_fields[0]->PhysDeriv(inarray[i],
819  m_wk2[i*nvel+0],
820  m_wk2[i*nvel+1],
821  m_wk2[i*nvel+2]);
822  }
823  for (int k=0; k< nvel; k++)
824  {
825  Vmath::Vadd(physTot,m_wk1[i*nvel+k],1,m_wk2[i*nvel+k],1,
826  m_wk1[i*nvel+k], 1);
827  }
828  }
829  // Calculate wk1 = A^(ij) = g^(jk)*u^i_(,k)
830  for (int i=0; i< nvel; i++)
831  {
832  for (int k=0; k< nvel; k++)
833  {
834  tmp[k] = m_wk1[i*nvel+k];
835  }
836  RaiseIndex(tmp, m_tmp);
837  for (int j=0; j<nvel; j++)
838  {
839  Vmath::Vcopy(physTot, m_tmp[j], 1, m_wk1[i*nvel+j], 1);
840  }
841  }
842  //
843  // Calculate L(U)^i = (A^(ij))_(,j) - alpha*d^2(u^i)/dx^jdx^j
844  //
845 
846  // Step 1 :
847  // d(A^(ij) - alpha*du^i/dx^j)/d(x^j)
848  for (int i=0; i< nvel; i++)
849  {
850  outarray[i] = Array<OneD, NekDouble>(physTot,0.0);
851  for (int j=0; j< nvel; j++)
852  {
853  Vmath::Smul(physTot, alpha, m_wk2[i*nvel+j], 1,
854  m_tmp[0], 1);
855  Vmath::Vsub(physTot, m_wk1[i*nvel+j], 1, m_tmp[0], 1,
856  m_tmp[0], 1);
857 
858  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[j],
859  m_tmp[0],
860  m_tmp[1]);
861  Vmath::Vadd(physTot,outarray[i],1,m_tmp[1],1,outarray[i], 1);
862  }
863  }
864 
865  // Step 2: d(A^(ij))/d(x^j) + {j,pj}*A^(ip)
866  for (int i=0; i< nvel; i++)
867  {
868  for (int p=0; p< nvel; p++)
869  {
870  tmp[p] = m_wk1[i*nvel+p];
871  }
873  for (int j=0; j< nvel; j++)
874  {
875  Vmath::Vadd(physTot,outarray[i],1,m_wk2[j*nvel+j],1,
876  outarray[i], 1);
877  }
878  }
879 
880  // Step 3: d(A^(ij))/d(x^j) + {j,pj}*A^(ip) + {i,pj} A^(pj)
881  for (int j=0; j< nvel; j++)
882  {
883  for (int p=0; p< nvel; p++)
884  {
885  tmp[p] = m_wk1[p*nvel+j];
886  }
888  for (int i=0; i< nvel; i++)
889  {
890  Vmath::Vadd(physTot,outarray[i], 1, m_wk2[i*nvel+j], 1,
891  outarray[i], 1);
892  }
893  }
894 
895  // Restore value of wavespace
896  m_fields[0]->SetWaveSpace(wavespace);
897 }
Array< OneD, Array< OneD, NekDouble > > m_wk2
Definition: Mapping.h:440
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
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:638
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:213
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:212
Array< OneD, Array< OneD, NekDouble > > m_tmp
Definition: Mapping.h:441
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:343
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
Array< OneD, Array< OneD, NekDouble > > m_wk1
Definition: Mapping.h:439
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
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:299
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 290 of file Mapping.h.

References v_VelocityLaplacian().

294  {
295  v_VelocityLaplacian( inarray, outarray, alpha);
296  }
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:792

Member Data Documentation

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

Array with the velocity of the coordinates.

Definition at line 413 of file Mapping.h.

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

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

Definition at line 409 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::MappingXYofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXYofXY::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingGeneral::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingTranslation::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXYofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXYofXY::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingGeneral::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingTranslation::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXYofZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXofZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingTranslation::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingGeneral::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXofZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingTranslation::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingGeneral::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXofZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingTranslation::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingGeneral::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXofZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingTranslation::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarToCartesian(), Nektar::GlobalMapping::MappingGeneral::v_CovarToCartesian(), v_CurlCurlField(), v_Divergence(), Nektar::GlobalMapping::MappingXofZ::v_DotGradJacobian(), Nektar::GlobalMapping::MappingXYofZ::v_DotGradJacobian(), Nektar::GlobalMapping::MappingXofXZ::v_DotGradJacobian(), Nektar::GlobalMapping::MappingTranslation::v_DotGradJacobian(), v_DotGradJacobian(), v_GetCartesianCoordinates(), v_GetCoordVelocity(), Nektar::GlobalMapping::MappingXofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingTranslation::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofXY::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingGeneral::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetJacobian(), Nektar::GlobalMapping::MappingXofXZ::v_GetJacobian(), Nektar::GlobalMapping::MappingXofZ::v_GetJacobian(), Nektar::GlobalMapping::MappingTranslation::v_GetJacobian(), Nektar::GlobalMapping::MappingXYofXY::v_GetJacobian(), Nektar::GlobalMapping::MappingGeneral::v_GetJacobian(), Nektar::GlobalMapping::MappingXYofZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingTranslation::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXYofXY::v_GetMetricTensor(), Nektar::GlobalMapping::MappingGeneral::v_GetMetricTensor(), v_gradgradU(), Nektar::GlobalMapping::MappingTranslation::v_InitObject(), v_InitObject(), Nektar::GlobalMapping::MappingXofXZ::v_LowerIndex(), Nektar::GlobalMapping::MappingXofZ::v_LowerIndex(), Nektar::GlobalMapping::MappingTranslation::v_LowerIndex(), v_LowerIndex(), Nektar::GlobalMapping::MappingXofZ::v_RaiseIndex(), Nektar::GlobalMapping::MappingXofXZ::v_RaiseIndex(), Nektar::GlobalMapping::MappingTranslation::v_RaiseIndex(), v_RaiseIndex(), v_UpdateBCs(), Nektar::GlobalMapping::MappingXYofZ::v_UpdateGeomInfo(), Nektar::GlobalMapping::MappingXYofXY::v_UpdateGeomInfo(), Nektar::GlobalMapping::MappingXofZ::v_UpdateGeomInfo(), Nektar::GlobalMapping::MappingXofXZ::v_UpdateGeomInfo(), v_UpdateMapping(), and v_VelocityLaplacian().

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

Definition at line 407 of file Mapping.h.

Referenced by Mapping(), and Output().

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

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

Name of the function containing the coordinates.

Definition at line 420 of file Mapping.h.

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

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

Array with metric terms of the mapping.

Definition at line 415 of file Mapping.h.

Referenced by Nektar::GlobalMapping::MappingXYofXY::CalculateMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXYofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXYofZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXofZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXofZ::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXofZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXYofZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXofZ::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXofXZ::v_DotGradJacobian(), Nektar::GlobalMapping::MappingXofXZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_GetJacobian(), Nektar::GlobalMapping::MappingXYofXY::v_GetJacobian(), Nektar::GlobalMapping::MappingXofZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::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::MappingXYofZ::v_UpdateGeomInfo(), Nektar::GlobalMapping::MappingXYofXY::v_UpdateGeomInfo(), Nektar::GlobalMapping::MappingXofXZ::v_UpdateGeomInfo(), and Nektar::GlobalMapping::MappingXofZ::v_UpdateGeomInfo().

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

Definition at line 435 of file Mapping.h.

Referenced by Load().

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

Definition at line 436 of file Mapping.h.

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

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

Definition at line 434 of file Mapping.h.

Referenced by Load().

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

Number of velocity components.

Definition at line 417 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::MappingXYofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXYofXY::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingGeneral::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingTranslation::v_ApplyChristoffelContravar(), Nektar::GlobalMapping::MappingXYofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXYofXY::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofXZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingXofZ::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingGeneral::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingTranslation::v_ApplyChristoffelCovar(), Nektar::GlobalMapping::MappingTranslation::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingGeneral::v_ContravarFromCartesian(), Nektar::GlobalMapping::MappingTranslation::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingGeneral::v_ContravarToCartesian(), Nektar::GlobalMapping::MappingTranslation::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingGeneral::v_CovarFromCartesian(), Nektar::GlobalMapping::MappingTranslation::v_CovarToCartesian(), Nektar::GlobalMapping::MappingXYofXY::v_CovarToCartesian(), Nektar::GlobalMapping::MappingGeneral::v_CovarToCartesian(), v_CurlCurlField(), v_Divergence(), v_DotGradJacobian(), v_GetCoordVelocity(), Nektar::GlobalMapping::MappingXofXZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingTranslation::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofXY::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingGeneral::v_GetInvMetricTensor(), Nektar::GlobalMapping::MappingXYofZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofXZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXofZ::v_GetMetricTensor(), Nektar::GlobalMapping::MappingTranslation::v_GetMetricTensor(), Nektar::GlobalMapping::MappingXYofXY::v_GetMetricTensor(), Nektar::GlobalMapping::MappingGeneral::v_GetMetricTensor(), v_gradgradU(), Nektar::GlobalMapping::MappingXofXZ::v_InitObject(), Nektar::GlobalMapping::MappingXofZ::v_InitObject(), Nektar::GlobalMapping::MappingXYofZ::v_InitObject(), Nektar::GlobalMapping::MappingTranslation::v_InitObject(), Nektar::GlobalMapping::MappingXYofXY::v_InitObject(), Nektar::GlobalMapping::MappingGeneral::v_InitObject(), v_InitObject(), Nektar::GlobalMapping::MappingTranslation::v_LowerIndex(), v_LowerIndex(), Nektar::GlobalMapping::MappingTranslation::v_RaiseIndex(), v_RaiseIndex(), v_UpdateBCs(), v_UpdateMapping(), and v_VelocityLaplacian().

LibUtilities::SessionReaderSharedPtr Nektar::GlobalMapping::Mapping::m_session
protected
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(), Nektar::GlobalMapping::MappingTranslation::v_InitObject(), and v_InitObject().

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

Name of the function containing the velocity of the coordinates.

Definition at line 422 of file Mapping.h.

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

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

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