43 #include "PanzerDiscFE_config.hpp"
47 #include "Kokkos_ViewFactory.hpp"
49 #include "Intrepid2_Utils.hpp"
50 #include "Intrepid2_FunctionSpaceTools.hpp"
51 #include "Intrepid2_Orientation.hpp"
52 #include "Intrepid2_OrientationTools.hpp"
57 template <
typename Scalar>
60 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> &
jac,
61 const PHX::MDField<Scalar,Cell,IP> & jac_det,
62 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
63 const int in_num_cells)
65 PHX::MDField<Scalar,Cell,IP> weighted_measure;
66 PHX::MDField<Scalar,Cell,NODE,Dim> vertex_coordinates;
67 build_weighted =
false;
68 evaluateValues(cub_points,jac,jac_det,jac_inv,weighted_measure,vertex_coordinates,
false,in_num_cells);
71 template <
typename Scalar>
74 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> &
jac,
75 const PHX::MDField<Scalar,Cell,IP> & jac_det,
76 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
77 const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
78 const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
79 bool use_vertex_coordinates,
80 const int in_num_cells)
84 int num_dim = basis_layout->dimension();
88 evaluateReferenceValues(cub_points,compute_derivatives,use_vertex_coordinates);
90 const int num_cells = in_num_cells < 0 ? jac.extent(0) : in_num_cells;
91 const std::pair<int,int> cell_range(0,num_cells);
96 auto s_basis_scalar = Kokkos::subview(basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
97 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
98 HGRADtransformVALUE(s_basis_scalar,
99 basis_ref_scalar.get_view());
102 auto s_weighted_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
103 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
104 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
105 multiplyMeasure(s_weighted_basis_scalar,
111 auto s_basis_vector = Kokkos::subview(basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
112 auto s_jac_inv = Kokkos::subview(jac_inv.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
113 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
114 HCURLtransformVALUE(s_basis_vector,
116 basis_ref_vector.get_view());
119 auto s_weighted_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
120 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
121 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
122 multiplyMeasure(s_weighted_basis_vector,
129 auto s_basis_vector = Kokkos::subview(basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
130 auto s_jac = Kokkos::subview(jac.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
131 auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
132 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
133 HDIVtransformVALUE(s_basis_vector,
136 basis_ref_vector.get_view());
139 auto s_weighted_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
140 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
141 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
142 multiplyMeasure(s_weighted_basis_vector,
149 auto s_basis_scalar = Kokkos::subview(basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
150 auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
151 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
152 HVOLtransformVALUE(s_basis_scalar,
154 basis_ref_scalar.get_view());
157 auto s_weighted_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
158 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
159 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
160 multiplyMeasure(s_weighted_basis_scalar,
168 auto s_grad_basis = Kokkos::subview(grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
169 auto s_jac_inv = Kokkos::subview(jac_inv.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
170 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
171 HGRADtransformGRAD(s_grad_basis,
173 grad_basis_ref.get_view());
176 auto s_weighted_grad_basis = Kokkos::subview(weighted_grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
177 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
178 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
179 multiplyMeasure(s_weighted_grad_basis,
185 auto s_curl_basis_scalar = Kokkos::subview(curl_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
186 auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
187 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
188 HDIVtransformDIV(s_curl_basis_scalar,
192 curl_basis_ref_scalar.get_view());
195 auto s_weighted_curl_basis_scalar = Kokkos::subview(weighted_curl_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
196 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
197 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
198 multiplyMeasure(s_weighted_curl_basis_scalar,
200 s_curl_basis_scalar);
204 auto s_curl_basis_vector = Kokkos::subview(curl_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
205 auto s_jac = Kokkos::subview(jac.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
206 auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
207 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
208 HCURLtransformCURL(s_curl_basis_vector,
211 curl_basis_ref_vector.get_view());
214 auto s_weighted_curl_basis_vector = Kokkos::subview(weighted_curl_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
215 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
216 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
217 multiplyMeasure(s_weighted_curl_basis_vector,
219 s_curl_basis_vector);
223 auto s_div_basis = Kokkos::subview(div_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
224 auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
225 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
226 HDIVtransformDIV(s_div_basis,
228 div_basis_ref.get_view());
231 auto s_weighted_div_basis = Kokkos::subview(weighted_div_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
232 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
233 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
234 multiplyMeasure(s_weighted_div_basis,
242 if(use_vertex_coordinates) {
256 auto s_basis_coordinates = Kokkos::subview(basis_coordinates.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
257 auto s_vertex_coordinates = Kokkos::subview(vertex_coordinates.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
258 Intrepid2::CellTools<PHX::Device::execution_space> cell_tools;
259 cell_tools.mapToPhysicalFrame(s_basis_coordinates,
260 basis_coordinates_ref.get_view(),
261 s_vertex_coordinates,
262 intrepid_basis->getBaseCellTopology());
272 template <
typename Scalar>
275 const int in_num_cells)
281 using coordsScalarType =
typename Intrepid2::Basis<PHX::Device::execution_space,Scalar,Scalar>::scalarType;
282 auto dyn_basis_coordinates_ref = af.
buildArray<coordsScalarType,
BASIS,
Dim>(
"basis_coordinates_ref",
283 basis_coordinates_ref.extent(0),
284 basis_coordinates_ref.extent(1));
285 intrepid_basis->getDofCoords(dyn_basis_coordinates_ref.get_view());
288 for (
int i = 0; i < basis_coordinates_ref.extent_int(0); ++i)
289 for (
int j = 0; j < basis_coordinates_ref.extent_int(1); ++j)
290 basis_coordinates_ref(i,j) = dyn_basis_coordinates_ref(i,j);
292 const int num_cells = in_num_cells < 0 ? vertex_coordinates.extent(0) : in_num_cells;
293 const std::pair<int,int> cell_range(0,num_cells);
294 const auto s_basis_coordinates = Kokkos::subview(basis_coordinates.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
295 const auto s_vertex_coordinates = Kokkos::subview(vertex_coordinates.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
297 Intrepid2::CellTools<PHX::Device::execution_space> cell_tools;
298 cell_tools.mapToPhysicalFrame(s_basis_coordinates,
299 basis_coordinates_ref.get_view(),
300 s_vertex_coordinates,
301 intrepid_basis->getBaseCellTopology());
308 template <
typename Scalar>
311 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> &
jac,
312 const PHX::MDField<Scalar,Cell,IP> & jac_det,
313 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
314 const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
315 const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
316 bool use_vertex_coordinates,
317 const int in_num_cells)
323 evaluateValues_Const(cub_points,jac_inv,weighted_measure,in_num_cells);
325 evaluateValues_HVol(cub_points,jac_det,jac_inv,weighted_measure,in_num_cells);
327 evaluateValues_HGrad(cub_points,jac_inv,weighted_measure,in_num_cells);
329 evaluateValues_HCurl(cub_points,jac,jac_det,jac_inv,weighted_measure,in_num_cells);
331 evaluateValues_HDiv(cub_points,jac,jac_det,weighted_measure,in_num_cells);
336 if(use_vertex_coordinates) {
338 evaluateBasisCoordinates(vertex_coordinates);
343 template <
typename Scalar>
346 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
347 const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
348 const int in_num_cells)
353 typedef Intrepid2::FunctionSpaceTools<PHX::Device::execution_space> fst;
358 const int num_points = basis_layout->numPoints();
360 const int num_dim = basis_layout->dimension();
361 const int num_cells = in_num_cells < 0 ? basis_layout->numCells() : in_num_cells;
364 auto cell_cub_points = af.
buildStaticArray<Scalar,IP,
Dim>(
"cell_cub_points",num_points,num_dim);
365 auto cell_grad_basis = af.
buildStaticArray<Scalar,Cell,BASIS,IP,Dim>(
"cell_grad_basis",1,num_basis,num_points,num_dim);
366 auto cell_jac_inv = af.
buildStaticArray<Scalar,Cell,IP,Dim,Dim>(
"cell_jac_inv",1,num_points,num_dim,num_dim);
368 auto cell_basis_ref_scalar = af.
buildStaticArray<Scalar,BASIS,IP>(
"cell_basis_ref_scalar",num_basis,num_points);
369 auto cell_grad_basis_ref = af.
buildStaticArray<Scalar,BASIS,IP,Dim>(
"cell_grad_basis_ref",num_basis,num_points,num_dim);
371 for(
int cell=0;cell<num_cells;++cell){
376 for(
int p=0;p<num_points;++p)
377 for(
int d=0;d<num_dim;++d)
378 for(
int d2=0;d2<num_dim;++d2)
379 cell_jac_inv(0,p,d,d2)=jac_inv(cell,p,d,d2);
380 for(
int p=0;p<num_points;++p)
381 for(
int d=0;d<num_dim;++d)
382 cell_cub_points(p,d)=cub_points(cell,p,d);
387 intrepid_basis->getValues(cell_basis_ref_scalar.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_VALUE);
389 if(compute_derivatives){
390 Kokkos::deep_copy(cell_grad_basis_ref.get_view(),0.0);
396 fst::HGRADtransformVALUE(cell_basis_scalar.get_view(),cell_basis_ref_scalar.get_view());
397 for(
int b=0;b<num_basis;++b)
398 for(
int p=0;p<num_points;++p)
399 basis_scalar(cell,b,p)=cell_basis_scalar(0,b,p);
401 if(compute_derivatives){
402 fst::HGRADtransformGRAD(cell_grad_basis.get_view(),cell_jac_inv.get_view(),cell_grad_basis_ref.get_view());
403 for(
int b=0;b<num_basis;++b)
404 for(
int p=0;p<num_points;++p)
405 for(
int d=0;d<num_dim;++d)
406 grad_basis(cell,b,p,d)=cell_grad_basis(0,b,p,d);
413 const std::pair<int,int> cell_range(0,num_cells);
414 auto s_weighted_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
415 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
416 auto s_basis_scalar = Kokkos::subview(basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
417 fst::multiplyMeasure(s_weighted_basis_scalar,s_weighted_measure,s_basis_scalar);
418 if(compute_derivatives){
419 auto s_weighted_grad_basis = Kokkos::subview(weighted_grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
420 auto s_grad_basis = Kokkos::subview(grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
421 fst::multiplyMeasure(s_weighted_grad_basis,s_weighted_measure,s_grad_basis);
428 template <
typename Scalar>
431 const PHX::MDField<Scalar,Cell,IP> & jac_det,
432 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & ,
433 const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
434 const int in_num_cells)
439 typedef Intrepid2::FunctionSpaceTools<PHX::Device::execution_space> fst;
444 const int num_points = basis_layout->numPoints();
446 const int num_dim = basis_layout->dimension();
447 const int num_cells = in_num_cells < 0 ? basis_layout->numCells() : in_num_cells;
450 auto cell_cub_points = af.
buildStaticArray<Scalar,IP,
Dim>(
"cell_cub_points",num_points,num_dim);
451 auto cell_jac_det = af.
buildStaticArray<Scalar,Cell,IP>(
"cell_jac_det",1,num_points);
453 auto cell_basis_ref_scalar = af.
buildStaticArray<Scalar,BASIS,IP>(
"cell_basis_ref_scalar",num_basis,num_points);
455 for(
int cell=0;cell<num_cells;++cell){
460 for(
int p=0;p<num_points;++p)
461 cell_jac_det(0,p)=jac_det(cell,p);
462 for(
int p=0;p<num_points;++p)
463 for(
int d=0;d<num_dim;++d)
464 cell_cub_points(p,d)=cub_points(cell,p,d);
469 intrepid_basis->getValues(cell_basis_ref_scalar.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_VALUE);
474 fst::HVOLtransformVALUE(cell_basis_scalar.get_view(),cell_jac_det.get_view(),cell_basis_ref_scalar.get_view());
475 for(
int b=0;b<num_basis;++b)
476 for(
int p=0;p<num_points;++p)
477 basis_scalar(cell,b,p)=cell_basis_scalar(0,b,p);
484 const std::pair<int,int> cell_range(0,num_cells);
485 auto s_weighted_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
486 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
487 auto s_basis_scalar = Kokkos::subview(basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
488 fst::multiplyMeasure(s_weighted_basis_scalar,s_weighted_measure,s_basis_scalar);
493 template <
typename Scalar>
496 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
497 const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
498 const int in_num_cells)
503 typedef Intrepid2::FunctionSpaceTools<PHX::Device::execution_space> fst;
508 const int num_points = basis_layout->numPoints();
510 const int num_dim = basis_layout->dimension();
511 const int num_cells = in_num_cells < 0 ? cub_points.extent(0) : in_num_cells;
514 auto cell_cub_points = af.
buildStaticArray<Scalar,IP,
Dim>(
"cell_cub_points",num_points,num_dim);
515 auto cell_grad_basis = af.
buildStaticArray<Scalar,Cell,BASIS,IP,Dim>(
"cell_grad_basis",1,num_basis,num_points,num_dim);
516 auto cell_jac_inv = af.
buildStaticArray<Scalar,Cell,IP,Dim,Dim>(
"cell_jac_inv",1,num_points,num_dim,num_dim);
518 auto cell_basis_ref_scalar = af.
buildStaticArray<Scalar,BASIS,IP>(
"cell_basis_ref_scalar",num_basis,num_points);
519 auto cell_grad_basis_ref = af.
buildStaticArray<Scalar,BASIS,IP,Dim>(
"cell_grad_basis_ref",num_basis,num_points,num_dim);
521 for(
int cell=0;cell<num_cells;++cell){
526 for(
int p=0;p<num_points;++p)
527 for(
int d=0;d<num_dim;++d)
528 for(
int d2=0;d2<num_dim;++d2)
529 cell_jac_inv(0,p,d,d2)=jac_inv(cell,p,d,d2);
530 for(
int p=0;p<num_points;++p)
531 for(
int d=0;d<num_dim;++d)
532 cell_cub_points(p,d)=cub_points(cell,p,d);
537 intrepid_basis->getValues(cell_basis_ref_scalar.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_VALUE);
539 if(compute_derivatives){
540 intrepid_basis->getValues(cell_grad_basis_ref.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_GRAD);
546 fst::HGRADtransformVALUE(cell_basis_scalar.get_view(),cell_basis_ref_scalar.get_view());
547 for(
int b=0;b<num_basis;++b)
548 for(
int p=0;p<num_points;++p)
549 basis_scalar(cell,b,p)=cell_basis_scalar(0,b,p);
551 if(compute_derivatives){
552 fst::HGRADtransformGRAD(cell_grad_basis.get_view(),cell_jac_inv.get_view(),cell_grad_basis_ref.get_view());
553 for(
int b=0;b<num_basis;++b)
554 for(
int p=0;p<num_points;++p)
555 for(
int d=0;d<num_dim;++d)
556 grad_basis(cell,b,p,d)=cell_grad_basis(0,b,p,d);
562 const std::pair<int,int> cell_range(0,num_cells);
563 auto s_weighted_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
564 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
565 auto s_basis_scalar = Kokkos::subview(basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
566 fst::multiplyMeasure(s_weighted_basis_scalar,s_weighted_measure,s_basis_scalar);
567 if(compute_derivatives){
568 auto s_weighted_grad_basis = Kokkos::subview(weighted_grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
569 auto s_grad_basis = Kokkos::subview(grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
570 fst::multiplyMeasure(s_weighted_grad_basis,s_weighted_measure,s_grad_basis);
577 template <
typename Scalar>
580 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> &
jac,
581 const PHX::MDField<Scalar,Cell,IP> & jac_det,
582 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
583 const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
584 const int in_num_cells)
590 typedef Intrepid2::FunctionSpaceTools<PHX::Device::execution_space> fst;
595 const int num_points = basis_layout->numPoints();
597 const int num_dim = basis_layout->dimension();
598 const int num_cells = in_num_cells < 0 ? basis_layout->numCells() : in_num_cells;
601 auto cell_jac = af.
buildStaticArray<Scalar,
Cell,IP,Dim,Dim>(
"cell_jac",1,num_points,num_dim,num_dim);
602 auto cell_jac_inv = af.
buildStaticArray<Scalar,Cell,IP,Dim,Dim>(
"cell_jac_inv",1,num_points,num_dim,num_dim);
603 auto cell_jac_det = af.
buildStaticArray<Scalar,Cell,IP>(
"cell_jac_det",1,num_points);
605 auto cell_basis_vector = af.
buildStaticArray<Scalar,Cell,
BASIS,IP,Dim>(
"cell_basis_vector",1,num_basis,num_points,num_dim);
606 auto cell_curl_basis_scalar = af.
buildStaticArray<Scalar,Cell,BASIS,IP>(
"cell_curl_basis_scalar",1,num_basis,num_points);
607 auto cell_curl_basis_vector = af.
buildStaticArray<Scalar,Cell,BASIS,IP,Dim>(
"cell_curl_basis_vector",1,num_basis,num_points,num_dim);
609 auto cell_curl_basis_ref = af.
buildArray<Scalar,BASIS,IP,Dim>(
"cell_curl_basis_ref",num_basis,num_points,num_dim);
610 auto cell_curl_basis_ref_scalar = af.
buildStaticArray<Scalar,BASIS,IP>(
"cell_curl_basis_ref_scalar",num_basis,num_points);
611 auto cell_basis_ref_vector = af.
buildArray<Scalar,BASIS,IP,Dim>(
"cell_basis_ref_vector",num_basis,num_points,num_dim);
613 for(
int cell=0;cell<num_cells;++cell){
618 for(
int p=0;p<num_points;++p)
619 for(
int d=0;d<num_dim;++d)
620 for(
int d2=0;d2<num_dim;++d2)
621 cell_jac(0,p,d,d2)=
jac(cell,p,d,d2);
622 for(
int p=0;p<num_points;++p)
623 for(
int d=0;d<num_dim;++d)
624 for(
int d2=0;d2<num_dim;++d2)
625 cell_jac_inv(0,p,d,d2)=jac_inv(cell,p,d,d2);
626 for(
int p=0;p<num_points;++p)
627 cell_jac_det(0,p)=jac_det(cell,p);
628 for(
int p=0;p<num_points;++p)
629 for(
int d=0;d<num_dim;++d)
630 cell_cub_points(p,d)=cub_points(cell,p,d);
635 intrepid_basis->getValues(cell_basis_ref_vector.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_VALUE);
637 if(compute_derivatives){
639 intrepid_basis->getValues(cell_curl_basis_ref_scalar.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_CURL);
640 }
else if(num_dim==3){
641 intrepid_basis->getValues(cell_curl_basis_ref.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_CURL);
648 fst::HCURLtransformVALUE(cell_basis_vector.get_view(),cell_jac_inv.get_view(),cell_basis_ref_vector.get_view());
649 for(
int b=0;b<num_basis;++b)
650 for(
int p=0;p<num_points;++p)
651 for(
int d=0;d<num_dim;++d)
652 basis_vector(cell,b,p,d)=cell_basis_vector(0,b,p,d);
654 if(compute_derivatives){
659 fst::HDIVtransformDIV(cell_curl_basis_scalar.get_view(),cell_jac_det.get_view(),cell_curl_basis_ref_scalar.get_view());
660 for(
int b=0;b<num_basis;++b)
661 for(
int p=0;p<num_points;++p)
662 curl_basis_scalar(cell,b,p)=cell_curl_basis_scalar(0,b,p);
663 }
else if(num_dim==3) {
664 fst::HCURLtransformCURL(cell_curl_basis_vector.get_view(),cell_jac.get_view(),cell_jac_det.get_view(),cell_curl_basis_ref.get_view());
665 for(
int b=0;b<num_basis;++b)
666 for(
int p=0;p<num_points;++p)
667 for(
int d=0;d<num_dim;++d)
668 curl_basis_vector(cell,b,p,d)=cell_curl_basis_vector(0,b,p,d);
676 const std::pair<int,int> cell_range(0,num_cells);
677 auto s_weighted_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
678 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
679 auto s_basis_vector = Kokkos::subview(basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
680 fst::multiplyMeasure(s_weighted_basis_vector,s_weighted_measure,s_basis_vector);
681 if(compute_derivatives){
683 auto s_weighted_curl_basis_scalar = Kokkos::subview(weighted_curl_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
684 auto s_curl_basis_scalar = Kokkos::subview(curl_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
685 fst::multiplyMeasure(s_weighted_curl_basis_scalar,s_weighted_measure,s_curl_basis_scalar);
686 }
else if(num_dim==3){
687 auto s_weighted_curl_basis_vector = Kokkos::subview(weighted_curl_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
688 auto s_curl_basis_vector = Kokkos::subview(curl_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
689 fst::multiplyMeasure(s_weighted_curl_basis_vector,s_weighted_measure,s_curl_basis_vector);
696 template <
typename Scalar>
699 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> &
jac,
700 const PHX::MDField<Scalar,Cell,IP> & jac_det,
701 const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
702 const int in_num_cells)
707 typedef Intrepid2::FunctionSpaceTools<PHX::Device::execution_space> fst;
712 const int num_points = basis_layout->numPoints();
714 const int num_dim = basis_layout->dimension();
715 const int num_cells = in_num_cells < 0 ? basis_layout->numCells() : in_num_cells;
718 auto cell_jac = af.
buildStaticArray<Scalar,
Cell,IP,Dim,Dim>(
"cell_jac",1,num_points,num_dim,num_dim);
719 auto cell_jac_det = af.
buildStaticArray<Scalar,Cell,IP>(
"cell_jac_det",1,num_points);
721 auto cell_basis_vector = af.
buildStaticArray<Scalar,Cell,
BASIS,IP,Dim>(
"cell_basis_vector",1,num_basis,num_points,num_dim);
722 auto cell_div_basis = af.
buildStaticArray<Scalar,Cell,BASIS,IP>(
"cell_div_basis",1,num_basis,num_points);
724 auto cell_basis_ref_vector = af.
buildArray<Scalar,BASIS,IP,Dim>(
"cell_basis_ref_vector",num_basis,num_points,num_dim);
725 auto cell_div_basis_ref = af.
buildStaticArray<Scalar,BASIS,IP>(
"cell_div_basis_ref",num_basis,num_points);
727 for(
int cell=0;cell<num_cells;++cell){
732 for(
int p=0;p<num_points;++p)
733 for(
int d=0;d<num_dim;++d)
734 for(
int d2=0;d2<num_dim;++d2)
735 cell_jac(0,p,d,d2)=
jac(cell,p,d,d2);
736 for(
int p=0;p<num_points;++p)
737 cell_jac_det(0,p)=jac_det(cell,p);
738 for(
int p=0;p<num_points;++p)
739 for(
int d=0;d<num_dim;++d)
740 cell_cub_points(p,d)=cub_points(cell,p,d);
744 intrepid_basis->getValues(cell_basis_ref_vector.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_VALUE);
746 if(compute_derivatives){
747 intrepid_basis->getValues(cell_div_basis_ref.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_DIV);
753 fst::HDIVtransformVALUE(cell_basis_vector.get_view(),cell_jac.get_view(),cell_jac_det.get_view(),cell_basis_ref_vector.get_view());
754 for(
int b=0;b<num_basis;++b)
755 for(
int p=0;p<num_points;++p)
756 for(
int d=0;d<num_dim;++d)
757 basis_vector(cell,b,p,d)=cell_basis_vector(0,b,p,d);
759 if(compute_derivatives){
760 fst::HDIVtransformDIV(cell_div_basis.get_view(),cell_jac_det.get_view(),cell_div_basis_ref.get_view());
761 for(
int b=0;b<num_basis;++b)
762 for(
int p=0;p<num_points;++p)
763 div_basis(cell,b,p)=cell_div_basis(0,b,p);
768 const std::pair<int,int> cell_range(0,num_cells);
769 auto s_weighted_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
770 auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
771 auto s_basis_vector = Kokkos::subview(basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
772 fst::multiplyMeasure(s_weighted_basis_vector,s_weighted_measure,s_basis_vector);
773 if(compute_derivatives){
774 auto s_weighted_div_basis = Kokkos::subview(weighted_div_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
775 auto s_div_basis = Kokkos::subview(div_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
776 fst::multiplyMeasure(s_weighted_div_basis,s_weighted_measure,s_div_basis);
782 template <
typename Scalar>
785 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> &
jac,
786 const PHX::MDField<Scalar,Cell,IP> & jac_det,
787 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv)
790 PHX::MDField<Scalar,Cell,NODE,Dim> vertex_coordinates;
791 const int in_num_cells = jac.extent(0);
792 evaluateValuesCV(cub_points,jac,jac_det,jac_inv,vertex_coordinates,
false,in_num_cells);
796 template <
typename Scalar>
799 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> &
jac,
800 const PHX::MDField<Scalar,Cell,IP> & jac_det,
801 const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
802 const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
803 bool use_vertex_coordinates,
804 const int in_num_cells)
808 int num_ip = basis_layout->numPoints();
809 int num_card = basis_layout->cardinality();
810 int num_dim = basis_layout->dimension();
819 for (
size_type icell = 0; icell < num_cells; ++icell)
821 for (
int ip = 0; ip < num_ip; ++ip)
822 for (
int d = 0; d < num_dim; ++d)
823 dyn_cub_points(ip,d) = cell_cub_points(icell,ip,d);
828 intrepid_basis->getValues(dyn_basis_ref_scalar.get_view(),
829 dyn_cub_points.get_view(),
830 Intrepid2::OPERATOR_VALUE);
833 for (
int b = 0; b < num_card; ++b)
834 for (
int ip = 0; ip < num_ip; ++ip)
835 basis_scalar(icell,b,ip) = dyn_basis_ref_scalar(b,ip);
841 intrepid_basis->getValues(dyn_basis_ref_scalar.get_view(),
842 dyn_cub_points.get_view(),
843 Intrepid2::OPERATOR_VALUE);
851 for (
int ip = 0; ip < num_ip; ++ip)
852 dyn_jac_det(cellInd,ip) = jac_det(icell,ip);
854 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HVOLtransformVALUE(dyn_basis_scalar.get_view(),
855 dyn_jac_det.get_view(),
856 dyn_basis_ref_scalar.get_view());
858 for (
int b = 0; b < num_card; ++b)
859 for (
int ip = 0; ip < num_ip; ++ip)
860 basis_scalar(icell,b,ip) = dyn_basis_scalar(0,b,ip);
866 intrepid_basis->getValues(dyn_basis_ref_scalar.get_view(),
867 dyn_cub_points.get_view(),
868 Intrepid2::OPERATOR_VALUE);
871 for (
int b = 0; b < num_card; ++b)
872 for (
int ip = 0; ip < num_ip; ++ip)
873 basis_scalar(icell,b,ip) = dyn_basis_ref_scalar(b,ip);
875 if(compute_derivatives) {
878 ArrayDynamic dyn_grad_basis_ref = af.
buildArray<Scalar,BASIS,IP,Dim>(
"dyn_grad_basis_ref",num_card,num_ip,num_dim);
880 ArrayDynamic dyn_jac_inv = af.
buildArray<Scalar,Cell,IP,Dim,Dim>(
"dyn_jac_inv",one_cell,num_ip,num_dim,num_dim);
882 intrepid_basis->getValues(dyn_grad_basis_ref.get_view(),
883 dyn_cub_points.get_view(),
884 Intrepid2::OPERATOR_GRAD);
887 for (
int ip = 0; ip < num_ip; ++ip)
888 for (
int d1 = 0; d1 < num_dim; ++d1)
889 for (
int d2 = 0; d2 < num_dim; ++d2)
890 dyn_jac_inv(cellInd,ip,d1,d2) = jac_inv(icell,ip,d1,d2);
892 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HGRADtransformGRAD<Scalar>(dyn_grad_basis.get_view(),
893 dyn_jac_inv.get_view(),
894 dyn_grad_basis_ref.get_view());
896 for (
int b = 0; b < num_card; ++b)
897 for (
int ip = 0; ip < num_ip; ++ip)
898 for (
int d = 0; d < num_dim; ++d)
899 grad_basis(icell,b,ip,d) = dyn_grad_basis(0,b,ip,d);
906 intrepid_basis->getValues(dyn_basis_ref_vector.get_view(),
907 dyn_cub_points.get_view(),
908 Intrepid2::OPERATOR_VALUE);
910 const int one_cell = 1;
912 ArrayDynamic dyn_jac_inv = af.
buildArray<Scalar,Cell,IP,Dim,Dim>(
"dyn_jac_inv",one_cell,num_ip,num_dim,num_dim);
914 const int cellInd = 0;
915 for (
int ip = 0; ip < num_ip; ++ip)
916 for (
int d1 = 0; d1 < num_dim; ++d1)
917 for (
int d2 = 0; d2 < num_dim; ++d2)
918 dyn_jac_inv(cellInd,ip,d1,d2) = jac_inv(icell,ip,d1,d2);
920 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HCURLtransformVALUE(dyn_basis_vector.get_view(),
921 dyn_jac_inv.get_view(),
922 dyn_basis_ref_vector.get_view());
924 for (
int b = 0; b < num_card; ++b)
925 for (
int ip = 0; ip < num_ip; ++ip)
926 for (
int d = 0; d < num_dim; ++d)
927 basis_vector(icell,b,ip,d) = dyn_basis_vector(0,b,ip,d);
929 if(compute_derivatives && num_dim ==2) {
931 ArrayDynamic dyn_curl_basis_ref_scalar = af.
buildArray<Scalar,BASIS,IP>(
"dyn_curl_basis_ref_scalar",num_card,num_ip);
932 ArrayDynamic dyn_curl_basis_scalar = af.
buildArray<Scalar,Cell,BASIS,IP>(
"dyn_curl_basis_scalar",one_cell,num_card,num_ip);
935 intrepid_basis->getValues(dyn_curl_basis_ref_scalar.get_view(),
936 dyn_cub_points.get_view(),
937 Intrepid2::OPERATOR_CURL);
939 for (
int ip = 0; ip < num_ip; ++ip)
940 dyn_jac_det(cellInd,ip) = jac_det(icell,ip);
942 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HDIVtransformDIV(dyn_curl_basis_scalar.get_view(),
943 dyn_jac_det.get_view(),
944 dyn_curl_basis_ref_scalar.get_view());
946 for (
int b = 0; b < num_card; ++b)
947 for (
int ip = 0; ip < num_ip; ++ip)
948 curl_basis_scalar(icell,b,ip) = dyn_curl_basis_scalar(0,b,ip);
951 if(compute_derivatives && num_dim ==3) {
953 ArrayDynamic dyn_curl_basis_ref = af.
buildArray<Scalar,BASIS,IP,Dim>(
"dyn_curl_basis_ref_vector",num_card,num_ip,num_dim);
954 ArrayDynamic dyn_curl_basis = af.
buildArray<Scalar,Cell,BASIS,IP,Dim>(
"dyn_curl_basis_vector",one_cell,num_card,num_ip,num_dim);
958 intrepid_basis->getValues(dyn_curl_basis_ref.get_view(),
959 dyn_cub_points.get_view(),
960 Intrepid2::OPERATOR_CURL);
962 for (
int ip = 0; ip < num_ip; ++ip)
964 dyn_jac_det(cellInd,ip) = jac_det(icell,ip);
965 for (
int d1 = 0; d1 < num_dim; ++d1)
966 for (
int d2 = 0; d2 < num_dim; ++d2)
967 dyn_jac(cellInd,ip,d1,d2) =
jac(icell,ip,d1,d2);
970 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HCURLtransformCURL(dyn_curl_basis.get_view(),
972 dyn_jac_det.get_view(),
973 dyn_curl_basis_ref.get_view());
975 for (
int b = 0; b < num_card; ++b)
976 for (
int ip = 0; ip < num_ip; ++ip)
977 for (
int d = 0; d < num_dim; ++d)
978 curl_basis_vector(icell,b,ip,d) = dyn_curl_basis(0,b,ip,d);
987 intrepid_basis->getValues(dyn_basis_ref_vector.get_view(),
988 dyn_cub_points.get_view(),
989 Intrepid2::OPERATOR_VALUE);
997 for (
int ip = 0; ip < num_ip; ++ip)
999 dyn_jac_det(cellInd,ip) = jac_det(icell,ip);
1000 for (
int d1 = 0; d1 < num_dim; ++d1)
1001 for (
int d2 = 0; d2 < num_dim; ++d2)
1002 dyn_jac(cellInd,ip,d1,d2) =
jac(icell,ip,d1,d2);
1005 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HDIVtransformVALUE(dyn_basis_vector.get_view(),
1007 dyn_jac_det.get_view(),
1008 dyn_basis_ref_vector.get_view());
1010 for (
int b = 0; b < num_card; ++b)
1011 for (
int ip = 0; ip < num_ip; ++ip)
1012 for (
int d = 0; d < num_dim; ++d)
1013 basis_vector(icell,b,ip,d) = dyn_basis_vector(0,b,ip,d);
1015 if(compute_derivatives) {
1017 ArrayDynamic dyn_div_basis_ref = af.
buildArray<Scalar,BASIS,IP>(
"dyn_div_basis_ref_scalar",num_card,num_ip);
1018 ArrayDynamic dyn_div_basis = af.
buildArray<Scalar,Cell,BASIS,IP>(
"dyn_div_basis_scalar",one_cell,num_card,num_ip);
1020 intrepid_basis->getValues(dyn_div_basis_ref.get_view(),
1021 dyn_cub_points.get_view(),
1022 Intrepid2::OPERATOR_DIV);
1024 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HDIVtransformDIV<Scalar>(dyn_div_basis.get_view(),
1025 dyn_jac_det.get_view(),
1026 dyn_div_basis_ref.get_view());
1028 for (
int b = 0; b < num_card; ++b)
1029 for (
int ip = 0; ip < num_ip; ++ip)
1030 div_basis(icell,b,ip) = dyn_div_basis(0,b,ip);
1039 if(use_vertex_coordinates) {
1041 evaluateBasisCoordinates(vertex_coordinates);
1046 template <
typename Scalar>
1052 int num_quad = basis_layout->numPoints();
1053 int num_dim = basis_layout->dimension();
1054 int num_card = basis_layout->cardinality();
1058 for (
int ip = 0; ip < num_quad; ++ip)
1059 for (
int d = 0; d < num_dim; ++d)
1060 dyn_cub_points(ip,d) = cub_points(ip,d);
1066 intrepid_basis->getValues(dyn_basis_ref_scalar.get_view(),
1067 dyn_cub_points.get_view(),
1068 Intrepid2::OPERATOR_VALUE);
1070 for (
int b = 0; b < num_card; ++b)
1071 for (
int ip = 0; ip < num_quad; ++ip)
1072 basis_ref_scalar(b,ip) = dyn_basis_ref_scalar(b,ip);
1077 intrepid_basis->getValues(dyn_basis_ref_vector.get_view(),
1078 dyn_cub_points.get_view(),
1079 Intrepid2::OPERATOR_VALUE);
1081 for (
int b = 0; b < num_card; ++b)
1082 for (
int ip = 0; ip < num_quad; ++ip)
1083 for (
int d = 0; d < num_dim; ++d)
1084 basis_ref_vector(b,ip,d) = dyn_basis_ref_vector(b,ip,d);
1091 intrepid_basis->getValues(dyn_grad_basis_ref.get_view(),
1092 dyn_cub_points.get_view(),
1093 Intrepid2::OPERATOR_GRAD);
1095 for (
int b = 0; b < num_card; ++b)
1096 for (
int ip = 0; ip < num_quad; ++ip)
1097 for (
int d = 0; d < num_dim; ++d)
1098 grad_basis_ref(b,ip,d) = dyn_grad_basis_ref(b,ip,d);
1100 else if(elmtspace==
PureBasis::HCURL && in_compute_derivatives && num_dim==2) {
1103 intrepid_basis->getValues(dyn_curl_basis_ref.get_view(),
1104 dyn_cub_points.get_view(),
1105 Intrepid2::OPERATOR_CURL);
1107 for (
int b = 0; b < num_card; ++b)
1108 for (
int ip = 0; ip < num_quad; ++ip)
1109 curl_basis_ref_scalar(b,ip) = dyn_curl_basis_ref(b,ip);
1111 else if(elmtspace==
PureBasis::HCURL && in_compute_derivatives && num_dim==3) {
1114 intrepid_basis->getValues(dyn_curl_basis_ref.get_view(),
1115 dyn_cub_points.get_view(),
1116 Intrepid2::OPERATOR_CURL);
1118 for (
int b = 0; b < num_card; ++b)
1119 for (
int ip = 0; ip < num_quad; ++ip)
1120 for (
int d = 0; d < num_dim; ++d)
1121 curl_basis_ref_vector(b,ip,d) = dyn_curl_basis_ref(b,ip,d);
1126 intrepid_basis->getValues(dyn_div_basis_ref.get_view(),
1127 dyn_cub_points.get_view(),
1128 Intrepid2::OPERATOR_DIV);
1130 for (
int b = 0; b < num_card; ++b)
1131 for (
int ip = 0; ip < num_quad; ++ip)
1132 div_basis_ref(b,ip) = dyn_div_basis_ref(b,ip);
1136 if(use_vertex_coordinates) {
1141 using coordsScalarType =
typename Intrepid2::Basis<PHX::Device::execution_space,Scalar,Scalar>::scalarType;
1142 auto dyn_basis_coordinates_ref = af.
buildArray<coordsScalarType,
BASIS,Dim>(
"basis_coordinates_ref",
1143 basis_coordinates_ref.extent(0),
1144 basis_coordinates_ref.extent(1));
1145 intrepid_basis->getDofCoords(dyn_basis_coordinates_ref.get_view());
1148 for (
int i = 0; i < basis_coordinates_ref.extent_int(0); ++i)
1149 for (
int j = 0; j < basis_coordinates_ref.extent_int(1); ++j)
1150 basis_coordinates_ref(i,j) = dyn_basis_coordinates_ref(i,j);
1154 references_evaluated =
true;
1158 template <
typename Scalar>
1161 const int in_num_cells)
1163 if (!intrepid_basis->requireOrientation())
1166 typedef Intrepid2::OrientationTools<PHX::Device> ots;
1173 const int num_cell_basis_layout = in_num_cells < 0 ? basis_layout->numCells() : in_num_cells;
1174 const int num_cell_orientation = orientations.size();
1175 const int num_cell = num_cell_basis_layout < num_cell_orientation ? num_cell_basis_layout : num_cell_orientation;
1176 const int num_dim = basis_layout->dimension();
1177 const Kokkos::pair<int,int> range_cell(0, num_cell);
1181 Kokkos::DynRankView<Intrepid2::Orientation,PHX::Device> drv_orts(
"drv_orts", num_cell);
1182 auto host_drv_orts = Kokkos::create_mirror_view(drv_orts);
1183 for (
size_t i=0; i < drv_orts.size(); ++i)
1184 host_drv_orts(i) = orientations[i];
1185 Kokkos::deep_copy(drv_orts,host_drv_orts);
1186 typename PHX::Device().fence();
1194 auto drv_basis_scalar = Kokkos::subview(basis_scalar.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1196 "drv_basis_scalar_tmp",
1197 drv_basis_scalar.extent(0),
1198 drv_basis_scalar.extent(1),
1199 drv_basis_scalar.extent(2));
1200 Kokkos::deep_copy(drv_basis_scalar_tmp, drv_basis_scalar);
1201 ots::modifyBasisByOrientation(drv_basis_scalar,
1202 drv_basis_scalar_tmp,
1204 intrepid_basis.getRawPtr());
1206 if(build_weighted) {
1207 auto drv_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1209 "drv_basis_scalar_tmp",
1210 drv_basis_scalar.extent(0),
1211 drv_basis_scalar.extent(1),
1212 drv_basis_scalar.extent(2));
1213 Kokkos::deep_copy(drv_basis_scalar_tmp, drv_basis_scalar);
1214 ots::modifyBasisByOrientation(drv_basis_scalar,
1215 drv_basis_scalar_tmp,
1217 intrepid_basis.getRawPtr());
1222 if (compute_derivatives) {
1224 auto drv_grad_basis = Kokkos::subview(grad_basis.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1226 "drv_grad_basis_tmp",
1227 drv_grad_basis.extent(0),
1228 drv_grad_basis.extent(1),
1229 drv_grad_basis.extent(2),
1230 drv_grad_basis.extent(3));
1231 Kokkos::deep_copy(drv_grad_basis_tmp, drv_grad_basis);
1232 ots::modifyBasisByOrientation(drv_grad_basis,
1235 intrepid_basis.getRawPtr());
1237 if(build_weighted) {
1238 auto drv_grad_basis = Kokkos::subview(weighted_grad_basis.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1240 "drv_grad_basis_tmp",
1241 drv_grad_basis.extent(0),
1242 drv_grad_basis.extent(1),
1243 drv_grad_basis.extent(2),
1244 drv_grad_basis.extent(3));
1245 Kokkos::deep_copy(drv_grad_basis_tmp, drv_grad_basis);
1246 ots::modifyBasisByOrientation(drv_grad_basis,
1249 intrepid_basis.getRawPtr());
1260 auto drv_basis_vector = Kokkos::subview(basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1262 "drv_basis_vector_tmp",
1263 drv_basis_vector.extent(0),
1264 drv_basis_vector.extent(1),
1265 drv_basis_vector.extent(2),
1266 drv_basis_vector.extent(3));
1267 Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1268 ots::modifyBasisByOrientation(drv_basis_vector,
1269 drv_basis_vector_tmp,
1271 intrepid_basis.getRawPtr());
1273 if(build_weighted) {
1274 auto drv_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1276 "drv_basis_vector_tmp",
1277 drv_basis_vector.extent(0),
1278 drv_basis_vector.extent(1),
1279 drv_basis_vector.extent(2),
1280 drv_basis_vector.extent(3));
1281 Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1282 ots::modifyBasisByOrientation(drv_basis_vector,
1283 drv_basis_vector_tmp,
1285 intrepid_basis.getRawPtr());
1289 if (compute_derivatives) {
1291 auto drv_curl_basis_scalar = Kokkos::subview(curl_basis_scalar.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1293 "drv_curl_basis_scalar_tmp",
1294 drv_curl_basis_scalar.extent(0),
1295 drv_curl_basis_scalar.extent(1),
1296 drv_curl_basis_scalar.extent(2));
1297 Kokkos::deep_copy(drv_curl_basis_scalar_tmp, drv_curl_basis_scalar);
1298 ots::modifyBasisByOrientation(drv_curl_basis_scalar,
1299 drv_curl_basis_scalar_tmp,
1301 intrepid_basis.getRawPtr());
1304 if(build_weighted) {
1305 auto drv_curl_basis_scalar = Kokkos::subview(weighted_curl_basis_scalar.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1307 "drv_curl_basis_scalar_tmp",
1308 drv_curl_basis_scalar.extent(0),
1309 drv_curl_basis_scalar.extent(1),
1310 drv_curl_basis_scalar.extent(2));
1311 Kokkos::deep_copy(drv_curl_basis_scalar_tmp, drv_curl_basis_scalar);
1312 ots::modifyBasisByOrientation(drv_curl_basis_scalar,
1313 drv_curl_basis_scalar_tmp,
1315 intrepid_basis.getRawPtr());
1326 auto drv_basis_vector = Kokkos::subview(basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1328 "drv_basis_vector_tmp",
1329 drv_basis_vector.extent(0),
1330 drv_basis_vector.extent(1),
1331 drv_basis_vector.extent(2),
1332 drv_basis_vector.extent(3));
1333 Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1334 ots::modifyBasisByOrientation(drv_basis_vector,
1335 drv_basis_vector_tmp,
1337 intrepid_basis.getRawPtr());
1339 if(build_weighted) {
1340 auto drv_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1342 "drv_basis_vector_tmp",
1343 drv_basis_vector.extent(0),
1344 drv_basis_vector.extent(1),
1345 drv_basis_vector.extent(2),
1346 drv_basis_vector.extent(3));
1347 Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1348 ots::modifyBasisByOrientation(drv_basis_vector,
1349 drv_basis_vector_tmp,
1351 intrepid_basis.getRawPtr());
1355 if (compute_derivatives) {
1357 auto drv_curl_basis_vector = Kokkos::subview(curl_basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1359 "drv_curl_basis_vector_tmp",
1360 drv_curl_basis_vector.extent(0),
1361 drv_curl_basis_vector.extent(1),
1362 drv_curl_basis_vector.extent(2),
1363 drv_curl_basis_vector.extent(3));
1364 Kokkos::deep_copy(drv_curl_basis_vector_tmp, drv_curl_basis_vector);
1365 ots::modifyBasisByOrientation(drv_curl_basis_vector,
1366 drv_curl_basis_vector_tmp,
1368 intrepid_basis.getRawPtr());
1370 if(build_weighted) {
1371 auto drv_curl_basis_vector = Kokkos::subview(weighted_curl_basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1373 "drv_curl_basis_vector_tmp",
1374 drv_curl_basis_vector.extent(0),
1375 drv_curl_basis_vector.extent(1),
1376 drv_curl_basis_vector.extent(2),
1377 drv_curl_basis_vector.extent(3));
1378 Kokkos::deep_copy(drv_curl_basis_vector_tmp, drv_curl_basis_vector);
1379 ots::modifyBasisByOrientation(drv_curl_basis_vector,
1380 drv_curl_basis_vector_tmp,
1382 intrepid_basis.getRawPtr());
1392 auto drv_basis_vector = Kokkos::subview(basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1394 "drv_basis_vector_tmp",
1395 drv_basis_vector.extent(0),
1396 drv_basis_vector.extent(1),
1397 drv_basis_vector.extent(2),
1398 drv_basis_vector.extent(3));
1399 Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1400 ots::modifyBasisByOrientation(drv_basis_vector,
1401 drv_basis_vector_tmp,
1403 intrepid_basis.getRawPtr());
1405 if(build_weighted) {
1406 auto drv_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1408 "drv_basis_vector_tmp",
1409 drv_basis_vector.extent(0),
1410 drv_basis_vector.extent(1),
1411 drv_basis_vector.extent(2),
1412 drv_basis_vector.extent(3));
1413 Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1414 ots::modifyBasisByOrientation(drv_basis_vector,
1415 drv_basis_vector_tmp,
1417 intrepid_basis.getRawPtr());
1420 if (compute_derivatives) {
1422 auto drv_div_basis = Kokkos::subview(div_basis.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1424 "drv_div_basis_tmp",
1425 drv_div_basis.extent(0),
1426 drv_div_basis.extent(1),
1427 drv_div_basis.extent(2));
1428 Kokkos::deep_copy(drv_div_basis_tmp, drv_div_basis);
1429 ots::modifyBasisByOrientation(drv_div_basis,
1432 intrepid_basis.getRawPtr());
1434 if(build_weighted) {
1435 auto drv_div_basis = Kokkos::subview(weighted_div_basis.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1437 "drv_div_basis_tmp",
1438 drv_div_basis.extent(0),
1439 drv_div_basis.extent(1),
1440 drv_div_basis.extent(2));
1441 Kokkos::deep_copy(drv_div_basis_tmp, drv_div_basis);
1442 ots::modifyBasisByOrientation(drv_div_basis,
1445 intrepid_basis.getRawPtr());
1452 template <
typename Scalar>
1459 int num_cell = orientations.extent(0);
1460 int num_basis = orientations.extent(1);
1461 int num_dim = basis_layout->dimension();
1462 int num_ip = basis_layout->numPoints();
1470 for (
int c=0; c<num_cell; c++)
1471 for (
int b=0; b<num_basis; b++)
1472 for (
int p=0; p<num_ip; p++)
1473 for (
int d=0; d<num_dim; d++)
1474 basis_vector(c, b, p, d) *= orientations(c, b);
1476 if(compute_derivatives) {
1478 for (
int c=0; c<num_cell; c++)
1479 for (
int b=0; b<num_basis; b++)
1480 for (
int p=0; p<num_ip; p++)
1481 curl_basis_scalar(c, b, p) *= orientations(c, b);
1485 if(build_weighted) {
1486 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_basis_vector.get_view(),orientations.get_view());
1488 if(compute_derivatives)
1489 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_curl_basis_scalar.get_view(),orientations.get_view());
1497 for (
int c=0; c<num_cell; c++)
1498 for (
int b=0; b<num_basis; b++)
1499 for (
int p=0; p<num_ip; p++)
1500 for (
int d=0; d<num_dim; d++)
1501 basis_vector(c, b, p, d) *= orientations(c, b);
1503 if(compute_derivatives) {
1505 for (
int c=0; c<num_cell; c++)
1506 for (
int b=0; b<num_basis; b++)
1507 for (
int p=0; p<num_ip; p++)
1508 for (
int d=0; d<num_dim; d++)
1509 curl_basis_vector(c, b, p,d) *= orientations(c, b);
1513 if(build_weighted) {
1514 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_basis_vector.get_view(),orientations.get_view());
1516 if(compute_derivatives)
1517 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_curl_basis_vector.get_view(),orientations.get_view());
1524 for (
int c=0; c<num_cell; c++)
1525 for (
int b=0; b<num_basis; b++)
1526 for (
int p=0; p<num_ip; p++)
1527 for (
int d=0; d<num_dim; d++)
1528 basis_vector(c, b, p, d) *= orientations(c, b);
1530 if(compute_derivatives) {
1533 for (
int c=0; c<num_cell; c++)
1534 for (
int b=0; b<num_basis; b++)
1535 for (
int p=0; p<num_ip; p++)
1536 div_basis(c, b, p) *= orientations(c, b);
1540 if(build_weighted) {
1541 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_basis_vector.get_view(),orientations.get_view());
1543 if(compute_derivatives)
1544 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_div_basis.get_view(),orientations.get_view());
1549 template <
typename Scalar>
1551 {
return basis_layout->getBasis()->getElementSpace(); }
1553 template <
typename Scalar>
1556 bool computeDerivatives)
1560 compute_derivatives = computeDerivatives;
1561 basis_layout = layout;
1568 int numcells = basisDesc->
numCells();
1572 intrepid_basis = basisDesc->
getIntrepid2Basis<PHX::Device::execution_space,Scalar,Scalar>();
1587 weighted_basis_scalar = af.
buildStaticArray<Scalar,Cell,BASIS,IP>(
"weighted_basis",numcells,card,num_quad);
1592 if(compute_derivatives) {
1593 grad_basis_ref = af.
buildStaticArray<Scalar,BASIS,IP,
Dim>(
"grad_basis_ref",card,num_quad,dim);
1594 grad_basis = af.
buildStaticArray<Scalar,Cell,BASIS,IP,Dim>(
"grad_basis",numcells,card,num_quad,dim);
1597 weighted_grad_basis = af.
buildStaticArray<Scalar,Cell,BASIS,IP,Dim>(
"weighted_grad_basis",numcells,card,num_quad,dim);
1612 basis_vector = af.
buildStaticArray<Scalar,
Cell,BASIS,IP,Dim>(
"basis",numcells,card,num_quad,dim);
1615 weighted_basis_vector = af.
buildStaticArray<Scalar,Cell,BASIS,IP,Dim>(
"weighted_basis",numcells,card,num_quad,dim);
1625 if(compute_derivatives) {
1628 curl_basis_ref_scalar = af.
buildStaticArray<Scalar,BASIS,IP>(
"curl_basis_ref",card,num_quad);
1629 curl_basis_scalar = af.
buildStaticArray<Scalar,Cell,BASIS,IP>(
"curl_basis",numcells,card,num_quad);
1632 weighted_curl_basis_scalar = af.
buildStaticArray<Scalar,Cell,BASIS,IP>(
"weighted_curl_basis",numcells,card,num_quad);
1635 curl_basis_ref_vector = af.
buildStaticArray<Scalar,BASIS,IP,Dim>(
"curl_basis_ref",card,num_quad,dim);
1636 curl_basis_vector = af.
buildStaticArray<Scalar,Cell,BASIS,IP,Dim>(
"curl_basis",numcells,card,num_quad,dim);
1639 weighted_curl_basis_vector = af.
buildStaticArray<Scalar,Cell,BASIS,IP,Dim>(
"weighted_curl_basis",numcells,card,num_quad,dim);
1651 basis_vector = af.
buildStaticArray<Scalar,
Cell,BASIS,IP,Dim>(
"basis",numcells,card,num_quad,dim);
1654 weighted_basis_vector = af.
buildStaticArray<Scalar,Cell,BASIS,IP,Dim>(
"weighted_basis",numcells,card,num_quad,dim);
1669 if(compute_derivatives) {
1670 div_basis_ref = af.
buildStaticArray<Scalar,BASIS,IP>(
"div_basis_ref",card,num_quad);
1671 div_basis = af.
buildStaticArray<Scalar,Cell,BASIS,IP>(
"div_basis",numcells,card,num_quad);
1674 weighted_div_basis = af.
buildStaticArray<Scalar,Cell,BASIS,IP>(
"weighted_div_basis",numcells,card,num_quad);
1686 weighted_basis_scalar = af.
buildStaticArray<Scalar,Cell,BASIS,IP>(
"weighted_basis",numcells,card,num_quad);
1706 basis_coordinates = af.
buildStaticArray<Scalar,
Cell,BASIS,Dim>(
"basis_coordinates",numcells,card,dim);
void evaluateValues_HCurl(const PHX::MDField< Scalar, Cell, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac, const PHX::MDField< Scalar, Cell, IP > &jac_det, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv, const PHX::MDField< Scalar, Cell, IP > &weighted_measure, const int in_num_cells)
Kokkos::DynRankView< typename InputArray::value_type, PHX::Device > createDynRankView(const InputArray &a, const std::string &name, const DimensionPack...dims)
Wrapper to simplify Panzer use of Sacado ViewFactory.
void evaluateValues_HDiv(const PHX::MDField< Scalar, Cell, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac, const PHX::MDField< Scalar, Cell, IP > &jac_det, const PHX::MDField< Scalar, Cell, IP > &weighted_measure, const int in_num_cells)
void evaluateValuesCV(const PHX::MDField< Scalar, Cell, IP, Dim > &cell_cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac, const PHX::MDField< Scalar, Cell, IP > &jac_det, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv)
void evaluateReferenceValues(const PHX::MDField< Scalar, IP, Dim > &cub_points, bool compute_derivatives, bool use_vertex_coordinates)
EElementSpace getElementSpace() const
void evaluateValues_HGrad(const PHX::MDField< Scalar, Cell, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv, const PHX::MDField< Scalar, Cell, IP > &weighted_measure, const int in_num_cells)
int cardinality() const
Returns the number of basis coefficients.
PHX::MDField< Scalar, T0 > buildStaticArray(const std::string &str, int d0) const
Teuchos::RCP< const PureBasis > getBasis() const
int dimension() const
Returns the dimension of the basis from the topology.
void applyOrientations(const PHX::MDField< const Scalar, Cell, BASIS > &orientations)
Method to apply orientations to a basis values container.
Teuchos::RCP< const shards::CellTopology > getCellTopology() const
#define TEUCHOS_TEST_FOR_EXCEPT_MSG(throw_exception_test, msg)
void evaluateBasisCoordinates(const PHX::MDField< Scalar, Cell, NODE, Dim > &vertex_coordinates, const int in_num_cells=-1)
ArrayTraits< Scalar, PHX::MDField< Scalar > >::size_type size_type
PureBasis::EElementSpace getElementSpace() const
Teuchos::RCP< Intrepid2::Basis< PHX::Device::execution_space, double, double > > getIntrepid2Basis() const
void evaluateValues_Const(const PHX::MDField< Scalar, Cell, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv, const PHX::MDField< Scalar, Cell, IP > &weighted_measure, const int in_num_cells)
void setupArrays(const Teuchos::RCP< const panzer::BasisIRLayout > &basis, bool computeDerivatives=true)
Sizes/allocates memory for arrays.
int numCells() const
Returns the number of cells in the data layouts.
void evaluateValues_HVol(const PHX::MDField< Scalar, Cell, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP > &jac_det, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv, const PHX::MDField< Scalar, Cell, IP > &weighted_measure, const int in_num_cells)
Description and data layouts associated with a particular basis.
#define TEUCHOS_ASSERT(assertion_test)
void evaluateValues(const PHX::MDField< Scalar, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac, const PHX::MDField< Scalar, Cell, IP > &jac_det, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv, const int in_num_cells=-1)
PHX::MDField< Scalar > buildArray(const std::string &str, int d0) const
PHX::MDField< Scalar > ArrayDynamic