Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator 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 (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const bool generalized)
 CurlCurl calculated on the whole field. More...
 
GLOBAL_MAPPING_EXPORT bool IsTimeDependent ()
 Get flag defining if mapping is time-dependent. More...
 
GLOBAL_MAPPING_EXPORT void SetTimeDependent (const bool value)
 Set flag defining if mapping is time-dependent. More...
 
GLOBAL_MAPPING_EXPORT bool IsFromFunction ()
 Get flag defining if mapping is defined by a function. More...
 
GLOBAL_MAPPING_EXPORT void SetFromFunction (const bool value)
 Set flag defining if mapping is defined by a function. More...
 
GLOBAL_MAPPING_EXPORT bool HasConstantJacobian ()
 Get flag defining if mapping has constant Jacobian. More...
 
GLOBAL_MAPPING_EXPORT bool IsDefined ()
 Get flag determining if the mapping was defined or is trivial. More...
 
GLOBAL_MAPPING_EXPORT void UpdateBCs (const NekDouble time)
 Update the Dirichlet Boundary Conditions when using Mappings. More...
 
GLOBAL_MAPPING_EXPORT void UpdateMapping (const NekDouble time, const Array< OneD, Array< OneD, NekDouble > > &coords=NullNekDoubleArrayofArray, const Array< OneD, Array< OneD, NekDouble > > &coordsVel=NullNekDoubleArrayofArray)
 Update the Mapping with new coordinates. More...
 
GLOBAL_MAPPING_EXPORT void UpdateGeomInfo ()
 Recompute the metric terms of the Mapping. More...
 

Static 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 (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const bool generalized)
 
virtual GLOBAL_MAPPING_EXPORT void v_UpdateMapping (const NekDouble time, const Array< OneD, Array< OneD, NekDouble > > &coords=NullNekDoubleArrayofArray, const Array< OneD, Array< OneD, NekDouble > > &coordsVel=NullNekDoubleArrayofArray)
 
virtual GLOBAL_MAPPING_EXPORT void v_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 Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, 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  (pSession->GetComm());
90 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
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

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

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

Referenced by v_UpdateBCs().

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

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

491 {
492  if(inarray == outarray)
493  {
494  int physTot = m_fields[0]->GetTotPoints();
495  int nvel = m_nConvectiveFields;
496 
497  for (int i=0; i< nvel; i++)
498  {
499  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
500  }
501  v_ContravarToCartesian( m_tmp, outarray);
502  }
503  else
504  {
505  v_ContravarToCartesian( inarray, outarray);
506  }
507 }
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:1047
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 578 of file Mapping.cpp.

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

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

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

