-
Notifications
You must be signed in to change notification settings - Fork 559
/
Intrepid2_CellTools.hpp
1637 lines (1376 loc) · 91.8 KB
/
Intrepid2_CellTools.hpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
// @HEADER
// *****************************************************************************
// Intrepid2 Package
//
// Copyright 2007 NTESS and the Intrepid2 contributors.
// SPDX-License-Identifier: BSD-3-Clause
// *****************************************************************************
// @HEADER
/** \file Intrepid2_CellTools.hpp
\brief Header file for the Intrepid2::CellTools class.
\author Created by P. Bochev and D. Ridzal.
Kokkorized by Kyungjoo Kim
*/
#ifndef __INTREPID2_CELLTOOLS_HPP__
#define __INTREPID2_CELLTOOLS_HPP__
#include "Intrepid2_ConfigDefs.hpp"
#include "Shards_CellTopology.hpp"
#include "Shards_BasicTopologies.hpp"
#include "Teuchos_RCP.hpp"
#include "Intrepid2_Types.hpp"
#include "Intrepid2_Utils.hpp"
#include "Intrepid2_RealSpaceTools.hpp"
#include "Intrepid2_Basis.hpp"
#include "Intrepid2_HGRAD_LINE_C1_FEM.hpp"
#include "Intrepid2_HGRAD_LINE_C2_FEM.hpp"
#include "Intrepid2_HGRAD_TRI_C1_FEM.hpp"
#include "Intrepid2_HGRAD_TRI_C2_FEM.hpp"
#include "Intrepid2_HGRAD_QUAD_C1_FEM.hpp"
#include "Intrepid2_HGRAD_QUAD_C2_FEM.hpp"
#include "Intrepid2_HGRAD_TET_C1_FEM.hpp"
#include "Intrepid2_HGRAD_TET_C2_FEM.hpp"
#include "Intrepid2_HGRAD_TET_COMP12_FEM.hpp"
#include "Intrepid2_HGRAD_HEX_C1_FEM.hpp"
#include "Intrepid2_HGRAD_HEX_C2_FEM.hpp"
#include "Intrepid2_HGRAD_WEDGE_C1_FEM.hpp"
#include "Intrepid2_HGRAD_WEDGE_C2_FEM.hpp"
#include "Intrepid2_HGRAD_PYR_C1_FEM.hpp"
#include "Intrepid2_HGRAD_PYR_I2_FEM.hpp"
#include "Intrepid2_Data.hpp"
#include "Intrepid2_CellData.hpp"
namespace Intrepid2 {
//============================================================================================//
// //
// CellTools //
// //
//============================================================================================//
/** \class Intrepid2::CellTools
\brief A stateless class for operations on cell data. Provides methods for:
\li computing Jacobians of reference-to-physical frame mappings, their inverses and determinants
\li application of the reference-to-physical frame mapping and its inverse
\li parametrizations of edges and faces of reference cells needed for edge and face integrals,
\li computation of edge and face tangents and face normals on both reference and physical frames
\li inclusion tests for point sets in reference and physical cells.
*/
template<typename DeviceType>
class CellTools {
using ExecSpaceType = typename DeviceType::execution_space;
using MemSpaceType = typename DeviceType::memory_space;
public:
/** \brief Checks if a cell topology has reference cell
\param cell [in] - cell topology
\return true if the cell topology has reference cell, false otherwise
*/
inline
static bool
hasReferenceCell( const shards::CellTopology cellTopo ) {
return RefSubcellParametrization<DeviceType>::isSupported(cellTopo.getKey());
}
private:
/** \brief Generates default HGrad basis based on cell topology
\param cellTopo [in] - cell topology
*/
template<typename outputValueType,
typename pointValueType>
static Teuchos::RCP<Basis<DeviceType,outputValueType,pointValueType> >
createHGradBasis( const shards::CellTopology cellTopo ) {
Teuchos::RCP<Basis<DeviceType,outputValueType,pointValueType> > r_val;
switch (cellTopo.getKey()) {
case shards::Line<2>::key: r_val = Teuchos::rcp(new Basis_HGRAD_LINE_C1_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Line<3>::key: r_val = Teuchos::rcp(new Basis_HGRAD_LINE_C2_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Triangle<3>::key: r_val = Teuchos::rcp(new Basis_HGRAD_TRI_C1_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Quadrilateral<4>::key: r_val = Teuchos::rcp(new Basis_HGRAD_QUAD_C1_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Tetrahedron<4>::key: r_val = Teuchos::rcp(new Basis_HGRAD_TET_C1_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Hexahedron<8>::key: r_val = Teuchos::rcp(new Basis_HGRAD_HEX_C1_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Wedge<6>::key: r_val = Teuchos::rcp(new Basis_HGRAD_WEDGE_C1_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Pyramid<5>::key: r_val = Teuchos::rcp(new Basis_HGRAD_PYR_C1_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Triangle<6>::key: r_val = Teuchos::rcp(new Basis_HGRAD_TRI_C2_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Quadrilateral<8>::key: r_val = Teuchos::rcp(new Basis_HGRAD_QUAD_I2_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Quadrilateral<9>::key: r_val = Teuchos::rcp(new Basis_HGRAD_QUAD_C2_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Tetrahedron<10>::key: r_val = Teuchos::rcp(new Basis_HGRAD_TET_C2_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Tetrahedron<11>::key: r_val = Teuchos::rcp(new Basis_HGRAD_TET_COMP12_FEM<DeviceType,outputValueType,pointValueType>()); break;
case shards::Hexahedron<20>::key: r_val = Teuchos::rcp(new Basis_HGRAD_HEX_I2_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Hexahedron<27>::key: r_val = Teuchos::rcp(new Basis_HGRAD_HEX_C2_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Wedge<15>::key: r_val = Teuchos::rcp(new Basis_HGRAD_WEDGE_I2_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Wedge<18>::key: r_val = Teuchos::rcp(new Basis_HGRAD_WEDGE_C2_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Pyramid<13>::key: r_val = Teuchos::rcp(new Basis_HGRAD_PYR_I2_FEM <DeviceType,outputValueType,pointValueType>()); break;
case shards::Beam<2>::key:
case shards::Beam<3>::key:
case shards::ShellLine<2>::key:
case shards::ShellLine<3>::key:
case shards::ShellTriangle<3>::key:
case shards::ShellTriangle<6>::key:
case shards::ShellQuadrilateral<4>::key:
case shards::ShellQuadrilateral<8>::key:
case shards::ShellQuadrilateral<9>::key:
default: {
INTREPID2_TEST_FOR_EXCEPTION( true, std::invalid_argument,
">>> ERROR (Intrepid2::CellTools::createHGradBasis): Cell topology not supported.");
}
}
return r_val;
}
public:
/** \brief Default constructor.
*/
CellTools() = default;
/** \brief Destructor
*/
~CellTools() = default;
//============================================================================================//
// //
// Jacobian, inverse Jacobian and Jacobian determinant //
// //
//============================================================================================//
/** \brief Computes the Jacobian matrix \e DF of the reference-to-physical frame map \e F.
There are two use cases:
\li Computes Jacobians \f$DF_{c}\f$ of the reference-to-physical map for \b all physical cells
in a cell workset on a \b single set of reference points stored in a rank-2 (P,D) array;
\li Computes Jacobians \f$DF_{c}\f$ of the reference-to-physical map for \b all physical cells
in a cell workset on \b multiple reference point sets having the same number of points,
indexed by cell ordinal, and stored in a rank-3 (C,P,D) array;
For a single point set in a rank-2 array (P,D) returns a rank-4 (C,P,D,D) array
such that
\f[
\mbox{jacobian}(c,p,i,j) = [DF_{c}(\mbox{points}(p))]_{ij} \quad c=0,\ldots, C
\f]
For multiple sets of reference points in a rank-3 (C,P,D) array returns
rank-4 (C,P,D,D) array such that
\f[
\mbox{jacobian}(c,p,i,j) = [DF_{c}(\mbox{points}(c,p))]_{ij} \quad c=0,\ldots, C
\f]
Requires pointer to HGrad basis that defines reference to physical cell mapping.
See Section \ref sec_cell_topology_ref_map_DF for definition of the Jacobian.
\warning
The points are not required to be in the reference cell associated with the specified
cell topology. CellTools provides several inclusion tests methods to check whether
or not the points are inside a reference cell.
\param jacobian [out] - rank-4 array with dimensions (C,P,D,D) with the Jacobians
\param points [in] - rank-2/3 array with dimensions (P,D)/(C,P,D) with the evaluation points
\param cellWorkset [in] - rank-3 container with logical dimensions (C,N,D) with the nodes of the cell workset
\param basis [in] - HGrad basis for reference to physical cell mapping
\param startCell [in] - first cell index in cellWorkset for which we should compute the Jacobian; corresponds to the 0 index in Jacobian and/or points container. Default: 0.
\param endCell [in] - first cell index in cellWorkset that we do not process; endCell - startCell must equal the extent of the Jacobian container in dimension 0. Default: -1, a special value that indicates the extent of the cellWorkset should be used.
*/
template<typename JacobianViewType,
typename PointViewType,
typename WorksetType,
typename HGradBasisType>
static void
setJacobian( JacobianViewType jacobian,
const PointViewType points,
const WorksetType worksetCell,
const Teuchos::RCP<HGradBasisType> basis,
const int startCell=0, const int endCell=-1);
/** \brief Computes the Jacobian matrix \e DF of the reference-to-physical frame map \e F.
There are two use cases:
\li Computes Jacobians \f$DF_{c}\f$ of the reference-to-physical map for \b all physical cells
in a cell workset on a \b single set of reference points stored in a rank-2 (P,D) array;
\li Computes Jacobians \f$DF_{c}\f$ of the reference-to-physical map for \b all physical cells
in a cell workset on \b multiple reference point sets having the same number of points,
indexed by cell ordinal, and stored in a rank-3 (C,P,D) array;
For a single point set in a rank-2 array (P,D) returns a rank-4 (C,P,D,D) array
such that
\f[
\mbox{jacobian}(c,p,i,j) = [DF_{c}(\mbox{points}(p))]_{ij} \quad c=0,\ldots, C
\f]
For multiple sets of reference points in a rank-3 (C,P,D) array returns
rank-4 (C,P,D,D) array such that
\f[
\mbox{jacobian}(c,p,i,j) = [DF_{c}(\mbox{points}(c,p))]_{ij} \quad c=0,\ldots, C
\f]
Requires pointer to HGrad basis that defines reference to physical cell mapping.
See Section \ref sec_cell_topology_ref_map_DF for definition of the Jacobian.
\warning
The points are not required to be in the reference cell associated with the specified
cell topology. CellTools provides several inclusion tests methods to check whether
or not the points are inside a reference cell.
\param jacobian [out] - rank-4 array with dimensions (C,P,D,D) with the Jacobians
\param cellWorkset [in] - rank-3 container with logical dimensions (C,N,D) with the nodes of the cell workset
\param gradients [in] - rank-3/4 array with dimensions (N,P,D)/(C,N,P,D) with the gradients of the physical-to-cell mapping
\param startCell [in] - first cell index in cellWorkset for which we should compute the Jacobian; corresponds to the 0 index in Jacobian and/or points container. Default: 0.
\param endCell [in] - first cell index in cellWorkset that we do not process; endCell - startCell must equal the extent of the Jacobian container in dimension 0. Default: -1, a special value that indicates the extent of the cellWorkset should be used.
*/
template<typename JacobianViewType,
typename BasisGradientsType,
typename WorksetType>
static void
setJacobian( JacobianViewType jacobian,
const WorksetType worksetCell,
const BasisGradientsType gradients,
const int startCell=0, const int endCell=-1);
/** \brief Computes the Jacobian matrix \e DF of the reference-to-physical frame map \e F.
There are two use cases:
\li Computes Jacobians \f$DF_{c}\f$ of the reference-to-physical map for \b all physical cells
in a cell workset on a \b single set of reference points stored in a rank-2 (P,D) array;
\li Computes Jacobians \f$DF_{c}\f$ of the reference-to-physical map for \b all physical cells
in a cell workset on \b multiple reference point sets having the same number of points,
indexed by cell ordinal, and stored in a rank-3 (C,P,D) array;
For a single point set in a rank-2 array (P,D) returns a rank-4 (C,P,D,D) array
such that
\f[
\mbox{jacobian}(c,p,i,j) = [DF_{c}(\mbox{points}(p))]_{ij} \quad c=0,\ldots, C
\f]
For multiple sets of reference points in a rank-3 (C,P,D) array returns
rank-4 (C,P,D,D) array such that
\f[
\mbox{jacobian}(c,p,i,j) = [DF_{c}(\mbox{points}(c,p))]_{ij} \quad c=0,\ldots, C
\f]
Requires cell topology with a reference cell. See Section \ref sec_cell_topology_ref_map_DF
for definition of the Jacobian.
\warning
The points are not required to be in the reference cell associated with the specified
cell topology. CellTools provides several inclusion tests methods to check whether
or not the points are inside a reference cell.
\param jacobian [out] - rank-4 array with dimensions (C,P,D,D) with the Jacobians
\param points [in] - rank-2/3 array with dimensions (P,D)/(C,P,D) with the evaluation points
\param cellWorkset [in] - rank-3 array with dimensions (C,N,D) with the nodes of the cell workset
\param cellTopo [in] - cell topology of the cells stored in \c cellWorkset
*/
template<typename JacobianViewType,
typename PointViewType,
typename WorksetCellViewType>
static void
setJacobian( JacobianViewType jacobian,
const PointViewType points,
const WorksetCellViewType worksetCell,
const shards::CellTopology cellTopo ) {
using nonConstPointValueType = typename PointViewType::non_const_value_type;
auto basis = createHGradBasis<nonConstPointValueType,nonConstPointValueType>(cellTopo);
setJacobian(jacobian,
points,
worksetCell,
basis);
}
/** \brief Computes the inverse of the Jacobian matrix \e DF of the reference-to-physical frame map \e F.
Returns rank-4 array with dimensions (C,P,D,D) such that
\f[
\mbox{jacobianInv}(c,p,*,*) = \mbox{jacobian}^{-1}(c,p,*,*) \quad c = 0,\ldots, C
\f]
\param jacobianInv [out] - rank-4 array with dimensions (C,P,D,D) with the inverse Jacobians
\param jacobian [in] - rank-4 array with dimensions (C,P,D,D) with the Jacobians
*/
template<typename JacobianInvViewType,
typename JacobianViewType>
static void
setJacobianInv( JacobianInvViewType jacobianInv,
const JacobianViewType jacobian );
/** \brief Computes the determinant of the Jacobian matrix \e DF of the reference-to-physical frame map \e F.
Returns rank-2 array with dimensions (C,P) such that
\f[
\mbox{jacobianDet}(c,p) = \mbox{det}(\mbox{jacobian}(c,p,*,*)) \quad c=0,\ldots, C
\f]
\param jacobianDet [out] - rank-2 array with dimensions (C,P) with Jacobian determinants
\param jacobian [in] - rank-4 array with dimensions (C,P,D,D) with the Jacobians
*/
template<typename JacobianDetViewType,
typename JacobianViewType>
static void
setJacobianDet( JacobianDetViewType jacobianDet,
const JacobianViewType jacobian );
/** \brief Allocates and returns a Data container suitable for storing determinants corresponding to the Jacobians in the Data container provided
\param jacobian [in] - data with shape (C,P,D,D), as returned by CellGeometry::allocateJacobianData()
\return Data container with shape (C,P)
*/
template<class PointScalar>
static Data<PointScalar,DeviceType> allocateJacobianDet( const Data<PointScalar,DeviceType> & jacobian );
/** \brief Allocates and returns a Data container suitable for storing inverses corresponding to the Jacobians in the Data container provided
\param jacobian [in] - data with shape (C,P,D,D), as returned by CellGeometry::allocateJacobianData()
\return Data container with shape (C,P,D,D)
*/
template<class PointScalar>
static Data<PointScalar,DeviceType> allocateJacobianInv( const Data<PointScalar,DeviceType> & jacobian );
/** \brief Computes determinants corresponding to the Jacobians in the Data container provided
\param jacobianDet [out] - data with shape (C,P), as returned by CellTools::allocateJacobianDet()
\param jacobian [in] - data with shape (C,P,D,D), as returned by CellGeometry::allocateJacobianData()
*/
template<class PointScalar>
static void setJacobianDet( Data<PointScalar,DeviceType> & jacobianDet,
const Data<PointScalar,DeviceType> & jacobian);
/** \brief Computes reciprocals of determinants corresponding to the Jacobians in the Data container provided
\param jacobianDet [out] - data with shape (C,P), as returned by CellTools::allocateJacobianDet()
\param jacobian [in] - data with shape (C,P,D,D), as returned by CellGeometry::allocateJacobianData()
*/
template<class PointScalar>
static void setJacobianDetInv( Data<PointScalar,DeviceType> & jacobianDet,
const Data<PointScalar,DeviceType> & jacobian);
/** \brief Computes determinants corresponding to the Jacobians in the Data container provided
\param jacobianInv [out] - data container with shape (C,P,D,D), as returned by CellTools::allocateJacobianInv()
\param jacobian [in] - data with shape (C,P,D,D), as returned by CellGeometry::allocateJacobianData()
*/
template<class PointScalar>
static void setJacobianInv( Data<PointScalar,DeviceType> & jacobianInv,
const Data<PointScalar,DeviceType> & jacobian);
/** \brief Multiplies the Jacobian with shape (C,P,D,D) by the reciprocals of the determinants, with shape (C,P), entrywise.
\param jacobianDividedByDet [out] - data container with shape (C,P,D,D), as returned by CellTools::allocateJacobianInv()
\param jacobian [in] - data container with shape (C,P,D,D), as returned by CellTools::allocateJacobianInv()
\param jacobianDetInv [in] - data with shape (C,P,D,D), as returned by CellGeometry::allocateJacobianData()
*/
template<class PointScalar>
static void setJacobianDividedByDet( Data<PointScalar,DeviceType> & jacobianDividedByDet,
const Data<PointScalar,DeviceType> & jacobian,
const Data<PointScalar,DeviceType> & jacobianDetInv);
//============================================================================================//
// //
// Node information //
// //
//============================================================================================//
// the node information can be used inside of kokkos functor and needs kokkos inline and
// exception should be an abort. for now, let's not decorate
/** \brief Computes the Cartesian coordinates of reference cell barycenter.
Requires cell topology with a reference cell.
\param cellCenter [out] - coordinates of the specified reference cell center
\param cell [in] - cell topology
*/
template<typename cellCenterValueType, class ...cellCenterProperties>
static void
getReferenceCellCenter( Kokkos::DynRankView<cellCenterValueType,cellCenterProperties...> cellCenter,
const shards::CellTopology cell );
/** \brief Retrieves the Cartesian coordinates of a reference cell vertex.
Requires cell topology with a reference cell.
\param cellVertex [out] - coordinates of the specified reference cell vertex
\param cell [in] - cell topology of the cell
\param vertexOrd [in] - vertex ordinal
*/
template<typename cellVertexValueType, class ...cellVertexProperties>
static void
getReferenceVertex( Kokkos::DynRankView<cellVertexValueType,cellVertexProperties...> cellVertex,
const shards::CellTopology cell,
const ordinal_type vertexOrd );
/** \brief Retrieves the Cartesian coordinates of all vertices of a reference subcell.
Returns rank-2 array with the Cartesian coordinates of the vertices of the
specified reference cell subcell. Requires cell topology with a reference cell.
\param subcellVertices [out] - rank-2 (V,D) array with the Cartesian coordinates of the reference subcell vertices
\param subcellDim [in] - dimension of the subcell; 0 <= subcellDim <= parentCell dimension
\param subcellOrd [in] - ordinal of the subcell
\param parentCell [in] - topology of the cell that owns the subcell
\remark When \c subcellDim = dimension of the \c parentCell this method returns the Cartesian
coordinates of the vertices of the reference cell itself.
Note that this requires \e subcellOrd=0.
*/
template<typename subcellVertexValueType, class ...subcellVertexProperties>
static void
getReferenceSubcellVertices( Kokkos::DynRankView<subcellVertexValueType,subcellVertexProperties...> subcellVertices,
const ordinal_type subcellDim,
const ordinal_type subcellOrd,
const shards::CellTopology parentCell );
/** \brief Retrieves the Cartesian coordinates of a reference cell node.
Returns Cartesian coordinates of a reference cell node. Requires cell topology
with a reference cell. Node coordinates are always returned as an (x,y,z)-triple
regardlesss of the actual topological cell dimension. The unused coordinates are
set to zero, e.g., node 0 of Line<2> is returned as {-1,0,0}.
\remark
Because the nodes of a cell with a base topology coincide with its vertices, for cells
with base topology this method is equivalent to Intrepid2::CellTools::getReferenceVertex.
\param cellNode [out] - coordinates of the specified reference vertex
\param cell [in] - cell topology of the cell
\param vertexOrd [in] - node ordinal
*/
template<typename cellNodeValueType, class ...cellNodeProperties>
static void
getReferenceNode( Kokkos::DynRankView<cellNodeValueType,cellNodeProperties...> cellNode,
const shards::CellTopology cell,
const ordinal_type nodeOrd );
/** \brief Retrieves the Cartesian coordinates of all nodes of a reference subcell.
Returns rank-2 array with the Cartesian coordinates of the nodes of the
specified reference cell subcell. Requires cell topology with a reference cell.
\param subcellNodes [out] - rank-2 (N,D) array with the Cartesian coordinates of the reference subcell nodes
\param subcellDim [in] - dimension of the subcell; 0 <= subcellDim <= parentCell dimension
\param subcellOrd [in] - ordinal of the subcell
\param parentCell [in] - topology of the cell that owns the subcell
\remark When \c subcellDim = dimension of the \c parentCell this method returns the Cartesian
coordinates of the nodes of the reference cell itself. Note that this requires \c subcellOrd=0.
*/
template<typename SubcellNodeViewType>
static void
getReferenceSubcellNodes( SubcellNodeViewType subcellNodes,
const ordinal_type subcellDim,
const ordinal_type subcellOrd,
const shards::CellTopology parentCell );
/** \brief Computes constant tangent vectors to edges of 2D or 3D reference cells.
Returns rank-1 array with dimension (D), D=2 or D=3; such that
\f[
{refEdgeTangent}(*) = \hat{\bf t}_i = {\partial\hat{\Phi}_i(t)\over\partial t}\,,
\f]
where \f$\hat{\Phi}_i : R =[-1,1]\mapsto \hat{\mathcal E}_i\f$ is the parametrization map
of the specified reference edge \f$\hat{\mathcal E}_i\f$, given by
\f[
\hat{\Phi}_i(t) = \left\{\begin{array}{ll}
(\hat{x}(t),\hat{y}(t),\hat{z}(t)) & \mbox{for 3D parent cells} \\[1ex]
(\hat{x}(t),\hat{y}(t)) & \mbox{for 2D parent cells} \\[1ex]
\end{array}\right.
\f]
The length of computed edge tangents is one-half the length of their associated edges:
\f[
|\hat{\bf t}_i | = {1\over 2} |\hat{\mathcal E}_i |\,.
\f]
Because the edges of all reference cells are always affine images of [-1,1],
the edge tangent is constant vector field.
\param refEdgeTangent [out] - rank-1 array (D) with the edge tangent; D = cell dimension
\param edgeOrd [in] - ordinal of the edge whose tangent is computed
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename RefEdgeTangentViewType>
static void
getReferenceEdgeTangent( RefEdgeTangentViewType refEdgeTangent,
const ordinal_type edgeOrd,
const shards::CellTopology parentCell );
/** \brief Computes pairs of constant tangent vectors to faces of a 3D reference cells.
Returns 2 rank-1 arrays with dimension (D), D=3, such that
\f[
{refFaceTanU}(*) = \hat{\bf t}_{i,u} = {\partial\hat{\Phi}_i(u,v)\over\partial u} =
\left({\partial\hat{x}(u,v)\over\partial u},
{\partial\hat{y}(u,v)\over\partial u},
{\partial\hat{z}(u,v)\over\partial u} \right) ;
\f]
\f[
{refFaceTanV}(*) = \hat{\bf t}_{i,v} = {\partial\hat{\Phi}_i(u,v)\over \partial v} =
\left({\partial\hat{x}(u,v)\over\partial v},
{\partial\hat{y}(u,v)\over\partial v},
{\partial\hat{z}(u,v)\over\partial v} \right)\,;
\f]
where \f$\hat{\Phi}_i: R \mapsto \hat{\mathcal F}_i\f$
is the parametrization map of the specified reference face \f$\hat{\mathcal F}_i\f$ given by
\f[
\hat{\Phi}_i(u,v) =(\hat{x}(u,v),\hat{y}(u,v),\hat{z}(u,v))
\f]
and
\f[
R = \left\{\begin{array}{rl}
\{(0,0),(1,0),(0,1)\} & \mbox{if $\hat{\mathcal F}_i$ is Triangle} \\[1ex]
[-1,1]\times [-1,1] & \mbox{if $\hat{\mathcal F}_i$ is Quadrilateral} \,.
\end{array}\right.
\f]
Because the faces of all reference cells are always affine images of \e R ,
the coordinate functions \f$\hat{x},\hat{y},\hat{z}\f$ of the parametrization map
are linear and the face tangents are constant vectors.
\param refFaceTanU [out] - rank-1 array (D) with (constant) tangent in u-direction
\param refFaceTanV [out] - rank-1 array (D) with (constant) tangent in v-direction
\param faceOrd [in] - ordinal of the face whose tangents are computed
\param parentCell [in] - cell topology of the parent 3D reference cell
*/
template<typename RefFaceTanViewType>
static void
getReferenceFaceTangents( RefFaceTanViewType refFaceTanU,
RefFaceTanViewType refFaceTanV,
const ordinal_type faceOrd,
const shards::CellTopology parentCell );
/** \brief Computes constant normal vectors to sides of 2D or 3D reference cells.
A side is defined as a subcell of dimension one less than that of its parent cell.
Therefore, sides of 2D cells are 1-subcells (edges) and sides of 3D cells
are 2-subcells (faces).
Returns rank-1 array with dimension (D), D = 2 or 3 such that
\f[
{refSideNormal}(*) = \hat{\bf n}_i =
\left\{\begin{array}{rl}
\displaystyle
\left({\partial\hat{\Phi}_i(t)\over\partial t}\right)^{\perp}
& \mbox{for 2D parent cells} \\[2ex]
\displaystyle
{\partial\hat{\Phi}_{i}\over\partial u} \times
{\partial\hat{\Phi}_{i}\over\partial v} & \mbox{for 3D parent cells}
\end{array}\right.
\f]
where \f$ (u_1,u_2)^\perp = (u_2, -u_1)\f$, and \f$\hat{\Phi}_i: R \mapsto \hat{\mathcal S}_i\f$
is the parametrization map of the specified reference side \f$\hat{\mathcal S}_i\f$ given by
\f[
\hat{\Phi}_i(u,v) =
\left\{\begin{array}{rl}
(\hat{x}(t),\hat{y}(t)) & \mbox{for 2D parent cells} \\[1ex]
(\hat{x}(u,v),\hat{y}(u,v),\hat{z}(u,v)) & \mbox{for 3D parent cells}
\end{array}\right.
\f]
For sides of 2D cells \e R=[-1,1] and for sides of 3D cells
\f[
R = \left\{\begin{array}{rl}
\{(0,0),(1,0),(0,1)\} & \mbox{if $\hat{\mathcal S}_i$ is Triangle} \\[1ex]
[-1,1]\times [-1,1] & \mbox{if $\hat{\mathcal S}_i$ is Quadrilateral} \,.
\end{array}\right.
\f]
For 3D cells the length of computed side normals is proportional to side area:
\f[
|\hat{\bf n}_i | = \left\{\begin{array}{rl}
2 \mbox{Area}(\hat{\mathcal F}_i) & \mbox{if $\hat{\mathcal F}_i$ is Triangle} \\[1ex]
\mbox{Area}(\hat{\mathcal F}_i) & \mbox{if $\hat{\mathcal F}_i$ is Quadrilateral} \,.
\end{array}\right.
\f]
For 2D cells the length of computed side normals is proportional to side length:
\f[
|\hat{\bf n}_i | = {1\over 2} |\hat{\mathcal F}_i |\,.
\f]
Because the sides of all reference cells are always affine images of \e R ,
the coordinate functions \f$\hat{x},\hat{y},\hat{z}\f$ of the parametrization maps
are linear and the side normal is a constant vector.
\remark
- For 3D cells the reference side normal coincides with the face normal computed by
Intrepid2::CellTools::getReferenceFaceNormal and these two methods are completely equivalent.
- For 2D cells the reference side normal is defined by \f$\hat{{\bf n}}= \hat{\bf t}^\perp = (t_2,-t_1)\f$
where \f$\hat{{\bf t}}=(t_1,t_2)\f$ is the tangent vector computed by
Intrepid2::CellTools::getReferenceEdgeTangent. Therefore, the pair
\f$(\hat{{\bf n}},\hat{{\bf t}})\f$ is positively oriented.
\param refSideNormal [out] - rank-1 array (D) with (constant) side normal
\param sideOrd [in] - ordinal of the side whose normal is computed
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename RefSideNormalViewType>
static void
getReferenceSideNormal( RefSideNormalViewType refSideNormal,
const ordinal_type sideOrd,
const shards::CellTopology parentCell );
/** \brief Computes constant normal vectors to faces of 3D reference cell.
Returns rank-1 array with dimension (D), D=3 such that
\f[
{refFaceNormal}(*) = \hat{\bf n}_i = {\partial\hat{\Phi}_{i}\over\partial u} \times
{\partial\hat{\Phi}_{i}\over\partial v}
\f]
where \f$\hat{\Phi}_i: R \mapsto \hat{\mathcal F}_i\f$
is the parametrization map of the specified reference face \f$\hat{\mathcal F}_i\f$ given by
\f[
\hat{\Phi}_i(u,v) =(\hat{x}(u,v),\hat{y}(u,v),\hat{z}(u,v))
\f]
and
\f[
R = \left\{\begin{array}{rl}
\{(0,0),(1,0),(0,1)\} & \mbox{if ${\mathcal F}$ is Triangle} \\[1ex]
[-1,1]\times [-1,1] & \mbox{if ${\mathcal F}$ is Quadrilateral} \,.
\end{array}\right.
\f]
The length of computed face normals is proportional to face area:
\f[
|\hat{\bf n}_i | = \left\{\begin{array}{rl}
2 \mbox{Area}(\hat{\mathcal F}_i) & \mbox{if $\hat{\mathcal F}_i$ is Triangle} \\[1ex]
\mbox{Area}(\hat{\mathcal F}_i) & \mbox{if $\hat{\mathcal F}_i$ is Quadrilateral} \,.
\end{array}\right.
\f]
Because the faces of all reference cells are always affine images of \e R ,
the coordinate functions \f$\hat{x},\hat{y},\hat{z}\f$ of the parametrization map
are linear and the face normal is a constant vector.
\remark
The method Intrepid2::CellTools::getReferenceFaceTangents computes the reference face tangents
\f${\partial\hat{\Phi}_{i}/\partial u}\f$ and \f${\partial\hat{\Phi}_{i}/\partial v}\f$.
\param refFaceNormal [out] - rank-1 array (D) with (constant) face normal
\param faceOrd [in] - ordinal of the face whose normal is computed
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename RefFaceNormalViewType>
static void
getReferenceFaceNormal( RefFaceNormalViewType refFaceNormal,
const ordinal_type faceOrd,
const shards::CellTopology parentCell );
/** \brief Computes non-normalized tangent vectors to physical edges in an edge workset
\f$\{\mathcal{E}_{c,i}\}_{c=0}^{N}\f$; (see \ref sec_cell_topology_subcell_wset for definition of edge worksets).
For every edge in the workset the tangents are computed at the points
\f${\bf x}_p = F_c(\hat{\Phi}_i(t_p))\in\mathcal{E}_{c,i}\f$ that are images of points
from <var>R=[-1,1]</var> on edge \f$\mathcal{E}_{c,i}\f$. Returns rank-3 array with
dimensions (C,P,D1), D1=2 or D1=3 such that
\f[
{edgeTangents}(c,p,d) =
DF_c(\hat{\Phi}_i(t_p))\cdot {\partial{\hat{\Phi}}_{i}(t_p)\over\partial t}\,; \qquad t_p \in R
\f]
In this formula:
\li \f$ DF_c \f$ is Jacobian of parent cell \f${\mathcal C}\f$ that owns physical edge \f${\mathcal E}_{c,i}\f$;
\li \f$ {\partial{\hat{\Phi}}_{i}/\partial t}\f$ is the (constant) tangent to reference edge
\f$\hat{\mathcal E}_i\f$; see Intrepid2::CellTools::getReferenceEdgeTangent that has the
same local ordinal as the edges in the workset;
\li \f$ \hat{\Phi}_i R\mapsto\hat{\mathcal E}_i \f$ is parametrization of reference edge \f$\hat{\mathcal E}_i\f$;
\warning
\c worksetJacobians must contain the values of \f$DF_c(\hat{\Phi}_i(t_p))\f$,
where \f$ t_p \in R=[-1,1] \f$, i.e., Jacobians of the parent cells evaluated at points
that are located on reference edge \f$\hat{\mathcal E}_i\f$ having the same local ordinal as
the edges in the workset.
\param edgeTangents [out] - rank-3 array (C,P,D1) with tangents on workset edges
\param worksetJacobians [in] - rank-4 array (C,P,D1,D1) with Jacobians evaluated at ref. edge points
\param worksetEdgeOrd [in] - edge ordinal, relative to ref. cell, of the edge workset
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename edgeTangentValueType, class ...edgeTangentProperties,
typename worksetJacobianValueType, class ...worksetJacobianProperties>
static void
getPhysicalEdgeTangents( Kokkos::DynRankView<edgeTangentValueType,edgeTangentProperties...> edgeTangents,
const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
const ordinal_type worksetEdgeOrd,
const shards::CellTopology parentCell );
/** \brief Computes non-normalized tangent vectors to physical edges in an edge workset
\f$\{\mathcal{E}_{c,i}\}_{c=0}^{N}\f$; (see \ref sec_cell_topology_subcell_wset for definition of edge worksets).
It is similar to the <var>CellTools::getPhysicalEdgeTangents</var> function above, with the difference that the edge ordinal can change from point to point,
and it is provided by the rank-2 input array <var><b>worksetEdgeOrds</b></var>, indexed by (C,P).
\param edgeTangents [out] - rank-3 array (C,P,D1) with tangents on workset edges
\param worksetJacobians [in] - rank-4 array (C,P,D1,D1) with Jacobians evaluated at ref. edge points
\param worksetEdgeOrds [in] - rank-2 array (C,P) with edge ordinals, relative to ref. cell, of the edge workset
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename edgeTangentValueType, class ...edgeTangentProperties,
typename worksetJacobianValueType, class ...worksetJacobianProperties,
typename edgeOrdValueType, class ...edgeOrdProperties>
static void
getPhysicalEdgeTangents( Kokkos::DynRankView<edgeTangentValueType,edgeTangentProperties...> edgeTangents,
const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
const Kokkos::DynRankView<edgeOrdValueType,edgeOrdProperties...> worksetEdgeOrds,
const shards::CellTopology parentCell );
/** \brief Computes non-normalized tangent vector pairs to physical faces in a face workset
\f$\{\mathcal{F}_{c,i}\}_{c=0}^{N}\f$; (see \ref sec_cell_topology_subcell_wset for definition of face worksets).
For every face in the workset the tangents are computed at the points
\f${\bf x}_p = F_c(\hat{\Phi}_i(u_p,v_p))\in\mathcal{F}_{c,i}\f$ that are images of points
from the parametrization domain \e R on face \f$\mathcal{F}_{c,i}\f$.
Returns 2 rank-3 arrays with dimensions (C,P,D), D=3 such that
\f[
{faceTanU}(c,p,d) = DF_c(\hat{\Phi}_i(u_p, v_p))\cdot {\partial\hat{\Phi}_i\over\partial u};\qquad
{faceTanV}(c,p,d) = DF_c(\hat{\Phi}_i(u_p, v_p))\cdot {\partial\hat{\Phi}_{i}\over\partial v}\,;
\qquad (u_p, v_p) \in R \,.
\f]
In this formula:
\li \f$ DF_c \f$ is Jacobian of parent cell \f${\mathcal C}\f$ that owns physical face \f${\mathcal F}_{c,i}\f$;
\li \f$ {\partial\hat{\Phi}_i/\partial u}, {\partial\hat{\Phi}_i/\partial v}\f$ are the (constant)
tangents on reference face \f$\hat{\mathcal F}_i\f$; see Intrepid2::CellTools::getReferenceFaceTangents;
that has the same local ordinal as the faces in the workset;
\li \f$ \hat{\Phi}_i : R\mapsto \hat{\mathcal F}_i\f$ is parametrization of reference face \f$\hat{\mathcal F}_i\f$;
\li \e R is the parametrization domain for reference face \f$\hat{\mathcal F}_i\f$:
\f[
R =
\left\{\begin{array}{rl}
\{(0,0),(1,0),(0,1)\} & \mbox{if $\hat{\mathcal F}_i$ is Triangle} \\[1ex]
[-1,1]\times [-1,1] & \mbox{if $\hat{\mathcal F}_i$ is Quadrilateral}
\end{array}\right.
\f]
\warning
\c worksetJacobians must contain the values of \f$DF_c(\hat{\Phi}_i(u_p,v_p))\f$,
where \f$(u_p,v_p)\in R\f$, i.e., Jacobians of the parent cells evaluated at points
that are located on reference face \f$\hat{\mathcal F}_i\f$ having the same local ordinal as
the faces in the workset.
\param faceTanU [out] - rank-3 array (C,P,D), image of ref. face u-tangent at workset faces
\param faceTanV [out] - rank-3 array (C,P,D), image of ref. face u-tangent at workset faces
\param worksetJacobians [in] - rank-4 array (C,P,D,D) with Jacobians at ref. face points
\param worksetFaceOrd [in] - face ordinal, relative to ref. cell, of the face workset
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename faceTanValueType, class ...faceTanProperties,
typename worksetJacobianValueType, class ...worksetJacobianProperties>
static void
getPhysicalFaceTangents( Kokkos::DynRankView<faceTanValueType,faceTanProperties...> faceTanU,
Kokkos::DynRankView<faceTanValueType,faceTanProperties...> faceTanV,
const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
const ordinal_type worksetFaceOrd,
const shards::CellTopology parentCell );
/** \brief Computes non-normalized tangent vector pairs to physical faces in a face workset
\f$\{\mathcal{F}_{c,i}\}_{c=0}^{N}\f$; (see \ref sec_cell_topology_subcell_wset for definition of face worksets).
It is similar to the <var>CellTools::getPhysicalFaceTangents</var> function above, with the difference that the face ordinal can change from point to point,
and it is provided by the rank-2 input array <var><b>worksetFaceOrds</b></var>, indexed by (C,P).
\param faceTanU [out] - rank-3 array (C,P,D), image of ref. face u-tangent at workset faces
\param faceTanV [out] - rank-3 array (C,P,D), image of ref. face u-tangent at workset faces
\param worksetJacobians [in] - rank-4 array (C,P,D,D) with Jacobians at ref. face points
\param worksetFaceOrds [in] - rank-2 array (C,P) with face ordinals, relative to ref. cell, of the face workset
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename faceTanValueType, class ...faceTanProperties,
typename worksetJacobianValueType, class ...worksetJacobianProperties,
typename faceOrdValueType, class ...faceOrdProperties>
static void
getPhysicalFaceTangents( Kokkos::DynRankView<faceTanValueType,faceTanProperties...> faceTanU,
Kokkos::DynRankView<faceTanValueType,faceTanProperties...> faceTanV,
const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
const Kokkos::DynRankView<faceOrdValueType,faceOrdProperties...> worksetFaceOrds,
const shards::CellTopology parentCell );
/** \brief Computes non-normalized normal vectors to physical sides in a side workset
\f$\{\mathcal{S}_{c,i}\}_{c=0}^{N}\f$.
For every side in the workset the normals are computed at the points
\f${\bf x}_p = F_c(\hat{\Phi}_i(P_p))\in\mathcal{S}_{c,i}\f$ that are images of points
from the parametrization domain \e R on side \f$\mathcal{S}_{c,i}\f$.
A side is defined as a subcell of dimension one less than that of its parent cell.
Therefore, sides of 2D cells are 1-subcells (edges) and sides of 3D cells are 2-subcells (faces).
Returns rank-3 array with dimensions (C,P,D), D = 2 or 3, such that
\f[
{sideNormals}(c,p,d) =
\left\{\begin{array}{crl}
\displaystyle
\left(DF_c(\hat{\Phi}_i(t_p))\cdot
{\partial{\hat{\Phi}}_{i}(t_p)\over\partial t}\right)^{\perp} & t_p\in R
& \mbox{for 2D parent cells} \\[2ex]
\displaystyle
\left( DF_c(\hat{\Phi}_i(u_p, v_p))\cdot {\partial\hat{\Phi}_i\over\partial u}\right) \times
\left( DF_c(\hat{\Phi}_i(u_p, v_p))\cdot {\partial\hat{\Phi}_i\over\partial v}\right) \,;
& (u_p, v_p) \in R & \mbox{for 3D parent cells}
\end{array}\right.
\f]
In this formula:
\li \f$ DF_c \f$ is Jacobian of parent cell \f${\mathcal C}\f$ that owns physical side \f${\mathcal S}_{c,i}\f$;
\li For 2D cells: \f$ {\partial{\hat{\Phi}}_{i}/\partial t}\f$ is the (constant) tangent to reference side (edge)
\f$\hat{\mathcal S}_i\f$; see Intrepid2::CellTools::getReferenceEdgeTangent, that has the
same local ordinal as the sides in the workset;
\li For 3D cells: \f$ {\partial\hat{\Phi}_i/\partial u}, {\partial\hat{\Phi}_i/\partial v}\f$ are the (constant)
tangents on reference side (face) \f$\hat{\mathcal S}_i\f$; see Intrepid2::CellTools::getReferenceFaceTangents,
that has the same local ordinal as the sides in the workset;
\li \f$ \hat{\Phi}_i : R\mapsto \hat{\mathcal S}_i\f$ is parametrization of reference side \f$\hat{\mathcal S}_i\f$;
\li \e R is the parametrization domain for reference side \f$\hat{\mathcal S}_i\f$. For
2D parent cells \e R=[-1,1] and for 3D parent cells
\f[
R = \left\{\begin{array}{rl}
\{(0,0),(1,0),(0,1)\} & \mbox{if $\hat{\mathcal S}_i$ is Triangle} \\[1ex]
[-1,1]\times [-1,1] & \mbox{if $\hat{\mathcal S}_i$ is Quadrilateral}
\end{array}\right.
\f]
\remark
- For 3D cells the physical side normals coincides with the face normals computed by
Intrepid2::CellTools::getPhysicalFaceNormals and these two methods are completely equivalent.
- For 2D cells the physical side normals are defined by \f${\bf n}=(t_2,-t_1)\f$
where \f${\bf t}=(t_1,t_2)\f$ are the physical edge tangents computed by
Intrepid2::CellTools::getPhysicalEdgeTangents. Therefore, the pairs \f$({\bf n},{\bf t})\f$ are positively oriented.
\warning
\c worksetJacobians must contain the values of \f$DF_c(\hat{\Phi}_i(P_p))\f$,
where \f$P_p\in R\f$, i.e., Jacobians of the parent cells evaluated at points
that are located on reference side \f$\hat{\mathcal S}_i\f$ having the same local ordinal as
the sides in the workset.
\param sideNormals [out] - rank-3 array (C,P,D), normals at workset sides
\param worksetJacobians [in] - rank-4 array (C,P,D,D) with Jacobians at ref. side points
\param worksetSideOrd [in] - side ordinal, relative to ref. cell, of the side workset
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename sideNormalValueType, class ...sideNormalProperties,
typename worksetJacobianValueType, class ...worksetJacobianProperties>
static void
getPhysicalSideNormals( Kokkos::DynRankView<sideNormalValueType,sideNormalProperties...> sideNormals,
const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
const ordinal_type worksetSideOrd,
const shards::CellTopology parentCell );
/** \brief Computes non-normalized normal vectors to physical sides in a side workset
\f$\{\mathcal{S}_{c,i}\}_{c=0}^{N}\f$.
It is similar to the <var>CellTools::getPhysicalSideNormals</var> function above, with the difference that the side ordinal can change from point to point,
and it is provided by the rank-2 input array <var><b>worksetSideOrds</b></var>, indexed by (C,P).
\param sideNormals [out] - rank-3 array (C,P,D), normals at workset sides
\param worksetJacobians [in] - rank-4 array (C,P,D,D) with Jacobians at ref. side points
\param worksetSideOrds [in] - rank-2 array (C,P) with side ordinals, relative to ref. cell, of the side workset
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename sideNormalValueType, class ...sideNormalProperties,
typename worksetJacobianValueType, class ...worksetJacobianProperties,
typename edgeOrdValueType, class ...edgeOrdProperties>
static void
getPhysicalSideNormals( Kokkos::DynRankView<sideNormalValueType,sideNormalProperties...> sideNormals,
const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
const Kokkos::DynRankView<edgeOrdValueType,edgeOrdProperties...> worksetSideOrds,
const shards::CellTopology parentCell );
/** \brief Computes non-normalized normal vectors to physical faces in a face workset
\f$\{\mathcal{F}_{c,i}\}_{c=0}^{N}\f$; (see \ref sec_cell_topology_subcell_wset for definition of face worksets).
For every face in the workset the normals are computed at the points
\f${\bf x}_p = F_c(\hat{\Phi}_i(u_p,v_p))\in\mathcal{F}_{c,i}\f$ that are images of points
from the parametrization domain \e R on face \f$\mathcal{F}_{c,i}\f$.
Returns rank-3 array with dimensions (C,P,D), D=3, such that
\f[
{faceNormals}(c,p,d) =
\left( DF_c(\hat{\Phi}_i(u_p, v_p))\cdot {\partial\hat{\Phi}_i\over\partial u}\right) \times
\left( DF_c(\hat{\Phi}_i(u_p, v_p))\cdot {\partial\hat{\Phi}_i\over\partial v}\right) \,;
\qquad (u_p, v_p) \in R \,.
\f]
In this formula:
\li \f$ DF_c \f$ is Jacobian of parent cell \f${\mathcal C}\f$ that owns physical face \f${\mathcal F}_{c,i}\f$;
\li \f$ {\partial\hat{\Phi}_i/\partial u}, {\partial\hat{\Phi}_i/\partial v}\f$ are the (constant)
tangents on reference face \f$\hat{\mathcal F}_i\f$; see Intrepid2::CellTools::getReferenceFaceTangents;
that has the same local ordinal as the faces in the workset;
\li \f$ \hat{\Phi}_i : R\mapsto \hat{\mathcal F}_i\f$ is parametrization of reference face \f$\hat{\mathcal F}_i\f$;
\li \e R is the parametrization domain for reference face \f$\hat{\mathcal F}_i\f$:
\f[
R = \left\{\begin{array}{rl}
\{(0,0),(1,0),(0,1)\} & \mbox{if $\hat{\mathcal F}_i$ is Triangle} \\[1ex]
[-1,1]\times [-1,1] & \mbox{if $\hat{\mathcal F}_i$ is Quadrilateral}
\end{array}\right.
\f]
\warning
\c worksetJacobians must contain the values of \f$DF_c(\hat{\Phi}_i(u_p,v_p))\f$,
where \f$(u_p,v_p)\in R\f$, i.e., Jacobians of the parent cells evaluated at points
that are located on reference face \f$\hat{\mathcal F}_i\f$ having the same local ordinal as
the faces in the workset.
\param faceNormals [out] - rank-3 array (C,P,D), normals at workset faces
\param worksetJacobians [in] - rank-4 array (C,P,D,D) with Jacobians at ref. face points
\param worksetFaceOrd [in] - face ordinal, relative to ref. cell, of the face workset
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename faceNormalValueType, class ...faceNormalProperties,
typename worksetJacobianValueType, class ...worksetJacobianProperties>
static void
getPhysicalFaceNormals( Kokkos::DynRankView<faceNormalValueType,faceNormalProperties...> faceNormals,
const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
const ordinal_type worksetFaceOrd,
const shards::CellTopology parentCell );
/** \brief Computes non-normalized normal vectors to physical faces in a face workset
\f$\{\mathcal{F}_{c,i}\}_{c=0}^{N}\f$; (see \ref sec_cell_topology_subcell_wset for definition of face worksets).
It is similar to the <var>CellTools::getPhysicalSideNormals</var> function above, with the difference that the side ordinal can change from point to point,
and it is provided by the rank-2 input array <var><b>worksetSideOrds</b></var>, indexed by (C,P).
\param faceNormals [out] - rank-3 array (C,P,D), normals at workset faces
\param worksetJacobians [in] - rank-4 array (C,P,D,D) with Jacobians at ref. face points
\param worksetFaceOrds [in] - rank-2 array (C,P) with face ordinals, relative to ref. cell, of the face workset
\param parentCell [in] - cell topology of the parent reference cell
*/
template<typename faceNormalValueType, class ...faceNormalProperties,
typename worksetJacobianValueType, class ...worksetJacobianProperties,
typename faceOrdValueType, class ...faceOrdProperties>
static void
getPhysicalFaceNormals( Kokkos::DynRankView<faceNormalValueType,faceNormalProperties...> faceNormals,
const Kokkos::DynRankView<worksetJacobianValueType,worksetJacobianProperties...> worksetJacobians,
const Kokkos::DynRankView<faceOrdValueType,faceOrdProperties...> worksetFaceOrds,
const shards::CellTopology parentCell );
//============================================================================================//
// //
// Reference-to-physical frame mapping and its inverse //
// //
//============================================================================================//
/** \brief Computes \e F, the reference-to-physical frame map.
There are 2 use cases:
\li Applies \f$ F_{c} \f$ for \b all cells in a cell workset to a \b single point set stored
in a rank-2 (P,D) array;
\li Applies \f$ F_{c} \f$ for \b all cells in a cell workset to \b multiple point sets having
the same number of points, indexed by cell ordinal, and stored in a rank-3 (C,P,D) array;
For a single point set in a rank-2 array (P,D) returns a rank-3 (C,P,D) array such that
\f[
\mbox{physPoints}(c,p,d) = \Big(F_c(\mbox{refPoint}(p,*)) \Big)_d \quad c=0,\ldots, C
\f]
For multiple point sets in a rank-3 (C,P,D) array returns a rank-3 (C,P,D) array such that
\f[
\mbox{physPoints}(c,p,d) = \Big(F_c(\mbox{refPoint}(c,p,*)) \Big)_d \quad c=0,\ldots, C
\f]
This corresponds to mapping multiple sets of reference points to a matching number of
physical cells.
Requires pointer to HGrad basis that defines reference to physical cell mapping.
See Section \ref sec_cell_topology_ref_map for definition of the mapping function.
\warning
The array \c refPoints represents an arbitrary set of points in the reference
frame that are not required to be in the reference cell corresponding to the
specified cell topology. As a result, the images of these points under a given
reference-to-physical map are not necessarily contained in the physical cell that
is the image of the reference cell under that map. CellTools provides several
inclusion tests methods to check whether or not the points are inside a reference cell.
\param physPoints [out] - rank-3 array with dimensions (C,P,D) with the images of the ref. points
\param refPoints [in] - rank-3/2 array with dimensions (C,P,D)/(P,D) with points in reference frame
\param cellWorkset [in] - rank-3 container with logical dimensions (C,N,D) with the nodes of the cell workset
\param basis [in] - pointer to HGrad basis used in reference-to-physical cell mapping
*/
template<typename PhysPointValueType,
typename RefPointValueType,
typename WorksetType,