Intrepid2
Intrepid2_OrientationTools.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Intrepid2 Package
4 //
5 // Copyright 2007 NTESS and the Intrepid2 contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
15 #ifndef __INTREPID2_ORIENTATIONTOOLS_HPP__
16 #define __INTREPID2_ORIENTATIONTOOLS_HPP__
17 
18 #include "Intrepid2_ConfigDefs.hpp"
19 #include "Intrepid2_Types.hpp"
20 #include "Intrepid2_Utils.hpp"
21 
22 #include "Shards_CellTopology.hpp"
23 #include "Shards_BasicTopologies.hpp"
24 
25 #include "Intrepid2_PointTools.hpp"
26 
27 #include "Intrepid2_Basis.hpp"
28 
29 // -- HGRAD family
33 
36 
37 // -- HCURL family
40 
44 
45 // -- HDIV family
48 
52 
53 // -- Lower order family
56 
59 
63 
67 
70 
71 #include "Teuchos_LAPACK.hpp"
72 
73 
74 namespace Intrepid2 {
75 
76  namespace Impl {
77 
82  public:
83 
84  // -----------------------------------------------------------------------------
85  // Point modification
86  //
87  //
88 
95  template<typename ValueType>
96  KOKKOS_INLINE_FUNCTION
97  static void
98  getModifiedLinePoint(ValueType &ot,
99  const ValueType pt,
100  const ordinal_type ort);
101 
110  template<typename ValueType>
111  KOKKOS_INLINE_FUNCTION
112  static void
113  getModifiedTrianglePoint(ValueType &ot0,
114  ValueType &ot1,
115  const ValueType pt0,
116  const ValueType pt1,
117  const ordinal_type ort);
118 
127  template<typename ValueType>
128  KOKKOS_INLINE_FUNCTION
129  static void
130  getModifiedQuadrilateralPoint(ValueType &ot0,
131  ValueType &ot1,
132  const ValueType pt0,
133  const ValueType pt1,
134  const ordinal_type ort);
135 
143  template<typename outPointViewType,
144  typename refPointViewType>
145  inline
146  static void
147  mapToModifiedReference(outPointViewType outPoints,
148  const refPointViewType refPoints,
149  const shards::CellTopology cellTopo,
150  const ordinal_type cellOrt = 0);
151 
159  template<typename outPointViewType,
160  typename refPointViewType>
161  KOKKOS_INLINE_FUNCTION
162  static void
163  mapToModifiedReference(outPointViewType outPoints,
164  const refPointViewType refPoints,
165  const unsigned cellTopoKey,
166  const ordinal_type cellOrt = 0);
167 
168 
174  template<typename JacobianViewType>
175  KOKKOS_INLINE_FUNCTION
176  static void
177  getLineJacobian(JacobianViewType jacobian, const ordinal_type ort);
178 
184  template<typename JacobianViewType>
185  KOKKOS_INLINE_FUNCTION
186  static void
187  getTriangleJacobian(JacobianViewType jacobian, const ordinal_type ort);
188 
194  template<typename JacobianViewType>
195  KOKKOS_INLINE_FUNCTION
196  static void
197  getQuadrilateralJacobian(JacobianViewType jacobian, const ordinal_type ort);
198 
199 
206  template<typename JacobianViewType>
207  inline
208  static void
209  getJacobianOfOrientationMap(JacobianViewType jacobian,
210  const shards::CellTopology cellTopo,
211  const ordinal_type cellOrt);
212 
219  template<typename JacobianViewType>
220  KOKKOS_INLINE_FUNCTION
221  static void
222  getJacobianOfOrientationMap(JacobianViewType jacobian,
223  const unsigned cellTopoKey,
224  const ordinal_type cellOrt);
225 
226 
234  template<typename TanViewType, typename ParamViewType>
235  KOKKOS_INLINE_FUNCTION
236  static void getRefSubcellTangents(TanViewType tangents,
237  const ParamViewType subCellParametrization,
238  const unsigned subcellTopoKey,
239  const ordinal_type subCellOrd,
240  const ordinal_type ort);
241 
242 
253  template<typename TanNormViewType, typename ParamViewType>
254  KOKKOS_INLINE_FUNCTION
255  static void getRefSideTangentsAndNormal(TanNormViewType tangentsAndNormal,
256  const ParamViewType subCellParametrization,
257  const unsigned subcellTopoKey,
258  const ordinal_type subCellOrd,
259  const ordinal_type ort);
260 
261 
271  template<typename coordsViewType, typename subcellCoordsViewType, typename ParamViewType>
272  KOKKOS_INLINE_FUNCTION
273  static void mapSubcellCoordsToRefCell(coordsViewType cellCoords,
274  const subcellCoordsViewType subCellCoords,
275  const ParamViewType subcellParametrization,
276  const unsigned subcellTopoKey,
277  const ordinal_type subCellOrd,
278  const ordinal_type ort);
279 
280  // -----------------------------------------------------------------------------
281  // Coefficient Matrix
282  //
283  //
284 
296  template<typename OutputViewType,
297  typename subcellBasisHostType,
298  typename cellBasisHostType>
299  inline
300  static void
301  getCoeffMatrix_HGRAD(OutputViewType &output,
302  const subcellBasisHostType& subcellBasis,
303  const cellBasisHostType& cellBasis,
304  const ordinal_type subcellId,
305  const ordinal_type subcellOrt,
306  const bool inverse = false);
307 
318  template<typename OutputViewType,
319  typename subcellBasisHostType,
320  typename cellBasisHostType>
321  inline
322  static void
323  getCoeffMatrix_HCURL(OutputViewType &output,
324  const subcellBasisHostType& subcellBasis,
325  const cellBasisHostType& cellBasis,
326  const ordinal_type subcellId,
327  const ordinal_type subcellOrt,
328  const bool inverse = false);
329 
330 
341  template<typename OutputViewType,
342  typename subcellBasisHostType,
343  typename cellBasisHostType>
344  inline
345  static void
346  getCoeffMatrix_HDIV(OutputViewType &output,
347  const subcellBasisHostType& subcellBasis,
348  const cellBasisHostType& cellBasis,
349  const ordinal_type subcellId,
350  const ordinal_type subcellOrt,
351  const bool inverse = false);
352 
353 
362  template<typename OutputViewType,
363  typename cellBasisHostType>
364  inline
365  static void
366  getCoeffMatrix_HVOL(OutputViewType &output,
367  const cellBasisHostType& cellBasis,
368  const ordinal_type cellOrt,
369  const bool inverse = false);
370 
371  };
372  }
373 
377  template<typename DeviceType>
379  public:
380 
383  using CoeffMatrixDataViewType = Kokkos::View<double****,DeviceType>;
384 
385  //
388  using OrtCoeffDataType = std::map<std::pair<std::string,ordinal_type>,Kokkos::View<double****,DeviceType> >;
389 
390  static OrtCoeffDataType ortCoeffData;
391  static OrtCoeffDataType ortInvCoeffData;
392 
393  private:
394 
395  template<typename BasisHostType>
396  inline
397  static CoeffMatrixDataViewType createCoeffMatrixInternal(const BasisHostType* basis, const bool invTrans = false);
398 
399 
402  template<typename BasisHostType>
403  inline
404  static void init_HGRAD(CoeffMatrixDataViewType matData,
405  BasisHostType const *cellBasis,
406  const bool inverse = false);
407 
410  template<typename BasisHostType>
411  inline
412  static void init_HCURL(CoeffMatrixDataViewType matData,
413  BasisHostType const *cellBasis,
414  const bool inverse = false);
415 
418  template<typename BasisHostType>
419  inline
420  static void init_HDIV(CoeffMatrixDataViewType matData,
421  BasisHostType const *cellBasis,
422  const bool inverse = false);
423 
426  template<typename BasisHostType>
427  inline
428  static void init_HVOL(CoeffMatrixDataViewType matData,
429  BasisHostType const *cellBasis,
430  const bool inverse = false);
431 
432  public:
433 
437  template<typename BasisType>
438  inline
439  static CoeffMatrixDataViewType createCoeffMatrix(const BasisType* basis);
440 
444  template<typename BasisType>
445  inline
446  static CoeffMatrixDataViewType createInvCoeffMatrix(const BasisType* basis);
447 
450  inline
451  static void clearCoeffMatrix();
452 
459  template<typename elemOrtValueType, class ...elemOrtProperties,
460  typename elemNodeValueType, class ...elemNodeProperties>
461  inline
462  static void
463  getOrientation(Kokkos::DynRankView<elemOrtValueType,elemOrtProperties...> elemOrts,
464  const Kokkos::DynRankView<elemNodeValueType,elemNodeProperties...> elemNodes,
465  const shards::CellTopology cellTopo,
466  bool isSide = false);
467 
468 
476  template<typename outputValueType, class ...outputProperties,
477  typename inputValueType, class ...inputProperties,
478  typename OrientationViewType,
479  typename BasisType>
480  inline
481  static void
482  modifyBasisByOrientation(Kokkos::DynRankView<outputValueType,outputProperties...> output,
483  const Kokkos::DynRankView<inputValueType, inputProperties...> input,
484  const OrientationViewType orts,
485  const BasisType * basis,
486  const bool transpose = false);
487 
494  template<typename outputValueType, class ...outputProperties,
495  typename inputValueType, class ...inputProperties,
496  typename OrientationViewType,
497  typename BasisType>
498  inline
499  static void
500  modifyBasisByOrientationTranspose(Kokkos::DynRankView<outputValueType,outputProperties...> output,
501  const Kokkos::DynRankView<inputValueType, inputProperties...> input,
502  const OrientationViewType orts,
503  const BasisType * basis);
504 
505 
513  template<typename outputValueType, class ...outputProperties,
514  typename inputValueType, class ...inputProperties,
515  typename OrientationViewType,
516  typename BasisType>
517  inline
518  static void
519  modifyBasisByOrientationInverse(Kokkos::DynRankView<outputValueType,outputProperties...> output,
520  const Kokkos::DynRankView<inputValueType, inputProperties...> input,
521  const OrientationViewType orts,
522  const BasisType * basis,
523  const bool transpose = false);
524 
525 
533  template<typename outputValueType, class ...outputProperties,
534  typename inputValueType, class ...inputProperties,
535  typename OrientationViewType,
536  typename BasisTypeLeft,
537  typename BasisTypeRight>
538  inline
539  static void
540  modifyMatrixByOrientation(Kokkos::DynRankView<outputValueType,outputProperties...> output,
541  const Kokkos::DynRankView<inputValueType, inputProperties...> input,
542  const OrientationViewType orts,
543  const BasisTypeLeft* basisLeft,
544  const BasisTypeRight* basisRight);
545  };
546 
547  //Definition of static members
548  template<typename T>
550 
551  template<typename T>
553 }
554 
555 // include templated function definitions
558 #include "Intrepid2_OrientationToolsDefCoeffMatrix_HCURL.hpp"
563 
564 #endif
static KOKKOS_INLINE_FUNCTION void getModifiedLinePoint(ValueType &ot, const ValueType pt, const ordinal_type ort)
Computes modified point for line segment.
Stateless classes that act as factories for two families of hierarchical bases. HierarchicalBasisFami...
Header file for the Intrepid2::Basis_HGRAD_LINE_Cn_FEM class.
static void getCoeffMatrix_HGRAD(OutputViewType &output, const subcellBasisHostType &subcellBasis, const cellBasisHostType &cellBasis, const ordinal_type subcellId, const ordinal_type subcellOrt, const bool inverse=false)
Compute orientation matrix for HGRAD basis for a given subcell and its reference basis.
Header file for the Intrepid2::Basis_HDIV_TET_I1_FEM class.
Header file for the Intrepid2::Basis_HCURL_HEX_In_FEM class.
static KOKKOS_INLINE_FUNCTION void getModifiedQuadrilateralPoint(ValueType &ot0, ValueType &ot1, const ValueType pt0, const ValueType pt1, const ordinal_type ort)
Computes modified point for quadrilateral.
Header file for the Intrepid2::Basis_HDIV_HEX_I1_FEM class.
Header file for the Intrepid2::Basis_HDIV_WEDGE_I1_FEM class.
static KOKKOS_INLINE_FUNCTION void mapSubcellCoordsToRefCell(coordsViewType cellCoords, const subcellCoordsViewType subCellCoords, const ParamViewType subcellParametrization, const unsigned subcellTopoKey, const ordinal_type subCellOrd, const ordinal_type ort)
Maps points defined on the subCell manifold into the parent Cell accounting for orientation.
Definition file for the Intrepid2::OrientationTools class.
static void getJacobianOfOrientationMap(JacobianViewType jacobian, const shards::CellTopology cellTopo, const ordinal_type cellOrt)
Computes jacobian of the parameterization maps of 1- and 2-subcells with orientation.
Header file for the Intrepid2::Basis_HDIV_HEX_In_FEM class.
static void getCoeffMatrix_HCURL(OutputViewType &output, const subcellBasisHostType &subcellBasis, const cellBasisHostType &cellBasis, const ordinal_type subcellId, const ordinal_type subcellOrt, const bool inverse=false)
Compute orientation matrix for HCURL basis for a given subcell and its reference basis.
Creation of orientation matrix A of a face or edge for HDIV elements.
static void clearCoeffMatrix()
Clear coefficient matrix.
Header function for Intrepid2::Util class and other utility functions.
static KOKKOS_INLINE_FUNCTION void getRefSubcellTangents(TanViewType tangents, const ParamViewType subCellParametrization, const unsigned subcellTopoKey, const ordinal_type subCellOrd, const ordinal_type ort)
Computes the (oriented) subCell tangents.
static void modifyMatrixByOrientation(Kokkos::DynRankView< outputValueType, outputProperties...> output, const Kokkos::DynRankView< inputValueType, inputProperties...> input, const OrientationViewType orts, const BasisTypeLeft *basisLeft, const BasisTypeRight *basisRight)
Modify an assembled (C,F1,F2) matrix according to orientation of the cells.
Header file for the Intrepid2::Basis_HGRAD_TRI_Cn_FEM class.
Header file for the Intrepid2::Basis_HVOL_LINE_Cn_FEM class.
Header file for the Intrepid2::Basis_HDIV_QUAD_I1_FEM class.
Header file for the Intrepid2::Basis_HCURL_TET_In_FEM class.
Header file for the Intrepid2::Basis_HGRAD_TET_Cn_FEM class.
static void init_HDIV(CoeffMatrixDataViewType matData, BasisHostType const *cellBasis, const bool inverse=false)
Compute orientation matrix for HDIV basis.
Kokkos::View< double ****, DeviceType > CoeffMatrixDataViewType
subcell ordinal, orientation, matrix m x n
Header file for the Intrepid2::Basis_HDIV_TRI_I1_FEM class.
static void modifyBasisByOrientation(Kokkos::DynRankView< outputValueType, outputProperties...> output, const Kokkos::DynRankView< inputValueType, inputProperties...> input, const OrientationViewType orts, const BasisType *basis, const bool transpose=false)
Modify basis due to orientation.
Header file for the Intrepid2::Basis_HDIV_QUAD_In_FEM class.
Header file for the Intrepid2::Basis_HCURL_QUAD_In_FEM class.
static KOKKOS_INLINE_FUNCTION void getModifiedTrianglePoint(ValueType &ot0, ValueType &ot1, const ValueType pt0, const ValueType pt1, const ordinal_type ort)
Computes modified point for triangle.
static KOKKOS_INLINE_FUNCTION void getTriangleJacobian(JacobianViewType jacobian, const ordinal_type ort)
Computes Jacobian of orientation map for triangle.
static KOKKOS_INLINE_FUNCTION void getLineJacobian(JacobianViewType jacobian, const ordinal_type ort)
Computes Jacobian of orientation map for line segment.
static void getCoeffMatrix_HVOL(OutputViewType &output, const cellBasisHostType &cellBasis, const ordinal_type cellOrt, const bool inverse=false)
Compute orientation matrix for HVOL basis for a given (2D or 1D) cell and its reference basis...
Header file for the Intrepid2::Basis_HDIV_TRI_In_FEM class.
Contains definitions of custom data types in Intrepid2.
Header file for the Intrepid2::Basis_HGRAD_QUAD_Cn_FEM class.
Creation of orientation matrix A of a face or edge for HGRAD elements.
Header file for the Intrepid2::Basis_HDIV_TET_In_FEM class.
Creation of orientation matrix A of a face or edge for HGRAD elements.
static void init_HGRAD(CoeffMatrixDataViewType matData, BasisHostType const *cellBasis, const bool inverse=false)
Compute orientation matrix for HGRAD basis.
Definition file for functions that modify points due to orientation in the Intrepid2::Impl::Orientati...
static void getCoeffMatrix_HDIV(OutputViewType &output, const subcellBasisHostType &subcellBasis, const cellBasisHostType &cellBasis, const ordinal_type subcellId, const ordinal_type subcellOrt, const bool inverse=false)
Compute orientation matrix for HDIV basis for a given subcell and its reference basis.
Header file for the Intrepid2::Basis_HCURL_QUAD_I1_FEM class.
Header file for the Intrepid2::Basis_HCURL_WEDGE_I1_FEM class.
static void getOrientation(Kokkos::DynRankView< elemOrtValueType, elemOrtProperties...> elemOrts, const Kokkos::DynRankView< elemNodeValueType, elemNodeProperties...> elemNodes, const shards::CellTopology cellTopo, bool isSide=false)
Compute orientations of cells in a workset.
static void modifyBasisByOrientationTranspose(Kokkos::DynRankView< outputValueType, outputProperties...> output, const Kokkos::DynRankView< inputValueType, inputProperties...> input, const OrientationViewType orts, const BasisType *basis)
Modify basis due to orientation, applying the transpose of the operator applied in modifyBasisByOrien...
Tools to compute orientations for degrees-of-freedom.
Header file for the Intrepid2::Basis_HCURL_TET_I1_FEM class.
static CoeffMatrixDataViewType createCoeffMatrix(const BasisType *basis)
Create coefficient matrix.
static KOKKOS_INLINE_FUNCTION void getQuadrilateralJacobian(JacobianViewType jacobian, const ordinal_type ort)
Computes Jacobian of orientation map for quadrilateral.
Header file for the Intrepid2::Basis_HCURL_HEX_I1_FEM class.
static void mapToModifiedReference(outPointViewType outPoints, const refPointViewType refPoints, const shards::CellTopology cellTopo, const ordinal_type cellOrt=0)
Computes modified parameterization maps of 1- and 2-subcells with orientation.
std::map< std::pair< std::string, ordinal_type >, Kokkos::View< double ****, DeviceType > > OrtCoeffDataType
key :: basis name, order, value :: matrix data view type
Definition file for matrix data in the Intrepid2::OrientationTools class.
static KOKKOS_INLINE_FUNCTION void getRefSideTangentsAndNormal(TanNormViewType tangentsAndNormal, const ParamViewType subCellParametrization, const unsigned subcellTopoKey, const ordinal_type subCellOrd, const ordinal_type ort)
Computes the (oriented) tangents and normal of the side subCell The normals are only defined for side...
Header file for the Intrepid2::Basis_HCURL_TRI_I1_FEM class.
static void init_HCURL(CoeffMatrixDataViewType matData, BasisHostType const *cellBasis, const bool inverse=false)
Compute orientation matrix for HCURL basis.
Header file for the Intrepid2::Basis_HVOL_TRI_Cn_FEM class.
static void modifyBasisByOrientationInverse(Kokkos::DynRankView< outputValueType, outputProperties...> output, const Kokkos::DynRankView< inputValueType, inputProperties...> input, const OrientationViewType orts, const BasisType *basis, const bool transpose=false)
Modify basis due to orientation, applying the inverse of the operator applied in modifyBasisByOrienta...
static void init_HVOL(CoeffMatrixDataViewType matData, BasisHostType const *cellBasis, const bool inverse=false)
Compute orientation matrix for HVOL basis.
static CoeffMatrixDataViewType createInvCoeffMatrix(const BasisType *basis)
Create inverse of coefficient matrix.
Stateless class that acts as a factory for a family of nodal bases (hypercube topologies only at this...
Header file for the abstract base class Intrepid2::Basis.
Header file for the Intrepid2::Basis_HCURL_TRI_In_FEM class.
Header file for Intrepid2::PointTools class to provide utilities for barycentric coordinates, equispaced lattices, and warp-blend point distrubtions.
Tools to compute orientations for degrees-of-freedom.
Header file for the Intrepid2::Basis_HGRAD_HEX_Cn_FEM class.