521 {
522  if(inarray == outarray)
523  {
524  int physTot = m_fields[0]->GetTotPoints();
525  int nvel = m_nConvectiveFields;
526 
527  for (int i=0; i< nvel; i++)
528  {
529  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
530  }
531  v_CovarToCartesian( m_tmp, outarray);
532  }
533  else
534  {
535  v_CovarToCartesian( inarray, outarray);
536  }
537 }
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:1047
GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::Mapping::CurlCurlField ( const 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(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const bool generalized)
Definition: Mapping.cpp:1011
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:765
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:689
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 404 of file Mapping.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::LibUtilities::eFunctionTypeExpression, Nektar::LibUtilities::eFunctionTypeFile, m_fields, m_session, and Vmath::Zero().

Referenced by v_InitObject(), and v_UpdateMapping().

411 {
412  ASSERTL0(pSession->DefinesFunction(pFunctionName),
413  "Function '" + pFunctionName + "' does not exist.");
414 
415  unsigned int nq = pFields[0]->GetNpoints();
416  if (pArray.num_elements() != nq)
417  {
418  pArray = Array<OneD, NekDouble> (nq);
419  }
420 
422  vType = pSession->GetFunctionType(pFunctionName, pFieldName);
424  {
425  Array<OneD, NekDouble> x0(nq);
426  Array<OneD, NekDouble> x1(nq);
427  Array<OneD, NekDouble> x2(nq);
428 
429  pFields[0]->GetCoords(x0, x1, x2);
431  pSession->GetFunction(pFunctionName, pFieldName);
432 
433  ffunc->Evaluate(x0, x1, x2, pTime, pArray);
434  }
435  else if (vType == LibUtilities::eFunctionTypeFile)
436  {
437  std::string filename = pSession->GetFunctionFilename(
438  pFunctionName,
439  pFieldName);
440 
441  std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
442  std::vector<std::vector<NekDouble> > FieldData;
443  Array<OneD, NekDouble> vCoeffs(pFields[0]->GetNcoeffs());
444  Vmath::Zero(vCoeffs.num_elements(), vCoeffs, 1);
445 
448  (m_session->GetComm());
449  fld->Import(filename, FieldDef, FieldData);
450 
451  int idx = -1;
452  for (int i = 0; i < FieldDef.size(); ++i)
453  {
454  for (int j = 0; j < FieldDef[i]->m_fields.size(); ++j)
455  {
456  if (FieldDef[i]->m_fields[j] == pFieldName)
457  {
458  idx = j;
459  }
460  }
461 
462  if (idx >= 0)
463  {
464  pFields[0]->ExtractDataToCoeffs(
465  FieldDef[i],
466  FieldData[i],
467  FieldDef[i]->m_fields[idx],
468  vCoeffs);
469  }
470  else
471  {
472  cout << "Field " + pFieldName + " not found." << endl;
473  }
474  }
475  pFields[0]->BwdTrans_IterPerExp(vCoeffs, pArray);
476  }
477 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
LibUtilities::SessionReaderSharedPtr m_session
Session reader.
Definition: Mapping.h:405
boost::shared_ptr< FieldIO > FieldIOSharedPtr
Definition: FieldIO.h:225
boost::shared_ptr< Equation > EquationSharedPtr
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:359
void Nektar::GlobalMapping::Mapping::EvaluateTimeFunction ( LibUtilities::SessionReaderSharedPtr  pSession,
std::string  pFieldName,
Array< OneD, NekDouble > &  pArray,
const std::string &  pFunctionName,
NekDouble  pTime = NekDouble(0) 
)
protected

Definition at line 383 of file Mapping.cpp.

References ASSERTL0.

389 {
390  ASSERTL0(pSession->DefinesFunction(pFunctionName),
391  "Function '" + pFunctionName + "' does not exist.");
392 
394  pSession->GetFunction(pFunctionName, pFieldName);
395 
396  Array<OneD, NekDouble> x0(1,0.0);
397  Array<OneD, NekDouble> x1(1,0.0);
398  Array<OneD, NekDouble> x2(1,0.0);
399 
400  ffunc->Evaluate(x0, x1, x2, pTime, pArray);
401 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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:661
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:677
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:901
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:100
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 266 of file Mapping.cpp.

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

Referenced by Nektar::Utilities::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().

269 {
270  if (!m_init)
271  {
272  TiXmlElement* vMapping = NULL;
273  string vType;
274  if (pSession->DefinesElement("Nektar/Mapping"))
275  {
276  vMapping = pSession->GetElement("Nektar/Mapping");
277  vType = vMapping->Attribute("TYPE");
278  m_isDefined = true;
279  }
280  else
281  {
282  vType = "Translation";
283  }
284 
286  vType, pSession, pFields,
287  vMapping);
288 
289  m_init = true;
290  }
291 
292  return m_mappingPtr;
293 }
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 609 of file Mapping.cpp.

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

612 {
613  if(inarray == outarray)
614  {
615  int physTot = m_fields[0]->GetTotPoints();
616  int nvel = m_nConvectiveFields;
617 
618  for (int i=0; i< nvel; i++)
619  {
620  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
621  }
622  v_LowerIndex( m_tmp, outarray);
623  }
624  else
625  {
626  v_LowerIndex( inarray, outarray);
627  }
628 }
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:719
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
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 303 of file Mapping.cpp.

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

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

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

Referenced by v_CurlCurlField(), and v_VelocityLaplacian().

643 {
644  if(inarray == outarray)
645  {
646  int physTot = m_fields[0]->GetTotPoints();
647  int nvel = m_nConvectiveFields;
648 
649  for (int i=0; i< nvel; i++)
650  {
651  Vmath::Vcopy(physTot, inarray[i], 1, m_tmp[i], 1);
652  }
653  v_RaiseIndex( m_tmp, outarray);
654  }
655  else
656  {
657  v_RaiseIndex( inarray, outarray);
658  }
659 }
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:742
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
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 243 of file Mapping.cpp.

References InitObject(), m_fields, and m_session.

245 {
246  m_fields = pFields;
247 
248  TiXmlElement* vMapping = NULL;
249 
250  if (m_session->DefinesElement("Nektar/Mapping"))
251  {
252  vMapping = m_session->GetElement("Nektar/Mapping");
253  }
254  InitObject(pFields, vMapping);
255 }
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:1187
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:1388
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 ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const bool  generalized 
)
protectedvirtual

Definition at line 1011 of file Mapping.cpp.

References Nektar::MultiRegions::DirCartesianMap, gradgradU(), m_fields, m_nConvectiveFields, m_tmp, RaiseIndex(), Vmath::Smul(), Vmath::Vadd(), Vmath::Vcopy(), and Vmath::Vsub().

Referenced by CurlCurlField().

1015 {
1016  int physTot = m_fields[0]->GetTotPoints();
1017  int nvel = m_nConvectiveFields;
1018 
1019  // Set wavespace to false and store current value
1020  bool wavespace = m_fields[0]->GetWaveSpace();
1021  m_fields[0]->SetWaveSpace(false);
1022 
1023  // For implicit treatment of viscous terms, we want the generalized
1024  // curlcurl and for explicit treatment, we want the cartesian one.
1025  if (generalized)
1026  {
1027  // Get the second derivatives u^i_{,jk}
1028  Array<OneD, Array<OneD, NekDouble> > tmp(nvel);
1029  Array<OneD, Array<OneD, NekDouble> > ddU(nvel*nvel*nvel);
1030  gradgradU(inarray, ddU);
1031 
1032  // Raise index to obtain A^{ip}_{k} = g^pj u^i_{,jk}
1033  for (int i = 0; i < nvel; ++i)
1034  {
1035  for (int k = 0; k < nvel; ++k)
1036  {
1037  // Copy to wk
1038  for (int j = 0; j < nvel; ++j)
1039  {
1040  tmp[j] = ddU[i*nvel*nvel+j*nvel+k];
1041  }
1042  RaiseIndex(tmp, m_tmp);
1043  for (int p=0; p<nvel; ++p)
1044  {
1045  Vmath::Vcopy(physTot, m_tmp[p], 1,
1046  ddU[i*nvel*nvel+p*nvel+k], 1);
1047  }
1048  }
1049  }
1050  // The curlcurl is g^ji u^k_{kj} - g^jk u^i_kj = A^{ki}_k - A^{ik}_k
1051  for (int i = 0; i < nvel; ++i)
1052  {
1053  outarray[i] = Array<OneD, NekDouble> (physTot, 0.0);
1054  for (int k = 0; k < nvel; ++k)
1055  {
1056  Vmath::Vadd(physTot, outarray[i], 1,
1057  ddU[k*nvel*nvel+i*nvel+k], 1,
1058  outarray[i], 1);
1059  Vmath::Vsub(physTot, outarray[i], 1,
1060  ddU[i*nvel*nvel+k*nvel+k], 1,
1061  outarray[i], 1);
1062  }
1063  }
1064  }
1065  else
1066  {
1067  switch(nvel)
1068  {
1069  case 2:
1070  {
1071  Array<OneD,NekDouble> Vx(physTot);
1072  Array<OneD,NekDouble> Uy(physTot);
1073  Array<OneD,NekDouble> Dummy(physTot);
1074 
1075  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
1076  inarray[1], Vx);
1077  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
1078  inarray[0], Uy);
1079 
1080  Vmath::Vsub(physTot, Vx, 1, Uy, 1, Dummy, 1);
1081 
1082  m_fields[0]->PhysDeriv(Dummy,outarray[1],outarray[0]);
1083 
1084  Vmath::Smul(physTot, -1.0, outarray[1], 1, outarray[1], 1);
1085  }
1086  break;
1087 
1088  case 3:
1089  {
1090  // Declare variables
1091  Array<OneD,NekDouble> Ux(physTot);
1092  Array<OneD,NekDouble> Uy(physTot);
1093  Array<OneD,NekDouble> Uz(physTot);
1094 
1095  Array<OneD,NekDouble> Vx(physTot);
1096  Array<OneD,NekDouble> Vy(physTot);
1097  Array<OneD,NekDouble> Vz(physTot);
1098 
1099  Array<OneD,NekDouble> Wx(physTot);
1100  Array<OneD,NekDouble> Wy(physTot);
1101  Array<OneD,NekDouble> Wz(physTot);
1102 
1103  Array<OneD,NekDouble> Dummy1(physTot);
1104  Array<OneD,NekDouble> Dummy2(physTot);
1105 
1106  // Calculate gradient
1107  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
1108  inarray[0], Ux);
1109  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
1110  inarray[0], Uy);
1111  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],
1112  inarray[0], Uz);
1113 
1114  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
1115  inarray[1], Vx);
1116  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
1117  inarray[1], Vy);
1118  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],
1119  inarray[1], Vz);
1120 
1121  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
1122  inarray[2], Wx);
1123  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
1124  inarray[2], Wy);
1125  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],
1126  inarray[2], Wz);
1127 
1128  // x-component
1129  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
1130  Vx, Dummy1); //Vxy
1131  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
1132  Uy, Dummy2); //Uyy
1133  Vmath::Vsub(physTot, Dummy1, 1, Dummy2, 1, outarray[0], 1);
1134 
1135  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],
1136  Uz, Dummy1); //Uzz
1137  Vmath::Vsub(physTot, outarray[0], 1, Dummy1, 1,
1138  outarray[0], 1);
1139 
1140  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
1141  Wz, Dummy2); //Wxz
1142  Vmath::Vadd(physTot, outarray[0], 1, Dummy2, 1,
1143  outarray[0], 1);
1144 
1145  // y-component
1146  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
1147  Wz, Dummy1); //Wzy
1148  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],
1149  Vz, Dummy2); //Vzz
1150  Vmath::Vsub(physTot, Dummy1, 1, Dummy2, 1, outarray[1], 1);
1151 
1152  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
1153  Vx, Dummy1); //Vxx
1154  Vmath::Vsub(physTot, outarray[1], 1, Dummy1, 1,
1155  outarray[1], 1);
1156 
1157  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
1158  Uy, Dummy2); //Uyx
1159  Vmath::Vadd(physTot, outarray[1], 1, Dummy2, 1,
1160  outarray[1], 1);
1161 
1162  // z-component
1163  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
1164  Uz, Dummy1); //Uxz
1165  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
1166  Wy, Dummy2); //Wyy
1167  Vmath::Vsub(physTot, Dummy1, 1, Dummy2, 1, outarray[2], 1);
1168 
1169  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
1170  Wx, Dummy1); //Wxx
1171  Vmath::Vsub(physTot, outarray[2], 1, Dummy1, 1,
1172  outarray[2], 1);
1173 
1174  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
1175  Vz, Dummy2); //Vyz
1176  Vmath::Vadd(physTot, outarray[2], 1, Dummy2, 1,
1177  outarray[2], 1);
1178  }
1179  break;
1180  }
1181  }
1182 
1183  // Restore value of wavespace
1184  m_fields[0]->SetWaveSpace(wavespace);
1185 }
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:640
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:199
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:329
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:285
void Nektar::GlobalMapping::Mapping::v_Divergence ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Definition at line 765 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().

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

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

