Intrepid2
Intrepid2_Types.hpp
Go to the documentation of this file.
1// @HEADER
2// ************************************************************************
3//
4// Intrepid2 Package
5// Copyright (2007) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Kyungjoo Kim (kyukim@sandia.gov), or
38// Mauro Perego (mperego@sandia.gov)
39//
40// ************************************************************************
41// @HEADER
42
48#ifndef __INTREPID2_TYPES_HPP__
49#define __INTREPID2_TYPES_HPP__
50
51#include <Kokkos_Core.hpp>
52#include <Kokkos_DynRankView.hpp>
53
54#include <stdexcept>
55
56namespace Intrepid2 {
57
58 // use ordinal_type and size_type everywhere (no index type)
59 typedef int ordinal_type;
60 typedef size_t size_type;
61
62 template<typename ValueType>
63 KOKKOS_FORCEINLINE_FUNCTION
64 ValueType epsilon() {
65 return 0;
66 }
67
68 template<>
69 KOKKOS_FORCEINLINE_FUNCTION
70 double epsilon<double>() {
71 typedef union {
72 long long i64;
73 double d64;
74 } dbl_64;
75
76 dbl_64 s;
77 s.d64 = 1;
78 s.i64++;
79 return (s.i64 < 0 ? 1 - s.d64 : s.d64 - 1);
80 }
81
82 template<>
83 KOKKOS_FORCEINLINE_FUNCTION
84 float epsilon<float>() {
85 typedef union {
86 int i32;
87 float f32;
88 } flt_32;
89
90 flt_32 s;
91 s.f32 = 1;
92 s.f32++;
93 return (s.i32 < 0 ? 1 - s.f32 : s.f32 - 1);
94 }
95
96 KOKKOS_FORCEINLINE_FUNCTION
97 double epsilon() {
98 return epsilon<double>();
99 }
100
101 KOKKOS_FORCEINLINE_FUNCTION
102 double tolerence() {
103 return 100.0*epsilon();
104 }
105
106 KOKKOS_FORCEINLINE_FUNCTION
107 double threshold() {
108 return 10.0*epsilon();
109 }
110
113 public:
114 // KK: do not chagne max num pts per basis eval bigger than 1.
115 // polylib point and order match needs to be first examined; now if it is set bigger than 1
116 // it creates silent error.
117 //
118 // MP: I tried setting max num pts per basis eval to 2 and everything seems working fine. Leaving it to 1 for now.
119
121 static constexpr ordinal_type MaxNumPtsPerBasisEval= 1;
123 static constexpr ordinal_type MaxOrder = 8;
125 static constexpr ordinal_type MaxIntegrationPoints = 1001;
127 static constexpr ordinal_type MaxCubatureDegreeEdge= 20;
129 static constexpr ordinal_type MaxCubatureDegreeTri = 20;
131 static constexpr ordinal_type MaxCubatureDegreeTet = 20;
133 static constexpr ordinal_type MaxCubatureDegreePyr = 11;
135 static constexpr ordinal_type MaxDimension = 3;
137 static constexpr ordinal_type MaxNewton = 15;
139 static constexpr ordinal_type MaxDerivative = 10;
141 static constexpr ordinal_type MaxTensorComponents = 7;
143 static constexpr ordinal_type MaxVectorComponents = 7;
144
145 // we do not want to use hard-wired epsilon, threshold and tolerence.
146 // static constexpr double Epsilon = 1.0e-16;
147 // static constexpr double Threshold = 1.0e-15;
148 // static constexpr double Tolerence = 1.0e-14;
149 };
150 // const double Parameters::Epsilon = epsilon<double>(); // Platform-dependent machine epsilon.
151 // const double Parameters::Threshold = 10.0*epsilon<double>(); // Tolerance for various cell inclusion tests
152 // const double Parameters::Tolerence = 100.0*epsilon<double>(); // General purpose tolerance in, e.g., internal Newton's method to invert ref to phys maps
153
154 // ===================================================================
155 // Enum classes
156 // - Enum, String (char*) helper, valid
157 // - can be used on device and inside of kernel for debugging purpose
158 // - let's decorate kokkos inline
159
163 enum EPolyType {
164 POLYTYPE_GAUSS=0,
165 POLYTYPE_GAUSS_RADAU_LEFT,
166 POLYTYPE_GAUSS_RADAU_RIGHT,
167 POLYTYPE_GAUSS_LOBATTO,
168 POLYTYPE_MAX
169 };
170
171 KOKKOS_INLINE_FUNCTION
172 const char* EPolyTypeToString(const EPolyType polytype) {
173 switch(polytype) {
174 case POLYTYPE_GAUSS: return "Gauss";
175 case POLYTYPE_GAUSS_RADAU_LEFT: return "GaussRadauLeft";
176 case POLYTYPE_GAUSS_RADAU_RIGHT: return "GaussRadauRight";
177 case POLYTYPE_GAUSS_LOBATTO: return "GaussRadauLobatto";
178 case POLYTYPE_MAX: return "Max PolyType";
179 }
180 return "INVALID EPolyType";
181 }
182
188 KOKKOS_FORCEINLINE_FUNCTION
189 bool isValidPolyType(const EPolyType polytype){
190 return( polytype == POLYTYPE_GAUSS ||
191 polytype == POLYTYPE_GAUSS_RADAU_LEFT ||
192 polytype == POLYTYPE_GAUSS_RADAU_RIGHT ||
193 polytype == POLYTYPE_GAUSS_LOBATTO );
194 }
195
196
200 enum ECoordinates{
201 COORDINATES_CARTESIAN=0,
202 COORDINATES_POLAR,
203 COORDINATES_CYLINDRICAL,
204 COORDINATES_SPHERICAL,
205 COORDINATES_MAX
206 };
207
208 KOKKOS_INLINE_FUNCTION
209 const char* ECoordinatesToString(const ECoordinates coords) {
210 switch(coords) {
211 case COORDINATES_CARTESIAN: return "Cartesian";
212 case COORDINATES_POLAR: return "Polar";
213 case COORDINATES_CYLINDRICAL: return "Cylindrical";
214 case COORDINATES_SPHERICAL: return "Spherical";
215 case COORDINATES_MAX: return "Max. Coordinates";
216 }
217 return "INVALID ECoordinates";
218 }
219
225 KOKKOS_FORCEINLINE_FUNCTION
226 bool isValidCoordinate(const ECoordinates coordinateType){
227 return( coordinateType == COORDINATES_CARTESIAN ||
228 coordinateType == COORDINATES_POLAR ||
229 coordinateType == COORDINATES_CYLINDRICAL ||
230 coordinateType == COORDINATES_SPHERICAL );
231 }
232
236 enum ENorm{
237 NORM_ONE = 0,
238 NORM_TWO,
239 NORM_INF,
240 NORM_FRO, // Frobenius matrix norm
241 NORM_MAX
242 };
243
244 KOKKOS_INLINE_FUNCTION
245 const char* ENormToString(const ENorm norm) {
246 switch(norm) {
247 case NORM_ONE: return "1-Norm";
248 case NORM_TWO: return "2-Norm";
249 case NORM_INF: return "Infinity Norm";
250 case NORM_FRO: return "Frobenius Norm";
251 case NORM_MAX: return "Max. Norm";
252 }
253 return "INVALID ENorm";
254 }
255
261 KOKKOS_FORCEINLINE_FUNCTION
262 bool isValidNorm(const ENorm normType){
263 return( normType == NORM_ONE ||
264 normType == NORM_TWO ||
265 normType == NORM_INF ||
266 normType == NORM_FRO ||
267 normType == NORM_MAX );
268 }
269
275 enum EOperator{
276 OPERATOR_VALUE = 0,
277 OPERATOR_GRAD, // 1
278 OPERATOR_CURL, // 2
279 OPERATOR_DIV, // 3
280 OPERATOR_D1, // 4
281 OPERATOR_D2, // 5
282 OPERATOR_D3, // 6
283 OPERATOR_D4, // 7
284 OPERATOR_D5, // 8
285 OPERATOR_D6, // 9
286 OPERATOR_D7, // 10
287 OPERATOR_D8, // 11
288 OPERATOR_D9, // 12
289 OPERATOR_D10, // 13
290 OPERATOR_Dn, // 14
291 OPERATOR_MAX = OPERATOR_Dn // 14
292 };
293
294 KOKKOS_INLINE_FUNCTION
295 const char* EOperatorToString(const EOperator op) {
296 switch(op) {
297 case OPERATOR_VALUE: return "Value";
298 case OPERATOR_GRAD: return "Grad";
299 case OPERATOR_CURL: return "Curl";
300 case OPERATOR_DIV: return "Div";
301 case OPERATOR_D1: return "D1";
302 case OPERATOR_D2: return "D2";
303 case OPERATOR_D3: return "D3";
304 case OPERATOR_D4: return "D4";
305 case OPERATOR_D5: return "D5";
306 case OPERATOR_D6: return "D6";
307 case OPERATOR_D7: return "D7";
308 case OPERATOR_D8: return "D8";
309 case OPERATOR_D9: return "D9";
310 case OPERATOR_D10: return "D10";
311 case OPERATOR_MAX: return "Dn Operator";
312 }
313 return "INVALID EOperator";
314 }
315
321 KOKKOS_FORCEINLINE_FUNCTION
322 bool isValidOperator(const EOperator operatorType){
323 return ( operatorType == OPERATOR_VALUE ||
324 operatorType == OPERATOR_GRAD ||
325 operatorType == OPERATOR_CURL ||
326 operatorType == OPERATOR_DIV ||
327 operatorType == OPERATOR_D1 ||
328 operatorType == OPERATOR_D2 ||
329 operatorType == OPERATOR_D3 ||
330 operatorType == OPERATOR_D4 ||
331 operatorType == OPERATOR_D5 ||
332 operatorType == OPERATOR_D6 ||
333 operatorType == OPERATOR_D7 ||
334 operatorType == OPERATOR_D8 ||
335 operatorType == OPERATOR_D9 ||
336 operatorType == OPERATOR_D10 );
337 }
338
339
343 enum EFunctionSpace {
344 FUNCTION_SPACE_HGRAD = 0,
345 FUNCTION_SPACE_HCURL = 1,
346 FUNCTION_SPACE_HDIV = 2,
347 FUNCTION_SPACE_HVOL = 3,
348 FUNCTION_SPACE_VECTOR_HGRAD = 4,
349 FUNCTION_SPACE_TENSOR_HGRAD = 5,
350 FUNCTION_SPACE_MAX
351 };
352
353 KOKKOS_INLINE_FUNCTION
354 const char* EFunctionSpaceToString(const EFunctionSpace space) {
355 switch(space) {
356 case FUNCTION_SPACE_HGRAD: return "H(grad)";
357 case FUNCTION_SPACE_HCURL: return "H(curl)";
358 case FUNCTION_SPACE_HDIV: return "H(div)";
359 case FUNCTION_SPACE_HVOL: return "H(vol)";
360 case FUNCTION_SPACE_VECTOR_HGRAD: return "Vector H(grad)";
361 case FUNCTION_SPACE_TENSOR_HGRAD: return "Tensor H(grad)";
362 case FUNCTION_SPACE_MAX: return "Max. Function space";
363 }
364 return "INVALID EFunctionSpace";
365 }
366
372 KOKKOS_FORCEINLINE_FUNCTION
373 bool isValidFunctionSpace(const EFunctionSpace spaceType){
374 return ( spaceType == FUNCTION_SPACE_HGRAD ||
375 spaceType == FUNCTION_SPACE_HCURL ||
376 spaceType == FUNCTION_SPACE_HDIV ||
377 spaceType == FUNCTION_SPACE_HVOL ||
378 spaceType == FUNCTION_SPACE_VECTOR_HGRAD ||
379 spaceType == FUNCTION_SPACE_TENSOR_HGRAD );
380 }
381
390 enum EDiscreteSpace {
391 DISCRETE_SPACE_COMPLETE = 0, // value = 0
392 DISCRETE_SPACE_INCOMPLETE, // value = 1
393 DISCRETE_SPACE_BROKEN, // value = 2
394 DISCRETE_SPACE_MAX // value = 3
395 };
396
397 KOKKOS_INLINE_FUNCTION
398 const char* EDiscreteSpaceToString(const EDiscreteSpace space) {
399 switch(space) {
400 case DISCRETE_SPACE_COMPLETE: return "Complete";
401 case DISCRETE_SPACE_INCOMPLETE: return "Incomplete";
402 case DISCRETE_SPACE_BROKEN: return "Broken";
403 case DISCRETE_SPACE_MAX: return "Max. Rec. Space";
404 }
405 return "INVALID EDiscreteSpace";
406 }
407
413 KOKKOS_FORCEINLINE_FUNCTION
414 bool isValidDiscreteSpace(const EDiscreteSpace spaceType){
415 return ( spaceType == DISCRETE_SPACE_COMPLETE ||
416 spaceType == DISCRETE_SPACE_INCOMPLETE ||
417 spaceType == DISCRETE_SPACE_BROKEN );
418 }
419
423 enum EPointType {
424 POINTTYPE_EQUISPACED = 0, // value = 0
425 POINTTYPE_WARPBLEND,
426 POINTTYPE_GAUSS,
427 POINTTYPE_DEFAULT,
428 };
429
430 KOKKOS_INLINE_FUNCTION
431 const char* EPointTypeToString(const EPointType pointType) {
432 switch (pointType) {
433 case POINTTYPE_EQUISPACED: return "Equispaced Points";
434 case POINTTYPE_WARPBLEND: return "WarpBlend Points";
435 case POINTTYPE_GAUSS: return "Gauss Points";
436 case POINTTYPE_DEFAULT: return "Default Points";
437 }
438 return "INVALID EPointType";
439 }
440
445 KOKKOS_FORCEINLINE_FUNCTION
446 bool isValidPointType(const EPointType pointType) {
447 return ( pointType == POINTTYPE_EQUISPACED ||
448 pointType == POINTTYPE_WARPBLEND ||
449 pointType == POINTTYPE_GAUSS );
450 }
451
455 enum EBasis {
456 BASIS_FEM_DEFAULT = 0, // value = 0
457 BASIS_FEM_HIERARCHICAL, // value = 1
458 BASIS_FEM_LAGRANGIAN, // value = 2
459 BASIS_FVD_DEFAULT, // value = 3
460 BASIS_FVD_COVOLUME, // value = 4
461 BASIS_FVD_MIMETIC, // value = 5
462 BASIS_MAX // value = 6
463 };
464
465 KOKKOS_INLINE_FUNCTION
466 const char* EBasisToString(const EBasis basis) {
467 switch(basis) {
468 case BASIS_FEM_DEFAULT: return "FEM Default";
469 case BASIS_FEM_HIERARCHICAL: return "FEM Hierarchical";
470 case BASIS_FEM_LAGRANGIAN: return "FEM FIAT";
471 case BASIS_FVD_DEFAULT: return "FVD Default";
472 case BASIS_FVD_COVOLUME: return "FVD Covolume";
473 case BASIS_FVD_MIMETIC: return "FVD Mimetic";
474 case BASIS_MAX: return "Max. Basis";
475 }
476 return "INVALID EBasis";
477 }
478
484 KOKKOS_FORCEINLINE_FUNCTION
485 bool isValidBasis(const EBasis basisType){
486 return ( basisType == BASIS_FEM_DEFAULT ||
487 basisType == BASIS_FEM_HIERARCHICAL ||
488 basisType == BASIS_FEM_LAGRANGIAN ||
489 basisType == BASIS_FVD_DEFAULT ||
490 basisType == BASIS_FVD_COVOLUME ||
491 basisType == BASIS_FVD_MIMETIC );
492 }
493
494 // /** \enum Intrepid2::ECompEngine
495 // \brief Specifies how operators and functionals are computed internally
496 // (COMP_MANUAL = native C++ implementation, COMP_BLAS = BLAS implementation, etc.).
497 // */
498 // enum ECompEngine {
499 // COMP_CPP = 0,
500 // COMP_BLAS,
501 // COMP_ENGINE_MAX
502 // };
503
504 // KOKKOS_INLINE_FUNCTION
505 // const char* ECompEngineToString(const ECompEngine cEngine) {
506 // switch(cEngine) {
507 // case COMP_CPP: return "Native C++";
508 // case COMP_BLAS: return "BLAS";
509 // case COMP_ENGINE_MAX: return "Max. Comp. Engine";
510 // default: return "INVALID ECompEngine";
511 // }
512 // return "Error";
513 // }
514
515
516 // /** \brief Verifies validity of a computational engine enum
517
518 // \param compEngType [in] - enum of the computational engine
519 // \return 1 if the argument is valid computational engine; 0 otherwise
520 // */
521 // KOKKOS_FORCEINLINE_FUNCTION
522 // bool isValidCompEngine(const ECompEngine compEngType){
523 // //at the moment COMP_BLAS is not a valid CompEngine.
524 // return (compEngType == COMP_CPP);
525 // }
526
527
528} //namespace Intrepid2
529
530#endif
KOKKOS_FORCEINLINE_FUNCTION bool isValidFunctionSpace(const EFunctionSpace spaceType)
Verifies validity of a function space enum.
KOKKOS_FORCEINLINE_FUNCTION bool isValidPointType(const EPointType pointType)
Verifies validity of a point type enum.
KOKKOS_FORCEINLINE_FUNCTION bool isValidPolyType(const EPolyType polytype)
Verifies validity of a PolyType enum.
KOKKOS_FORCEINLINE_FUNCTION bool isValidNorm(const ENorm normType)
Verifies validity of a Norm enum.
KOKKOS_FORCEINLINE_FUNCTION bool isValidBasis(const EBasis basisType)
Verifies validity of a basis enum.
KOKKOS_FORCEINLINE_FUNCTION bool isValidCoordinate(const ECoordinates coordinateType)
Verifies validity of a Coordinate enum.
KOKKOS_FORCEINLINE_FUNCTION bool isValidDiscreteSpace(const EDiscreteSpace spaceType)
Verifies validity of a discrete space enum.
KOKKOS_FORCEINLINE_FUNCTION bool isValidOperator(const EOperator operatorType)
Verifies validity of an operator enum.
static constexpr ordinal_type MaxCubatureDegreeTri
The maximum degree of the polynomial that can be integrated exactly by a direct triangle rule.
static constexpr ordinal_type MaxCubatureDegreePyr
The maximum degree of the polynomial that can be integrated exactly by a direct pyramid rule.
static constexpr ordinal_type MaxCubatureDegreeTet
The maximum degree of the polynomial that can be integrated exactly by a direct tetrahedron rule.
static constexpr ordinal_type MaxCubatureDegreeEdge
The maximum degree of the polynomial that can be integrated exactly by a direct edge rule.
static constexpr ordinal_type MaxOrder
The maximum reconstruction order.
static constexpr ordinal_type MaxIntegrationPoints
The maximum number of integration points for direct cubature rules.
static constexpr ordinal_type MaxNewton
Maximum number of Newton iterations used internally in methods such as computing the action of the in...
static constexpr ordinal_type MaxTensorComponents
Maximum number of tensor/Cartesian products that can be taken: this allows hypercube basis in 7D to b...
static constexpr ordinal_type MaxNumPtsPerBasisEval
The maximum number of points to eval in serial mode.
static constexpr ordinal_type MaxVectorComponents
Maximum number of components that a VectorData object will store – 66 corresponds to OPERATOR_D10 on ...
static constexpr ordinal_type MaxDimension
The maximum ambient space dimension.
static constexpr ordinal_type MaxDerivative
Maximum order of derivatives allowed in intrepid.