72 using ExecutionSpace =
typename DeviceType::execution_space;
73 using ScratchSpace =
typename ExecutionSpace::scratch_memory_space;
74 using OutputScratchView = Kokkos::View<OutputScalar*,ScratchSpace,Kokkos::MemoryTraits<Kokkos::Unmanaged>>;
75 using PointScratchView = Kokkos::View<PointScalar*, ScratchSpace,Kokkos::MemoryTraits<Kokkos::Unmanaged>>;
77 using TeamPolicy = Kokkos::TeamPolicy<ExecutionSpace>;
78 using TeamMember =
typename TeamPolicy::member_type;
82 OutputFieldType output_;
83 InputPointsType inputPoints_;
85 ordinal_type polyOrder_;
86 ordinal_type numFields_, numPoints_;
88 size_t fad_size_output_;
90 static constexpr ordinal_type numVertices = 4;
91 static constexpr ordinal_type numEdges = 6;
92 static constexpr ordinal_type numEdgesPerFace = 3;
93 static constexpr ordinal_type numFaceFamilies = 2;
94 static constexpr ordinal_type numFaces = 4;
95 static constexpr ordinal_type numVerticesPerFace = 3;
96 static constexpr ordinal_type numInteriorFamilies = 3;
99 const ordinal_type face_vertices[numFaces*numVerticesPerFace] = {0,1,2,
108 const ordinal_type face_edges[numFaces * numEdgesPerFace] = {0,1,2,
114 const ordinal_type edge_start_[numEdges] = {0,1,0,0,1,2};
115 const ordinal_type edge_end_[numEdges] = {1,2,2,3,3,3};
116 const ordinal_type face_family_start_ [numFaceFamilies] = {0,1};
117 const ordinal_type face_family_middle_[numFaceFamilies] = {1,2};
118 const ordinal_type face_family_end_ [numFaceFamilies] = {2,0};
120 const ordinal_type numEdgeFunctions_;
121 const ordinal_type numFaceFunctionsPerFace_;
122 const ordinal_type numFaceFunctions_;
123 const ordinal_type numInteriorFunctionsPerFamily_;
124 const ordinal_type numInteriorFunctions_;
127 const ordinal_type faceOrdinalForInterior_[numInteriorFamilies] = {0,2,3};
128 const ordinal_type faceFamilyForInterior_[numInteriorFamilies] = {0,0,1};
129 const ordinal_type interiorCoordinateOrdinal_[numInteriorFamilies] = {3,0,1};
131 KOKKOS_INLINE_FUNCTION
132 ordinal_type dofOrdinalForFace(
const ordinal_type &faceOrdinal,
133 const ordinal_type &zeroBasedFaceFamily,
134 const ordinal_type &i,
135 const ordinal_type &j)
const
138 const ordinal_type faceDofOffset = numEdgeFunctions_ + faceOrdinal * numFaceFunctionsPerFace_;
144 const ordinal_type max_ij_sum = polyOrder_ - 1;
146 ordinal_type fieldOrdinal = faceDofOffset + zeroBasedFaceFamily;
148 for (ordinal_type ij_sum=1; ij_sum <= max_ij_sum; ij_sum++)
150 for (ordinal_type ii=0; ii<ij_sum; ii++)
153 const ordinal_type jj = ij_sum - ii;
154 if ( (ii == i) && (jj == j))
159 fieldOrdinal += numFaceFamilies;
166 : opType_(opType), output_(output), inputPoints_(inputPoints),
167 polyOrder_(polyOrder),
169 numEdgeFunctions_(polyOrder * numEdges),
170 numFaceFunctionsPerFace_(polyOrder * (polyOrder-1)),
171 numFaceFunctions_(numFaceFunctionsPerFace_*numFaces),
172 numInteriorFunctionsPerFamily_((polyOrder > 2) ? (polyOrder-2)*(polyOrder-1)*polyOrder/6 : 0),
173 numInteriorFunctions_(numInteriorFunctionsPerFamily_ * numInteriorFamilies)
175 numFields_ = output.extent_int(0);
176 numPoints_ = output.extent_int(1);
178 const ordinal_type expectedCardinality = numEdgeFunctions_ + numFaceFunctions_ + numInteriorFunctions_;
184 INTREPID2_TEST_FOR_EXCEPTION(numPoints_ != inputPoints.extent_int(0), std::invalid_argument,
"point counts need to match!");
185 INTREPID2_TEST_FOR_EXCEPTION(numFields_ != expectedCardinality, std::invalid_argument,
"output field size does not match basis cardinality");
188 KOKKOS_INLINE_FUNCTION
189 void computeEdgeLegendre(OutputScratchView &P,
190 const ordinal_type &edgeOrdinal,
191 const PointScalar* lambda)
const
193 const auto & s0 = lambda[edge_start_[edgeOrdinal]];
194 const auto & s1 = lambda[ edge_end_[edgeOrdinal]];
196 Polynomials::shiftedScaledLegendreValues(P, polyOrder_-1, PointScalar(s1), PointScalar(s0+s1));
199 KOKKOS_INLINE_FUNCTION
200 void edgeFunctionValue(OutputScalar &edgeValue_x,
201 OutputScalar &edgeValue_y,
202 OutputScalar &edgeValue_z,
203 const ordinal_type &edgeOrdinal,
204 OutputScratchView &P,
205 const ordinal_type &i,
206 const PointScalar* lambda,
207 const PointScalar* lambda_dx,
208 const PointScalar* lambda_dy,
209 const PointScalar* lambda_dz
212 const auto & s0 = lambda [edge_start_[edgeOrdinal]];
213 const auto & s0_dx = lambda_dx[edge_start_[edgeOrdinal]];
214 const auto & s0_dy = lambda_dy[edge_start_[edgeOrdinal]];
215 const auto & s0_dz = lambda_dz[edge_start_[edgeOrdinal]];
217 const auto & s1 = lambda [ edge_end_[edgeOrdinal]];
218 const auto & s1_dx = lambda_dx[ edge_end_[edgeOrdinal]];
219 const auto & s1_dy = lambda_dy[ edge_end_[edgeOrdinal]];
220 const auto & s1_dz = lambda_dz[ edge_end_[edgeOrdinal]];
222 const auto & P_i = P(i);
223 const PointScalar xWeight = s0 * s1_dx - s1 * s0_dx;
224 const PointScalar yWeight = s0 * s1_dy - s1 * s0_dy;
225 const PointScalar zWeight = s0 * s1_dz - s1 * s0_dz;
226 edgeValue_x = P_i * xWeight;
227 edgeValue_y = P_i * yWeight;
228 edgeValue_z = P_i * zWeight;
231 KOKKOS_INLINE_FUNCTION
232 void computeFaceIntegratedJacobi(OutputScratchView &L_2ip1,
233 const ordinal_type &zeroBasedFaceOrdinal,
234 const ordinal_type &zeroBasedFamilyOrdinal,
235 const ordinal_type &i,
236 const PointScalar* lambda)
const
238 const auto &s0_vertex_number = face_family_start_ [zeroBasedFamilyOrdinal];
239 const auto &s1_vertex_number = face_family_middle_[zeroBasedFamilyOrdinal];
240 const auto &s2_vertex_number = face_family_end_ [zeroBasedFamilyOrdinal];
243 const auto &s0_index = face_vertices[zeroBasedFaceOrdinal * numVerticesPerFace + s0_vertex_number];
244 const auto &s1_index = face_vertices[zeroBasedFaceOrdinal * numVerticesPerFace + s1_vertex_number];
245 const auto &s2_index = face_vertices[zeroBasedFaceOrdinal * numVerticesPerFace + s2_vertex_number];
247 const auto & s0 = lambda[s0_index];
248 const auto & s1 = lambda[s1_index];
249 const auto & s2 = lambda[s2_index];
250 const PointScalar jacobiScaling = s0 + s1 + s2;
252 const double alpha = i*2.0 + 1;
253 Polynomials::shiftedScaledIntegratedJacobiValues(L_2ip1, alpha, polyOrder_-1, s2, jacobiScaling);
256 KOKKOS_INLINE_FUNCTION
257 void faceFunctionValue(OutputScalar &value_x,
258 OutputScalar &value_y,
259 OutputScalar &value_z,
260 const ordinal_type &j,
261 const OutputScratchView &L_2ip1,
262 const OutputScalar &edgeValue_x,
263 const OutputScalar &edgeValue_y,
264 const OutputScalar &edgeValue_z,
265 const PointScalar* lambda)
const
267 const auto & L_2ip1_j = L_2ip1(j);
268 value_x = edgeValue_x * L_2ip1_j;
269 value_y = edgeValue_y * L_2ip1_j;
270 value_z = edgeValue_z * L_2ip1_j;
273 KOKKOS_INLINE_FUNCTION
274 void operator()(
const TeamMember & teamMember )
const
276 const ordinal_type numFunctionsPerFace = polyOrder_ * (polyOrder_ - 1);
277 auto pointOrdinal = teamMember.league_rank();
278 OutputScratchView edge_field_values_at_point, jacobi_values_at_point, other_values_at_point, other_values2_at_point;
279 if (fad_size_output_ > 0) {
280 edge_field_values_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_, fad_size_output_);
281 jacobi_values_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_, fad_size_output_);
282 other_values_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_, fad_size_output_);
283 other_values2_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_, fad_size_output_);
286 edge_field_values_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_);
287 jacobi_values_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_);
288 other_values_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_);
289 other_values2_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_);
292 const auto & x = inputPoints_(pointOrdinal,0);
293 const auto & y = inputPoints_(pointOrdinal,1);
294 const auto & z = inputPoints_(pointOrdinal,2);
297 const PointScalar lambda[4] = {1. - x - y - z, x, y, z};
298 const PointScalar lambda_dx[4] = {-1., 1., 0., 0.};
299 const PointScalar lambda_dy[4] = {-1., 0., 1., 0.};
300 const PointScalar lambda_dz[4] = {-1., 0., 0., 1.};
302 const int num1DEdgeFunctions = polyOrder_;
311 auto & P = edge_field_values_at_point;
313 int fieldOrdinalOffset = 0;
314 for (
int edgeOrdinal=0; edgeOrdinal<numEdges; edgeOrdinal++)
316 computeEdgeLegendre(P, edgeOrdinal, lambda);
318 for (
int i=0; i<num1DEdgeFunctions; i++)
320 auto &output_x = output_(i+fieldOrdinalOffset,pointOrdinal,0);
321 auto &output_y = output_(i+fieldOrdinalOffset,pointOrdinal,1);
322 auto &output_z = output_(i+fieldOrdinalOffset,pointOrdinal,2);
324 edgeFunctionValue(output_x, output_y, output_z,
326 lambda, lambda_dx, lambda_dy, lambda_dz);
328 fieldOrdinalOffset += num1DEdgeFunctions;
334 auto & L_2ip1 = jacobi_values_at_point;
337 const int max_ij_sum = polyOrder_ - 1;
340 int faceFieldOrdinalOffset = fieldOrdinalOffset;
341 for (
int faceOrdinal = 0; faceOrdinal < numFaces; faceOrdinal++) {
342 for (
int familyOrdinal=1; familyOrdinal<=numFaceFamilies; familyOrdinal++)
344 int fieldOrdinal = faceFieldOrdinalOffset + familyOrdinal - 1;
346 for (
int ij_sum=1; ij_sum <= max_ij_sum; ij_sum++)
348 for (
int i=0; i<ij_sum; i++)
350 computeFaceIntegratedJacobi(L_2ip1, faceOrdinal, familyOrdinal-1, i, lambda);
352 const int j = ij_sum - i;
354 const int faceEdgeOrdinal = familyOrdinal-1;
355 const int volumeEdgeOrdinal = face_edges[faceOrdinal * numEdgesPerFace + faceEdgeOrdinal];
356 const int edgeBasisOrdinal = i + volumeEdgeOrdinal*num1DEdgeFunctions;
357 const auto & edgeValue_x = output_(edgeBasisOrdinal,pointOrdinal,0);
358 const auto & edgeValue_y = output_(edgeBasisOrdinal,pointOrdinal,1);
359 const auto & edgeValue_z = output_(edgeBasisOrdinal,pointOrdinal,2);
361 auto & output_x = output_(fieldOrdinal,pointOrdinal,0);
362 auto & output_y = output_(fieldOrdinal,pointOrdinal,1);
363 auto & output_z = output_(fieldOrdinal,pointOrdinal,2);
365 faceFunctionValue(output_x, output_y, output_z, j, L_2ip1, edgeValue_x, edgeValue_y, edgeValue_z, lambda);
367 fieldOrdinal += numFaceFamilies;
370 fieldOrdinalOffset = fieldOrdinal - numFaceFamilies + 1;
372 faceFieldOrdinalOffset += numFunctionsPerFace;
378 const int interiorFieldOrdinalOffset = fieldOrdinalOffset;
379 const int min_ijk_sum = 2;
380 const int max_ijk_sum = polyOrder_-1;
383 const auto & L_2ipj = jacobi_values_at_point;
384 for (
int interiorFamilyOrdinal=1; interiorFamilyOrdinal<=numInteriorFamilies; interiorFamilyOrdinal++)
389 const auto & lambda_m = lambda[interiorCoordinateOrdinal_[interiorFamilyOrdinal-1]];
390 const PointScalar jacobiScaling = 1.0;
392 ordinal_type fieldOrdinal = interiorFieldOrdinalOffset + interiorFamilyOrdinal - 1;
394 const ordinal_type relatedFaceOrdinal = faceOrdinalForInterior_[interiorFamilyOrdinal-1];
395 const ordinal_type relatedFaceFamily = faceFamilyForInterior_ [interiorFamilyOrdinal-1];
397 for (
int ijk_sum=min_ijk_sum; ijk_sum <= max_ijk_sum; ijk_sum++)
399 for (
int i=0; i<ijk_sum-1; i++)
401 for (
int j=1; j<ijk_sum-i; j++)
404 const ordinal_type faceDofOrdinal = dofOrdinalForFace(relatedFaceOrdinal, relatedFaceFamily, i, j);
406 const double alpha = 2 * (i + j);
408 Polynomials::shiftedScaledIntegratedJacobiValues(L_2ipj, alpha, polyOrder_-1, lambda_m, jacobiScaling);
410 const int k = ijk_sum - i - j;
411 const auto & L_k = L_2ipj(k);
412 for (
int d=0; d<3; d++)
414 const auto & E_face_d = output_(faceDofOrdinal,pointOrdinal,d);
415 output_(fieldOrdinal,pointOrdinal,d) = L_k * E_face_d;
417 fieldOrdinal += numInteriorFamilies;
421 fieldOrdinalOffset = fieldOrdinal - numInteriorFamilies + 1;
430 int fieldOrdinalOffset = 0;
437 auto & P_i = edge_field_values_at_point;
438 auto & L_2ip1_j = jacobi_values_at_point;
439 for (
int edgeOrdinal=0; edgeOrdinal<numEdges; edgeOrdinal++)
441 const auto & s0 = lambda[edge_start_[edgeOrdinal]];
442 const auto & s1 = lambda[ edge_end_[edgeOrdinal]];
444 const auto & s0_dx = lambda_dx[edge_start_[edgeOrdinal]];
445 const auto & s0_dy = lambda_dy[edge_start_[edgeOrdinal]];
446 const auto & s0_dz = lambda_dz[edge_start_[edgeOrdinal]];
447 const auto & s1_dx = lambda_dx[ edge_end_[edgeOrdinal]];
448 const auto & s1_dy = lambda_dy[ edge_end_[edgeOrdinal]];
449 const auto & s1_dz = lambda_dz[ edge_end_[edgeOrdinal]];
451 const OutputScalar grad_s0_cross_grad_s1[3] = {s0_dy * s1_dz - s0_dz * s1_dy,
452 s0_dz * s1_dx - s0_dx * s1_dz,
453 s0_dx * s1_dy - s0_dy * s1_dx};
455 Polynomials::shiftedScaledLegendreValues(P_i, polyOrder_-1, PointScalar(s1), PointScalar(s0+s1));
456 for (
int i=0; i<num1DEdgeFunctions; i++)
458 for (
int d=0; d<3; d++)
460 output_(i+fieldOrdinalOffset,pointOrdinal,d) = (i+2) * P_i(i) * grad_s0_cross_grad_s1[d];
463 fieldOrdinalOffset += num1DEdgeFunctions;
484 auto & P_2ip1_j = other_values_at_point;
485 auto & L_2ip1_j_dt = other_values2_at_point;
488 int faceFieldOrdinalOffset = fieldOrdinalOffset;
489 for (
int faceOrdinal=0; faceOrdinal<numFaces; faceOrdinal++)
491 for (
int familyOrdinal=1; familyOrdinal<=numFaceFamilies; familyOrdinal++)
493 int fieldOrdinal = faceFieldOrdinalOffset + familyOrdinal - 1;
495 const auto &s0_vertex_number = face_family_start_ [familyOrdinal-1];
496 const auto &s1_vertex_number = face_family_middle_[familyOrdinal-1];
497 const auto &s2_vertex_number = face_family_end_ [familyOrdinal-1];
500 const auto &s0_index = face_vertices[faceOrdinal * numVerticesPerFace + s0_vertex_number];
501 const auto &s1_index = face_vertices[faceOrdinal * numVerticesPerFace + s1_vertex_number];
502 const auto &s2_index = face_vertices[faceOrdinal * numVerticesPerFace + s2_vertex_number];
504 const auto & s0 = lambda[s0_index];
505 const auto & s1 = lambda[s1_index];
506 const auto & s2 = lambda[s2_index];
507 const PointScalar jacobiScaling = s0 + s1 + s2;
509 const auto & s0_dx = lambda_dx[s0_index];
510 const auto & s0_dy = lambda_dy[s0_index];
511 const auto & s0_dz = lambda_dz[s0_index];
512 const auto & s1_dx = lambda_dx[s1_index];
513 const auto & s1_dy = lambda_dy[s1_index];
514 const auto & s1_dz = lambda_dz[s1_index];
515 const auto & s2_dx = lambda_dx[s2_index];
516 const auto & s2_dy = lambda_dy[s2_index];
517 const auto & s2_dz = lambda_dz[s2_index];
519 const PointScalar grad_s2[3] = {s2_dx, s2_dy, s2_dz};
520 const PointScalar gradJacobiScaling[3] = {s0_dx + s1_dx + s2_dx,
521 s0_dy + s1_dy + s2_dy,
522 s0_dz + s1_dz + s2_dz};
524 const PointScalar grad_s0_cross_grad_s1[3] = {s0_dy * s1_dz - s0_dz * s1_dy,
525 s0_dz * s1_dx - s0_dx * s1_dz,
526 s0_dx * s1_dy - s0_dy * s1_dx};
528 const PointScalar s0_grad_s1_minus_s1_grad_s0[3] = {s0 * s1_dx - s1 * s0_dx,
529 s0 * s1_dy - s1 * s0_dy,
530 s0 * s1_dz - s1 * s0_dz};
532 Polynomials::shiftedScaledLegendreValues (P_i, polyOrder_-1, PointScalar(s1), PointScalar(s0+s1));
539 const int max_ij_sum = polyOrder_ - 1;
540 for (
int ij_sum=1; ij_sum <= max_ij_sum; ij_sum++)
542 for (
int i=0; i<ij_sum; i++)
544 const int j = ij_sum - i;
546 const double alpha = i*2.0 + 1;
548 Polynomials::shiftedScaledJacobiValues (P_2ip1_j, alpha, polyOrder_-1, PointScalar(s2), jacobiScaling);
549 Polynomials::shiftedScaledIntegratedJacobiValues (L_2ip1_j, alpha, polyOrder_-1, PointScalar(s2), jacobiScaling);
550 Polynomials::shiftedScaledIntegratedJacobiValues_dt(L_2ip1_j_dt, alpha, polyOrder_-1, PointScalar(s2), jacobiScaling);
552 const PointScalar & edgeValue = P_i(i);
554 PointScalar grad_L_2ip1_j[3];
555 for (
int d=0; d<3; d++)
557 grad_L_2ip1_j[d] = P_2ip1_j(j-1) * grad_s2[d]
558 + L_2ip1_j_dt(j) * gradJacobiScaling[d];
561 const PointScalar grad_L_2ip1_j_cross_E_i[3] = { grad_L_2ip1_j[1] * edgeValue * s0_grad_s1_minus_s1_grad_s0[2] - grad_L_2ip1_j[2] * edgeValue * s0_grad_s1_minus_s1_grad_s0[1],
562 grad_L_2ip1_j[2] * edgeValue * s0_grad_s1_minus_s1_grad_s0[0] - grad_L_2ip1_j[0] * edgeValue * s0_grad_s1_minus_s1_grad_s0[2],
563 grad_L_2ip1_j[0] * edgeValue * s0_grad_s1_minus_s1_grad_s0[1] - grad_L_2ip1_j[1] * edgeValue * s0_grad_s1_minus_s1_grad_s0[0] };
565 for (
int d=0; d<3; d++)
567 const OutputScalar edgeCurl_d = (i+2.) * P_i(i) * grad_s0_cross_grad_s1[d];
568 output_(fieldOrdinal,pointOrdinal,d) = L_2ip1_j(j) * edgeCurl_d
569 + grad_L_2ip1_j_cross_E_i[d];
572 fieldOrdinal += numFaceFamilies;
575 fieldOrdinalOffset = fieldOrdinal - numFaceFamilies + 1;
577 faceFieldOrdinalOffset += numFunctionsPerFace;
583 auto & L_2ipj = jacobi_values_at_point;
584 auto & P_2ipj = other_values_at_point;
585 auto & L_2ip1 = edge_field_values_at_point;
586 auto & P = other_values2_at_point;
588 const int interiorFieldOrdinalOffset = fieldOrdinalOffset;
589 const int min_ijk_sum = 2;
590 const int max_ijk_sum = polyOrder_-1;
591 for (
int interiorFamilyOrdinal=1; interiorFamilyOrdinal<=numInteriorFamilies; interiorFamilyOrdinal++)
595 const ordinal_type relatedFaceOrdinal = faceOrdinalForInterior_[interiorFamilyOrdinal-1];
596 const ordinal_type relatedFaceFamily = faceFamilyForInterior_ [interiorFamilyOrdinal-1];
599 const auto & m = interiorCoordinateOrdinal_[interiorFamilyOrdinal-1];
600 const auto & lambda_m = lambda[m];
601 const PointScalar jacobiScaling = 1.0;
603 ordinal_type fieldOrdinal = interiorFieldOrdinalOffset + interiorFamilyOrdinal - 1;
605 for (
int ijk_sum=min_ijk_sum; ijk_sum <= max_ijk_sum; ijk_sum++)
607 for (
int i=0; i<ijk_sum-1; i++)
609 computeFaceIntegratedJacobi(L_2ip1, relatedFaceOrdinal, relatedFaceFamily, i, lambda);
611 const ordinal_type faceEdgeOrdinal = relatedFaceFamily;
612 const int volumeEdgeOrdinal = face_edges[relatedFaceOrdinal * numEdgesPerFace + faceEdgeOrdinal];
613 computeEdgeLegendre(P, volumeEdgeOrdinal, lambda);
615 OutputScalar edgeValue[3];
616 edgeFunctionValue(edgeValue[0], edgeValue[1], edgeValue[2], volumeEdgeOrdinal, P, i, lambda, lambda_dx, lambda_dy, lambda_dz);
618 for (
int j=1; j<ijk_sum-i; j++)
621 const ordinal_type faceDofOrdinal = dofOrdinalForFace(relatedFaceOrdinal, relatedFaceFamily, i, j);
623 const double alpha = 2 * (i + j);
625 Polynomials::shiftedScaledIntegratedJacobiValues(L_2ipj, alpha, polyOrder_-1, lambda_m, jacobiScaling);
626 Polynomials::shiftedScaledJacobiValues (P_2ipj, alpha, polyOrder_-1, lambda_m, jacobiScaling);
631 const int k = ijk_sum - i - j;
632 const auto & L_k = L_2ipj(k);
633 const auto & P_km1 = P_2ipj(k-1);
635 const PointScalar grad_L_k[3] = {P_km1 * lambda_dx[m],
636 P_km1 * lambda_dy[m],
637 P_km1 * lambda_dz[m]};
640 OutputScalar E_face[3];
641 faceFunctionValue(E_face[0], E_face[1], E_face[2], j, L_2ip1, edgeValue[0], edgeValue[1], edgeValue[2], lambda);
643 PointScalar grad_L_k_cross_E_face[3] = {grad_L_k[1] * E_face[2] - grad_L_k[2] * E_face[1],
644 grad_L_k[2] * E_face[0] - grad_L_k[0] * E_face[2],
645 grad_L_k[0] * E_face[1] - grad_L_k[1] * E_face[0]};
646 for (
int d=0; d<3; d++)
648 const auto & curl_E_face_d = output_(faceDofOrdinal,pointOrdinal,d);
649 output_(fieldOrdinal,pointOrdinal,d) = L_k * curl_E_face_d + grad_L_k_cross_E_face[d];
652 fieldOrdinal += numInteriorFamilies;
656 fieldOrdinalOffset = fieldOrdinal - numInteriorFamilies + 1;
672 INTREPID2_TEST_FOR_ABORT(
true,
673 ">>> ERROR: (Intrepid2::Hierarchical_HCURL_TET_Functor) Unsupported differential operator");
683 size_t team_shmem_size (
int team_size)
const
686 size_t shmem_size = 0;
687 if (fad_size_output_ > 0)
688 shmem_size += 4 * OutputScratchView::shmem_size(polyOrder_ + 1, fad_size_output_);
690 shmem_size += 4 * OutputScratchView::shmem_size(polyOrder_ + 1);