Nektar++
Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Friends | List of all members
Nektar::GlobalMapping::MappingGeneral Class Reference

#include <MappingGeneral.h>

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

Static Public Member Functions

static GLOBAL_MAPPING_EXPORT MappingSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
 Creates an instance of this class. More...
 
- Static Public Member Functions inherited from Nektar::GlobalMapping::Mapping
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...
 

Static Public Attributes

static std::string className
 Name of the class. More...
 

Protected Member Functions

 MappingGeneral (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
 
GLOBAL_MAPPING_EXPORT void v_InitObject (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping) override
 
GLOBAL_MAPPING_EXPORT void v_ContravarToCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
GLOBAL_MAPPING_EXPORT void v_CovarToCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
GLOBAL_MAPPING_EXPORT void v_ContravarFromCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
GLOBAL_MAPPING_EXPORT void v_CovarFromCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
GLOBAL_MAPPING_EXPORT void v_GetJacobian (Array< OneD, NekDouble > &outarray) override
 
GLOBAL_MAPPING_EXPORT void v_GetMetricTensor (Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
GLOBAL_MAPPING_EXPORT void v_GetInvMetricTensor (Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelContravar (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelCovar (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
GLOBAL_MAPPING_EXPORT void v_UpdateGeomInfo () override
 
void CalculateMetricTerms ()
 
void CalculateChristoffel ()
 
- Protected Member Functions inherited from Nektar::GlobalMapping::Mapping
GLOBAL_MAPPING_EXPORT Mapping (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
 Constructor. More...
 
GLOBAL_MAPPING_EXPORT void EvaluateFunction (Array< OneD, MultiRegions::ExpListSharedPtr > pFields, LibUtilities::SessionReaderSharedPtr pSession, std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, NekDouble pTime=NekDouble(0))
 
GLOBAL_MAPPING_EXPORT void EvaluateTimeFunction (LibUtilities::SessionReaderSharedPtr pSession, std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, NekDouble pTime=NekDouble(0))
 
virtual GLOBAL_MAPPING_EXPORT void v_InitObject (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
 
virtual GLOBAL_MAPPING_EXPORT void v_ContravarToCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_CovarToCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_ContravarFromCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_CovarFromCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_GetCartesianCoordinates (Array< OneD, NekDouble > &out0, Array< OneD, NekDouble > &out1, Array< OneD, NekDouble > &out2)
 
virtual GLOBAL_MAPPING_EXPORT void v_GetCoordVelocity (Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_GetJacobian (Array< OneD, NekDouble > &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_DotGradJacobian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_GetMetricTensor (Array< OneD, Array< OneD, NekDouble > > &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_GetInvMetricTensor (Array< OneD, Array< OneD, NekDouble > > &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_LowerIndex (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_RaiseIndex (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelContravar (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_ApplyChristoffelCovar (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)=0
 
virtual GLOBAL_MAPPING_EXPORT void v_Divergence (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_VelocityLaplacian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble alpha)
 
virtual GLOBAL_MAPPING_EXPORT void v_gradgradU (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual GLOBAL_MAPPING_EXPORT void v_CurlCurlField (Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const bool generalized)
 
virtual GLOBAL_MAPPING_EXPORT void v_UpdateMapping (const NekDouble time, const Array< OneD, Array< OneD, NekDouble > > &coords=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble > > &coordsVel=NullNekDoubleArrayOfArray)
 
virtual GLOBAL_MAPPING_EXPORT void v_UpdateGeomInfo ()=0
 
virtual GLOBAL_MAPPING_EXPORT void v_UpdateBCs (const NekDouble time)
 

Protected Attributes

Array< OneD, Array< OneD, NekDouble > > m_metricTensor
 
Array< OneD, Array< OneD, NekDouble > > m_invMetricTensor
 
Array< OneD, Array< OneD, NekDouble > > m_deriv
 
Array< OneD, Array< OneD, NekDouble > > m_invDeriv
 
Array< OneD, Array< OneD, NekDouble > > m_Christoffel
 
Array< OneD, NekDoublem_jac
 
- Protected Attributes inherited from Nektar::GlobalMapping::Mapping
LibUtilities::SessionReaderSharedPtr m_session
 Session reader. More...
 
LibUtilities::FieldIOSharedPtr m_fld
 
Array< OneD, MultiRegions::ExpListSharedPtrm_fields
 
Array< OneD, Array< OneD, NekDouble > > m_coords
 Array with the Cartesian coordinates. More...
 
Array< OneD, Array< OneD, NekDouble > > m_coordsVel
 Array with the velocity of the coordinates. More...
 
Array< OneD, Array< OneD, NekDouble > > m_GeometricInfo
 Array with metric terms of the mapping. More...
 
int m_nConvectiveFields
 Number of velocity components. More...
 
std::string m_funcName
 Name of the function containing the coordinates. More...
 
std::string m_velFuncName
 Name of the function containing the velocity of the coordinates. More...
 
bool m_constantJacobian
 Flag defining if the Jacobian is constant. More...
 
bool m_timeDependent
 Flag defining if the Mapping is time-dependent. More...
 
bool m_fromFunction
 Flag defining if the Mapping is defined by a function. More...
 
Array< OneD, Array< OneD, NekDouble > > m_wk1
 
Array< OneD, Array< OneD, NekDouble > > m_wk2
 
Array< OneD, Array< OneD, NekDouble > > m_tmp
 

Friends

class MemoryManager< MappingGeneral >
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::GlobalMapping::Mapping
virtual GLOBAL_MAPPING_EXPORT ~Mapping ()
 Destructor. More...
 
GLOBAL_MAPPING_EXPORT void InitObject (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
 Initialise the mapping object. More...
 
GLOBAL_MAPPING_EXPORT void ReplaceField (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
 Replace the Expansion List used by the mapping. More...
 
GLOBAL_MAPPING_EXPORT void Output (LibUtilities::FieldMetaDataMap &fieldMetaDataMap, const std::string &outname)
 Output function called when a chk or fld file is written. More...
 
GLOBAL_MAPPING_EXPORT void ContravarToCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Convert a contravariant vector to the Cartesian system. More...
 
GLOBAL_MAPPING_EXPORT void CovarToCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Convert a covariant vector to the Cartesian system. More...
 
GLOBAL_MAPPING_EXPORT void ContravarFromCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Convert a contravariant vector to the transformed system. More...
 
GLOBAL_MAPPING_EXPORT void CovarFromCartesian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Convert a covariant vector to the transformed system. More...
 
GLOBAL_MAPPING_EXPORT void GetCartesianCoordinates (Array< OneD, NekDouble > &out0, Array< OneD, NekDouble > &out1, Array< OneD, NekDouble > &out2)
 Get the Cartesian coordinates in the field. More...
 
GLOBAL_MAPPING_EXPORT void GetJacobian (Array< OneD, NekDouble > &outarray)
 Get the Jacobian of the transformation. More...
 
GLOBAL_MAPPING_EXPORT void DotGradJacobian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
 Calculate the dot product with the gradient of the Jacobian. More...
 
GLOBAL_MAPPING_EXPORT void GetMetricTensor (Array< OneD, Array< OneD, NekDouble > > &outarray)
 Get the metric tensor \(g_{ij}\). More...
 
GLOBAL_MAPPING_EXPORT void GetInvMetricTensor (Array< OneD, Array< OneD, NekDouble > > &outarray)
 Get the inverse of metric tensor \(g^{ij}\). More...
 
GLOBAL_MAPPING_EXPORT void LowerIndex (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Lower index of vector: \(v_{i} = g_{ij}*v^{j}\). More...
 
GLOBAL_MAPPING_EXPORT void RaiseIndex (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Raise index of vector: \(v^{i} = g^{ij}*v_{j}\). More...
 
GLOBAL_MAPPING_EXPORT void ApplyChristoffelContravar (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Apply the Christoffel symbols to a contravariant vector. More...
 
GLOBAL_MAPPING_EXPORT void ApplyChristoffelCovar (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Apply the Christoffel symbols to a covariant vector. More...
 
GLOBAL_MAPPING_EXPORT void GetCoordVelocity (Array< OneD, Array< OneD, NekDouble > > &outarray)
 Obtain the velocity of the coordinates. More...
 
GLOBAL_MAPPING_EXPORT void Divergence (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
 Calculate the generalised divergence operator. More...
 
GLOBAL_MAPPING_EXPORT void VelocityLaplacian (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble alpha=0.0)
 Generalised (correction to the) velocity Laplacian operator. More...
 
GLOBAL_MAPPING_EXPORT void gradgradU (const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Second order covariant derivatives of a contravariant vector. More...
 
GLOBAL_MAPPING_EXPORT void CurlCurlField (Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const bool generalized)
 CurlCurl calculated on the whole field. More...
 
GLOBAL_MAPPING_EXPORT bool IsTimeDependent ()
 Get flag defining if mapping is time-dependent. More...
 
GLOBAL_MAPPING_EXPORT void SetTimeDependent (const bool value)
 Set flag defining if mapping is time-dependent. More...
 
GLOBAL_MAPPING_EXPORT bool IsFromFunction ()
 Get flag defining if mapping is defined by a function. More...
 
GLOBAL_MAPPING_EXPORT void SetFromFunction (const bool value)
 Set flag defining if mapping is defined by a function. More...
 
GLOBAL_MAPPING_EXPORT bool HasConstantJacobian ()
 Get flag defining if mapping has constant Jacobian. More...
 
GLOBAL_MAPPING_EXPORT bool IsDefined ()
 Get flag determining if the mapping was defined or is trivial. More...
 
GLOBAL_MAPPING_EXPORT void UpdateBCs (const NekDouble time)
 Update the Dirichlet Boundary Conditions when using Mappings. More...
 
GLOBAL_MAPPING_EXPORT void UpdateMapping (const NekDouble time, const Array< OneD, Array< OneD, NekDouble > > &coords=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble > > &coordsVel=NullNekDoubleArrayOfArray)
 Update the Mapping with new coordinates. More...
 
GLOBAL_MAPPING_EXPORT void UpdateGeomInfo ()
 Recompute the metric terms of the Mapping. More...
 
- Static Protected Attributes inherited from Nektar::GlobalMapping::Mapping
static MappingSharedPtr m_mappingPtr = MappingSharedPtr()
 
static bool m_init = false
 
static bool m_isDefined = false
 

Detailed Description

This class implements the most general mapping, defined by the transformation

\[ \bar{x} = \bar{x}(x,y,z) \]

\[ \bar{y} = \bar{y}(x,y,z) \]

\[ \bar{z} = \bar{z}(x,y,z) \]

where \((\bar{x},\bar{y},\bar{z})\) are the Cartesian (physical) coordinates and \((x,y,z)\) are the transformed (computational) coordinates.

Definition at line 43 of file MappingGeneral.h.

Constructor & Destructor Documentation

◆ MappingGeneral()

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

Definition at line 56 of file MappingGeneral.cpp.

59 : Mapping(pSession, pFields)
60{
61}
GLOBAL_MAPPING_EXPORT Mapping(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
Constructor.
Definition: Mapping.cpp:53

Member Function Documentation

◆ CalculateChristoffel()

void Nektar::GlobalMapping::MappingGeneral::CalculateChristoffel ( )
protected

Definition at line 375 of file MappingGeneral.cpp.

376{
377 int physTot = m_fields[0]->GetTotPoints();
378 int nvel = m_nConvectiveFields;
379
380 Array<OneD, Array<OneD, NekDouble>> gradG(nvel * nvel * nvel);
381 Array<OneD, Array<OneD, NekDouble>> tmp(nvel * nvel * nvel);
382 m_Christoffel = Array<OneD, Array<OneD, NekDouble>>(nvel * nvel * nvel);
383 // Allocate memory
384 for (int i = 0; i < gradG.size(); i++)
385 {
386 gradG[i] = Array<OneD, NekDouble>(physTot, 0.0);
387 tmp[i] = Array<OneD, NekDouble>(physTot, 0.0);
388 m_Christoffel[i] = Array<OneD, NekDouble>(physTot, 0.0);
389 }
390
391 // Set wavespace to false and store current value
392 bool waveSpace = m_fields[0]->GetWaveSpace();
393 m_fields[0]->SetWaveSpace(false);
394
395 // Calculate gradients of g_ij
396 for (int i = 0; i < nvel; i++)
397 {
398 for (int j = 0; j < nvel; j++)
399 {
400 for (int k = 0; k < nvel; k++)
401 {
403 m_metricTensor[i * nvel + j],
404 gradG[i * nvel * nvel + j * nvel + k]);
405 }
406 }
407 }
408
409 // Calculate tmp[p,j,k] = 1/2( gradG[pj,k]+ gradG[pk,j]-gradG[jk,p])
410 for (int p = 0; p < nvel; p++)
411 {
412 for (int j = 0; j < nvel; j++)
413 {
414 for (int k = 0; k < nvel; k++)
415 {
416 Vmath::Vadd(physTot, gradG[p * nvel * nvel + j * nvel + k], 1,
417 gradG[p * nvel * nvel + k * nvel + j], 1,
418 tmp[p * nvel * nvel + j * nvel + k], 1);
419 Vmath::Vsub(physTot, tmp[p * nvel * nvel + j * nvel + k], 1,
420 gradG[j * nvel * nvel + k * nvel + p], 1,
421 tmp[p * nvel * nvel + j * nvel + k], 1);
422 Vmath::Smul(physTot, 0.5, tmp[p * nvel * nvel + j * nvel + k],
423 1, tmp[p * nvel * nvel + j * nvel + k], 1);
424 }
425 }
426 }
427
428 // Calculate Christoffel symbols = g^ip tmp[p,j,k]
429 for (int i = 0; i < nvel; i++)
430 {
431 for (int j = 0; j < nvel; j++)
432 {
433 for (int k = 0; k < nvel; k++)
434 {
435 for (int p = 0; p < nvel; p++)
436 {
438 physTot, m_invMetricTensor[i * nvel + p], 1,
439 tmp[p * nvel * nvel + j * nvel + k], 1,
440 m_Christoffel[i * nvel * nvel + j * nvel + k], 1,
441 m_Christoffel[i * nvel * nvel + j * nvel + k], 1);
442 }
443 }
444 }
445 }
446 // Restore wavespace
447 m_fields[0]->SetWaveSpace(waveSpace);
448}
Array< OneD, Array< OneD, NekDouble > > m_metricTensor
Array< OneD, Array< OneD, NekDouble > > m_invMetricTensor
Array< OneD, Array< OneD, NekDouble > > m_Christoffel
int m_nConvectiveFields
Number of velocity components.
Definition: Mapping.h:418
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Mapping.h:410
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:87
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.hpp:366
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.hpp:180
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.hpp:100
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
Definition: Vmath.hpp:220

References Nektar::MultiRegions::DirCartesianMap, m_Christoffel, Nektar::GlobalMapping::Mapping::m_fields, m_invMetricTensor, m_metricTensor, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, CellMLToNektar.cellml_metadata::p, Vmath::Smul(), Vmath::Vadd(), Vmath::Vsub(), and Vmath::Vvtvp().

Referenced by v_UpdateGeomInfo().

◆ CalculateMetricTerms()

void Nektar::GlobalMapping::MappingGeneral::CalculateMetricTerms ( )
protected

Definition at line 241 of file MappingGeneral.cpp.

242{
243 int physTot = m_fields[0]->GetTotPoints();
244 int nvel = m_nConvectiveFields;
245
246 // Set wavespace to false and store current value
247 bool wavespace = m_fields[0]->GetWaveSpace();
248 m_fields[0]->SetWaveSpace(false);
249
250 // Allocate memory
251 m_metricTensor = Array<OneD, Array<OneD, NekDouble>>(nvel * nvel);
252 m_invMetricTensor = Array<OneD, Array<OneD, NekDouble>>(nvel * nvel);
253 m_deriv = Array<OneD, Array<OneD, NekDouble>>(nvel * nvel);
254 m_invDeriv = Array<OneD, Array<OneD, NekDouble>>(nvel * nvel);
255 for (int i = 0; i < m_metricTensor.size(); i++)
256 {
257 m_metricTensor[i] = Array<OneD, NekDouble>(physTot, 0.0);
258 m_invMetricTensor[i] = Array<OneD, NekDouble>(physTot, 0.0);
259 m_deriv[i] = Array<OneD, NekDouble>(physTot, 0.0);
260 m_invDeriv[i] = Array<OneD, NekDouble>(physTot, 0.0);
261 }
262 m_jac = Array<OneD, NekDouble>(physTot, 0.0);
263
264 // First, calculate derivatives of the mapping -> dX^i/dx^j = c^i_j
265 for (int i = 0; i < nvel; i++)
266 {
267 for (int j = 0; j < nvel; j++)
268 {
270 m_coords[i], m_deriv[i * nvel + j]);
271 }
272 }
273 // In Homogeneous case, m_deriv(2,2) needs to be set to 1
274 // because differentiation in wavespace is not valid for non-periodic field
275 if (m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
276 {
277 Vmath::Fill(physTot, 1.0, m_deriv[2 * nvel + 2], 1);
278 }
279
280 // Now calculate the metric tensor --> g_ij = sum_k { c^k_i c^k_j }
281 for (int i = 0; i < nvel; i++)
282 {
283 for (int j = 0; j < nvel; j++)
284 {
285 for (int k = 0; k < nvel; k++)
286 {
287 Vmath::Vvtvp(physTot, m_deriv[k * nvel + i], 1,
288 m_deriv[k * nvel + j], 1,
289 m_metricTensor[i * nvel + j], 1,
290 m_metricTensor[i * nvel + j], 1);
291 }
292 }
293 }
294
295 // Put the adjoint of g in m_invMetricTensor
296 switch (nvel)
297 {
298 case 1:
299 Vmath::Fill(physTot, 1.0, m_invMetricTensor[0], 1);
300 break;
301 case 2:
302 Vmath::Vcopy(physTot, m_metricTensor[1 * nvel + 1], 1,
303 m_invMetricTensor[0 * nvel + 0], 1);
304 Vmath::Smul(physTot, -1.0, m_metricTensor[0 * nvel + 1], 1,
305 m_invMetricTensor[1 * nvel + 0], 1);
306 Vmath::Smul(physTot, -1.0, m_metricTensor[1 * nvel + 0], 1,
307 m_invMetricTensor[0 * nvel + 1], 1);
308 Vmath::Vcopy(physTot, m_metricTensor[0 * nvel + 0], 1,
309 m_invMetricTensor[1 * nvel + 1], 1);
310 break;
311 case 3:
312 {
313 int a, b, c, d, e, i, j;
314
315 // Compute g^{ij} by computing Cofactors(g_ij)^T
316 for (i = 0; i < nvel; ++i)
317 {
318 for (j = 0; j < nvel; ++j)
319 {
320 a = ((i + 1) % nvel) * nvel + ((j + 1) % nvel);
321 b = ((i + 1) % nvel) * nvel + ((j + 2) % nvel);
322 c = ((i + 2) % nvel) * nvel + ((j + 1) % nvel);
323 d = ((i + 2) % nvel) * nvel + ((j + 2) % nvel);
324 e = i * nvel + j;
325 // a*d - b*c
326 Vmath::Vmul(physTot, m_metricTensor[b], 1,
327 m_metricTensor[c], 1, m_invMetricTensor[e], 1);
328 Vmath::Vvtvm(physTot, m_metricTensor[a], 1,
330 m_invMetricTensor[e], 1);
331 }
332 }
333 break;
334 }
335 }
336
337 // Compute g = det(g_{ij}) (= Jacobian squared) and store
338 // temporarily in m_jac.
339 for (int i = 0; i < nvel; ++i)
340 {
341 Vmath::Vvtvp(physTot, m_metricTensor[i], 1, m_invMetricTensor[i * nvel],
342 1, m_jac, 1, m_jac, 1);
343 }
344
345 // Calculate g^ij (the inverse of g_ij) by dividing by jac
346 for (int i = 0; i < nvel * nvel; ++i)
347 {
348 Vmath::Vdiv(physTot, m_invMetricTensor[i], 1, m_jac, 1,
349 m_invMetricTensor[i], 1);
350 }
351
352 // Compute the Jacobian = sqrt(g)
353 Vmath::Vsqrt(physTot, m_jac, 1, m_jac, 1);
354
355 // Calculate the derivatives of the inverse transformation
356 // c'^j_i = dx^j/dX^i = sum_k {g^jk c^i_k}
357 for (int i = 0; i < nvel; ++i)
358 {
359 for (int j = 0; j < nvel; ++j)
360 {
361 for (int k = 0; k < nvel; ++k)
362 {
363 Vmath::Vvtvp(physTot, m_deriv[i * nvel + k], 1,
364 m_invMetricTensor[j * nvel + k], 1,
365 m_invDeriv[i * nvel + j], 1,
366 m_invDeriv[i * nvel + j], 1);
367 }
368 }
369 }
370
371 // Restore value of wavespace
372 m_fields[0]->SetWaveSpace(wavespace);
373}
Array< OneD, Array< OneD, NekDouble > > m_invDeriv
Array< OneD, Array< OneD, NekDouble > > m_deriv
Array< OneD, NekDouble > m_jac
Array< OneD, Array< OneD, NekDouble > > m_coords
Array with the Cartesian coordinates.
Definition: Mapping.h:412
std::vector< double > d(NPUPPER *NPUPPER)
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
Definition: Vmath.hpp:340
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.hpp:72
void Vvtvm(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)
vvtvm (vector times vector minus vector): z = w*x - y
Definition: Vmath.hpp:381
void Vdiv(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x/y.
Definition: Vmath.hpp:126
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.hpp:54
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

References Nektar::UnitTests::d(), Nektar::MultiRegions::DirCartesianMap, Nektar::MultiRegions::e3DH1D, Vmath::Fill(), Nektar::GlobalMapping::Mapping::m_coords, m_deriv, Nektar::GlobalMapping::Mapping::m_fields, m_invDeriv, m_invMetricTensor, m_jac, m_metricTensor, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, Vmath::Smul(), Vmath::Vcopy(), Vmath::Vdiv(), Vmath::Vmul(), Vmath::Vsqrt(), Vmath::Vvtvm(), and Vmath::Vvtvp().

Referenced by v_UpdateGeomInfo().

◆ create()

static GLOBAL_MAPPING_EXPORT MappingSharedPtr Nektar::GlobalMapping::MappingGeneral::create ( const LibUtilities::SessionReaderSharedPtr pSession,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const TiXmlElement *  pMapping 
)
inlinestatic

Creates an instance of this class.

Definition at line 50 of file MappingGeneral.h.

54 {
57 p->InitObject(pFields, pMapping);
58 return p;
59 }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
GLOBAL_MAPPING_EXPORT typedef std::shared_ptr< Mapping > MappingSharedPtr
A shared pointer to a Mapping object.
Definition: Mapping.h:57

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::GlobalMapping::MappingSharedPtr, and CellMLToNektar.cellml_metadata::p.

◆ v_ApplyChristoffelContravar()

void Nektar::GlobalMapping::MappingGeneral::v_ApplyChristoffelContravar ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
overrideprotectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 187 of file MappingGeneral.cpp.

190{
191 int physTot = m_fields[0]->GetTotPoints();
192 int nvel = m_nConvectiveFields;
193
194 // Calculate {i,jk} u^j
195 for (int i = 0; i < nvel; i++)
196 {
197 for (int k = 0; k < nvel; k++)
198 {
199 outarray[i * nvel + k] = Array<OneD, NekDouble>(physTot, 0.0);
200 for (int j = 0; j < nvel; j++)
201 {
202 Vmath::Vvtvp(physTot, inarray[j], 1,
203 m_Christoffel[i * nvel * nvel + j * nvel + k], 1,
204 outarray[i * nvel + k], 1, outarray[i * nvel + k],
205 1);
206 }
207 }
208 }
209}

References m_Christoffel, Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, and Vmath::Vvtvp().

◆ v_ApplyChristoffelCovar()

void Nektar::GlobalMapping::MappingGeneral::v_ApplyChristoffelCovar ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
overrideprotectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 211 of file MappingGeneral.cpp.

214{
215 int physTot = m_fields[0]->GetTotPoints();
216 int nvel = m_nConvectiveFields;
217
218 // Calculate {i,jk} u_i
219 for (int j = 0; j < nvel; j++)
220 {
221 for (int k = 0; k < nvel; k++)
222 {
223 outarray[j * nvel + k] = Array<OneD, NekDouble>(physTot, 0.0);
224 for (int i = 0; i < nvel; i++)
225 {
226 Vmath::Vvtvp(physTot, inarray[i], 1,
227 m_Christoffel[i * nvel * nvel + j * nvel + k], 1,
228 outarray[j * nvel + k], 1, outarray[j * nvel + k],
229 1);
230 }
231 }
232 }
233}

References m_Christoffel, Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, and Vmath::Vvtvp().

◆ v_ContravarFromCartesian()

void Nektar::GlobalMapping::MappingGeneral::v_ContravarFromCartesian ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
overrideprotectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 111 of file MappingGeneral.cpp.

114{
115 int physTot = m_fields[0]->GetTotPoints();
116 int nvel = m_nConvectiveFields;
117
118 for (int i = 0; i < nvel; i++)
119 {
120 outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
121 for (int j = 0; j < nvel; j++)
122 {
123 Vmath::Vvtvp(physTot, inarray[j], 1, m_invDeriv[j * nvel + i], 1,
124 outarray[i], 1, outarray[i], 1);
125 }
126 }
127}

References Nektar::GlobalMapping::Mapping::m_fields, m_invDeriv, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, and Vmath::Vvtvp().

◆ v_ContravarToCartesian()

void Nektar::GlobalMapping::MappingGeneral::v_ContravarToCartesian ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
overrideprotectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 75 of file MappingGeneral.cpp.

78{
79 int physTot = m_fields[0]->GetTotPoints();
80 int nvel = m_nConvectiveFields;
81
82 for (int i = 0; i < nvel; i++)
83 {
84 outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
85 for (int j = 0; j < nvel; j++)
86 {
87 Vmath::Vvtvp(physTot, inarray[j], 1, m_deriv[i * nvel + j], 1,
88 outarray[i], 1, outarray[i], 1);
89 }
90 }
91}

References m_deriv, Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, and Vmath::Vvtvp().

◆ v_CovarFromCartesian()

void Nektar::GlobalMapping::MappingGeneral::v_CovarFromCartesian ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
overrideprotectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 129 of file MappingGeneral.cpp.

132{
133 int physTot = m_fields[0]->GetTotPoints();
134 int nvel = m_nConvectiveFields;
135
136 for (int i = 0; i < nvel; i++)
137 {
138 outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
139 for (int j = 0; j < nvel; j++)
140 {
141 Vmath::Vvtvp(physTot, inarray[j], 1, m_deriv[j * nvel + i], 1,
142 outarray[i], 1, outarray[i], 1);
143 }
144 }
145}

References m_deriv, Nektar::GlobalMapping::Mapping::m_fields, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, and Vmath::Vvtvp().

◆ v_CovarToCartesian()

void Nektar::GlobalMapping::MappingGeneral::v_CovarToCartesian ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
overrideprotectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 93 of file MappingGeneral.cpp.

96{
97 int physTot = m_fields[0]->GetTotPoints();
98 int nvel = m_nConvectiveFields;
99
100 for (int i = 0; i < nvel; i++)
101 {
102 outarray[i] = Array<OneD, NekDouble>(physTot, 0.0);
103 for (int j = 0; j < nvel; j++)
104 {
105 Vmath::Vvtvp(physTot, inarray[j], 1, m_invDeriv[i * nvel + j], 1,
106 outarray[i], 1, outarray[i], 1);
107 }
108 }
109}

References Nektar::GlobalMapping::Mapping::m_fields, m_invDeriv, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, and Vmath::Vvtvp().

◆ v_GetInvMetricTensor()

void Nektar::GlobalMapping::MappingGeneral::v_GetInvMetricTensor ( Array< OneD, Array< OneD, NekDouble > > &  outarray)
overrideprotectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 170 of file MappingGeneral.cpp.

172{
173 int physTot = m_fields[0]->GetTotPoints();
174 int nvel = m_nConvectiveFields;
175
176 for (int i = 0; i < nvel; i++)
177 {
178 for (int j = 0; j < nvel; j++)
179 {
180 outarray[i * nvel + j] = Array<OneD, NekDouble>(physTot, 0.0);
181 Vmath::Vcopy(physTot, m_invMetricTensor[i * nvel + j], 1,
182 outarray[i * nvel + j], 1);
183 }
184 }
185}

References Nektar::GlobalMapping::Mapping::m_fields, m_invMetricTensor, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, and Vmath::Vcopy().

◆ v_GetJacobian()

void Nektar::GlobalMapping::MappingGeneral::v_GetJacobian ( Array< OneD, NekDouble > &  outarray)
overrideprotectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 147 of file MappingGeneral.cpp.

148{
149 int physTot = m_fields[0]->GetTotPoints();
150 Vmath::Vcopy(physTot, m_jac, 1, outarray, 1);
151}

References Nektar::GlobalMapping::Mapping::m_fields, m_jac, and Vmath::Vcopy().

◆ v_GetMetricTensor()

void Nektar::GlobalMapping::MappingGeneral::v_GetMetricTensor ( Array< OneD, Array< OneD, NekDouble > > &  outarray)
overrideprotectedvirtual

Implements Nektar::GlobalMapping::Mapping.

Definition at line 153 of file MappingGeneral.cpp.

155{
156 int physTot = m_fields[0]->GetTotPoints();
157 int nvel = m_nConvectiveFields;
158
159 for (int i = 0; i < nvel; i++)
160 {
161 for (int j = 0; j < nvel; j++)
162 {
163 outarray[i * nvel + j] = Array<OneD, NekDouble>(physTot, 0.0);
164 Vmath::Vcopy(physTot, m_metricTensor[i * nvel + j], 1,
165 outarray[i * nvel + j], 1);
166 }
167 }
168}

References Nektar::GlobalMapping::Mapping::m_fields, m_metricTensor, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, and Vmath::Vcopy().

◆ v_InitObject()

void Nektar::GlobalMapping::MappingGeneral::v_InitObject ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const TiXmlElement *  pMapping 
)
overrideprotectedvirtual

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 from Nektar::GlobalMapping::Mapping.

Definition at line 63 of file MappingGeneral.cpp.

66{
67 Mapping::v_InitObject(pFields, pMapping);
68
69 m_constantJacobian = false;
70
72 "General Mapping needs at least 2 velocity components.");
73}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
virtual GLOBAL_MAPPING_EXPORT void v_InitObject(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
Definition: Mapping.cpp:95
bool m_constantJacobian
Flag defining if the Jacobian is constant.
Definition: Mapping.h:427

References ASSERTL0, Nektar::GlobalMapping::Mapping::m_constantJacobian, Nektar::GlobalMapping::Mapping::m_nConvectiveFields, and Nektar::GlobalMapping::Mapping::v_InitObject().

◆ v_UpdateGeomInfo()

void Nektar::GlobalMapping::MappingGeneral::v_UpdateGeomInfo ( )
overrideprotectedvirtual

Friends And Related Function Documentation

◆ MemoryManager< MappingGeneral >

friend class MemoryManager< MappingGeneral >
friend

Definition at line 1 of file MappingGeneral.h.

Member Data Documentation

◆ className

std::string Nektar::GlobalMapping::MappingGeneral::className
static
Initial value:
=
"X = X(x,y,z), Y = Y(x,y,z), Z=Z(x,y,z)")
static GLOBAL_MAPPING_EXPORT MappingSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pMapping)
Creates an instance of this class.
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
MappingFactory & GetMappingFactory()
Declaration of the mapping factory singleton.
Definition: Mapping.cpp:47

Name of the class.

Definition at line 62 of file MappingGeneral.h.

◆ m_Christoffel

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::MappingGeneral::m_Christoffel
protected

◆ m_deriv

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::MappingGeneral::m_deriv
protected

◆ m_invDeriv

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::MappingGeneral::m_invDeriv
protected

◆ m_invMetricTensor

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::MappingGeneral::m_invMetricTensor
protected

◆ m_jac

Array<OneD, NekDouble> Nektar::GlobalMapping::MappingGeneral::m_jac
protected

Definition at line 70 of file MappingGeneral.h.

Referenced by CalculateMetricTerms(), and v_GetJacobian().

◆ m_metricTensor

Array<OneD, Array<OneD, NekDouble> > Nektar::GlobalMapping::MappingGeneral::m_metricTensor
protected

Definition at line 65 of file MappingGeneral.h.

Referenced by CalculateChristoffel(), CalculateMetricTerms(), and v_GetMetricTensor().