Referenced by DotGradJacobian().

692 {
693  int physTot = m_fields[0]->GetTotPoints();
694 
695  outarray = Array<OneD, NekDouble>(physTot, 0.0);
696  if ( !HasConstantJacobian() )
697  {
698  // Set wavespace to false and store current value
699  bool wavespace = m_fields[0]->GetWaveSpace();
700  m_fields[0]->SetWaveSpace(false);
701 
702  // Get Mapping Jacobian
703  Array<OneD, NekDouble> Jac(physTot, 0.0);
704  GetJacobian(Jac);
705 
706  // Calculate inarray . grad(Jac)
707  Array<OneD, NekDouble> wk(physTot, 0.0);
708  for(int i = 0; i < m_nConvectiveFields; ++i)
709  {
710  m_fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[i],
711  Jac, wk);
712  Vmath::Vvtvp(physTot, inarray[i], 1, wk, 1,
713  outarray, 1, outarray, 1);
714  }
715  m_fields[0]->SetWaveSpace(wavespace);
716  }
717 }
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:428
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 661 of file Mapping.cpp.

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

Referenced by GetCartesianCoordinates().

665 {
666  int physTot = m_fields[0]->GetTotPoints();
667 
668  out0 = Array<OneD, NekDouble>(physTot, 0.0);
669  out1 = Array<OneD, NekDouble>(physTot, 0.0);
670  out2 = Array<OneD, NekDouble>(physTot, 0.0);
671 
672  Vmath::Vcopy(physTot, m_coords[0], 1, out0, 1);
673  Vmath::Vcopy(physTot, m_coords[1], 1, out1, 1);
674  Vmath::Vcopy(physTot, m_coords[2], 1, out2, 1);
675 }
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:1047
void Nektar::GlobalMapping::Mapping::v_GetCoordVelocity ( Array< OneD, Array< OneD, NekDouble > > &  outarray)
protectedvirtual

Definition at line 677 of file Mapping.cpp.

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

Referenced by GetCoordVelocity().

679 {
680  int physTot = m_fields[0]->GetTotPoints();
681 
682  for(int i = 0; i < m_nConvectiveFields; ++i)
683  {
684  outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
685  Vmath::Vcopy(physTot, m_coordsVel[i], 1, outarray[i], 1);
686  }
687 }
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:1047
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 901 of file Mapping.cpp.

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

Referenced by gradgradU().

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

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

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

Referenced by LowerIndex().

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

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

Referenced by RaiseIndex().

745 {
746  int physTot = m_fields[0]->GetTotPoints();
747  int nvel = m_nConvectiveFields;
748 
749  Array<OneD, Array<OneD, NekDouble> > g(nvel*nvel);
750 
752 
753  for (int i=0; i< nvel; i++)
754  {
755  outarray[i] = Array<OneD, NekDouble> (physTot, 0.0);
756  for (int j=0; j< nvel; j++)
757  {
758  Vmath::Vvtvp(physTot, g[i*nvel+j], 1, inarray[j], 1,
759  outarray[i], 1,
760  outarray[i], 1);
761  }
762  }
763 }
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:428
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:417
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:409
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 1187 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().

1188 {
1189  int physTot = m_fields[0]->GetTotPoints();
1190  int nvel = m_nConvectiveFields;
1191  int nfields = m_fields.num_elements();
1192  int nbnds = m_fields[0]->GetBndConditions().num_elements();
1193 
1194  // Declare variables
1195  Array<OneD, int> BCtoElmtID;
1196  Array<OneD, int> BCtoTraceID;
1197  Array<OneD, const SpatialDomains::BoundaryConditionShPtr> BndConds;
1198  Array<OneD, MultiRegions::ExpListSharedPtr> BndExp;
1201 
1202  Array<OneD, NekDouble> ElmtVal(physTot, 0.0);
1203  Array<OneD, NekDouble> BndVal(physTot, 0.0);
1204  Array<OneD, NekDouble> coordVelElmt(physTot, 0.0);
1205  Array<OneD, NekDouble> coordVelBnd(physTot, 0.0);
1206  Array<OneD, NekDouble> Vals(physTot, 0.0);
1207 
1208  Array<OneD, bool> isDirichlet(nfields);
1209 
1210  Array<OneD, Array<OneD, NekDouble> > values(nfields);
1211  for (int i=0; i < nfields; i++)
1212  {
1213  values[i] = Array<OneD, NekDouble> (physTot, 0.0);
1214  }
1215 
1216  Array<OneD, Array<OneD, NekDouble> > tmp(nvel);
1217  Array<OneD, Array<OneD, NekDouble> > tmp2(nvel);
1218  Array<OneD, Array<OneD, NekDouble> > coordVel(nvel);
1219  for (int i = 0; i< nvel; i++)
1220  {
1221  tmp[i] = Array<OneD, NekDouble> (physTot, 0.0);
1222  tmp2[i] = Array<OneD, NekDouble> (physTot, 0.0);
1223  coordVel[i] = Array<OneD, NekDouble> (physTot, 0.0);
1224  }
1225 
1226  // Get coordinates velocity in transformed system (for MovingBody regions)
1227  GetCoordVelocity(tmp);
1228  ContravarFromCartesian(tmp, coordVel);
1229 
1230  // Get Cartesian coordinates for evaluating boundary conditions
1231  Array<OneD, Array<OneD, NekDouble> > coords(3);
1232  for (int dir=0; dir < 3; dir++)
1233  {
1234  coords[dir] = Array<OneD, NekDouble> (physTot, 0.0);
1235  }
1236  GetCartesianCoordinates(coords[0],coords[1],coords[2]);
1237 
1238  // Loop boundary conditions looking for Dirichlet bc's
1239  for(int n = 0 ; n < nbnds ; ++n)
1240  {
1241  // Evaluate original Dirichlet boundary conditions in whole domain
1242  for (int i = 0; i < nfields; ++i)
1243  {
1244  BndConds = m_fields[i]->GetBndConditions();
1245  BndExp = m_fields[i]->GetBndCondExpansions();
1246  if ( BndConds[n]->GetBoundaryConditionType() ==
1248  {
1249  isDirichlet[i] = true;
1250  // If we have the a velocity component
1251  // check if all vel bc's are also Dirichlet
1252  if ( i<nvel )
1253  {
1254  for (int j = 0; j < nvel; ++j)
1255  {
1256  ASSERTL0(m_fields[j]->GetBndConditions()[n]->
1257  GetBoundaryConditionType() ==
1259  "Mapping only supported when all velocity components have the same type of boundary conditions");
1260  }
1261  }
1262  // Check if bc is time-dependent
1263  ASSERTL0( !BndConds[n]->IsTimeDependent(),
1264  "Time-dependent Dirichlet boundary conditions not supported with mapping yet.");
1265 
1266  // Get boundary condition
1267  LibUtilities::Equation condition =
1268  boost::static_pointer_cast<
1269  SpatialDomains::DirichletBoundaryCondition>
1270  (BndConds[n])->
1271  m_dirichletCondition;
1272  // Evaluate
1273  condition.Evaluate(coords[0], coords[1], coords[2],
1274  time, values[i]);
1275  }
1276  else
1277  {
1278  isDirichlet[i] = false;
1279  }
1280  }
1281  // Convert velocity vector to transformed system
1282  if ( isDirichlet[0])
1283  {
1284  for (int i = 0; i < nvel; ++i)
1285  {
1286  Vmath::Vcopy(physTot, values[i], 1, tmp[i], 1);
1287  }
1288  ContravarFromCartesian(tmp, tmp2);
1289  for (int i = 0; i < nvel; ++i)
1290  {
1291  Vmath::Vcopy(physTot, tmp2[i], 1, values[i], 1);
1292  }
1293  }
1294 
1295  // Now, project result to boundary
1296  for (int i = 0; i < nfields; ++i)
1297  {
1298  BndConds = m_fields[i]->GetBndConditions();
1299  BndExp = m_fields[i]->GetBndCondExpansions();
1300 
1301  // Loop boundary conditions again to get correct
1302  // values for cnt
1303  int cnt = 0;
1304  for(int m = 0 ; m < nbnds; ++m)
1305  {
1306  int exp_size = BndExp[m]->GetExpSize();
1307  if (m==n && isDirichlet[i])
1308  {
1309  for (int j = 0; j < exp_size; ++j, cnt++)
1310  {
1311  m_fields[i]->GetBoundaryToElmtMap(BCtoElmtID,
1312  BCtoTraceID);
1313  /// Casting the bnd exp to the specific case
1314  Bc = boost::dynamic_pointer_cast<
1315  StdRegions::StdExpansion>
1316  (BndExp[n]->GetExp(j));
1317  // Get element expansion
1318  elmt = m_fields[i]->GetExp(BCtoElmtID[cnt]);
1319  // Get values on the element
1320  ElmtVal = values[i] +
1321  m_fields[i]->GetPhys_Offset(
1322  BCtoElmtID[cnt]);
1323  // Get values on boundary
1324  elmt->GetTracePhysVals(BCtoTraceID[cnt],
1325  Bc, ElmtVal, BndVal);
1326 
1327  // Pointer to value that should be updated
1328  Vals = BndExp[n]->UpdatePhys()
1329  + BndExp[n]->GetPhys_Offset(j);
1330 
1331  // Copy result
1332  Vmath::Vcopy(Bc->GetTotPoints(),
1333  BndVal, 1, Vals, 1);
1334 
1335  // Apply MovingBody correction
1336  if ( (i<nvel) &&
1337  BndConds[n]->GetUserDefined() ==
1338  "MovingBody" )
1339  {
1340  // get coordVel in the element
1341  coordVelElmt = coordVel[i] +
1342  m_fields[i]->GetPhys_Offset(
1343  BCtoElmtID[cnt]);
1344 
1345  // Get values on boundary
1346  elmt->GetTracePhysVals(
1347  BCtoTraceID[cnt], Bc,
1348  coordVelElmt, coordVelBnd);
1349 
1350  // Apply correction
1351  Vmath::Vadd(Bc->GetTotPoints(),
1352  coordVelBnd, 1,
1353  Vals, 1, Vals, 1);
1354  }
1355  }
1356  }
1357  else // setting if m!=n
1358  {
1359  cnt += exp_size;
1360  }
1361  }
1362  }
1363  }
1364 
1365  // Finally, perform FwdTrans in all fields
1366  for (int i = 0; i < m_fields.num_elements(); ++i)
1367  {
1368  // Get boundary condition information
1369  BndConds = m_fields[i]->GetBndConditions();
1370  BndExp = m_fields[i]->GetBndCondExpansions();
1371  for(int n = 0 ; n < BndConds.num_elements(); ++n)
1372  {
1373  if ( BndConds[n]->GetBoundaryConditionType() ==
1375  {
1376  BndExp[n]->FwdTrans_BndConstrained(BndExp[n]->GetPhys(),
1377  BndExp[n]->UpdateCoeffs());
1378  if (m_fields[i]->GetExpType() == MultiRegions::e3DH1D)
1379  {
1380  BndExp[n]->HomogeneousFwdTrans(BndExp[n]->GetCoeffs(),
1381  BndExp[n]->UpdateCoeffs());
1382  }
1383  }
1384  }
1385  }
1386 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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:548
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:285
GLOBAL_MAPPING_EXPORT 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 1388 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().

1392 {
1393  if (m_fromFunction)
1394  {
1395  std::string s_FieldStr;
1396  string fieldNames[3] = {"x", "y", "z"};
1397  string velFieldNames[3] = {"vx", "vy", "vz"};
1398  // Check if function from session file defines each component
1399  // and evaluate them, otherwise there is no need to update
1400  // coords
1401  for(int i = 0; i < 3; i++)
1402  {
1403  s_FieldStr = fieldNames[i];
1404  if ( m_session->DefinesFunction(m_funcName, s_FieldStr))
1405  {
1406  EvaluateFunction(m_fields, m_session, s_FieldStr, m_coords[i],
1407  m_funcName, time);
1408  if ( i==2 && m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
1409  {
1410  ASSERTL0 (false,
1411  "3DH1D does not support mapping in the z-direction.");
1412  }
1413  }
1414  s_FieldStr = velFieldNames[i];
1415  if ( m_session->DefinesFunction(m_velFuncName, s_FieldStr))
1416  {
1417  EvaluateFunction(m_fields, m_session, s_FieldStr,
1418  m_coordsVel[i], m_velFuncName, time);
1419  if ( i==2 && m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
1420  {
1421  ASSERTL0 (false,
1422  "3DH1D does not support mapping in the z-direction.");
1423  }
1424  }
1425  }
1426  }
1427  else
1428  {
1429  int physTot = m_fields[0]->GetTotPoints();
1430  int nvel = m_nConvectiveFields;
1431  // Copy coordinates
1432  for(int i = 0; i < 3; i++)
1433  {
1434  Vmath::Vcopy(physTot, coords[i], 1, m_coords[i], 1);
1435  }
1436 
1437  for(int i = 0; i < nvel; i++)
1438  {
1439  Vmath::Vcopy(physTot, coordsVel[i], 1, m_coordsVel[i], 1);
1440  }
1441  }
1442 
1443  // Update the information required by the specific mapping
1444  UpdateGeomInfo();
1445 }
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:161
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:404
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:1047
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 794 of file Mapping.cpp.

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

Referenced by VelocityLaplacian().

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

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