49#ifndef __INTREPID2_PROJECTIONTOOLSDEFL2_HPP__
50#define __INTREPID2_PROJECTIONTOOLSDEFL2_HPP__
58namespace Experimental {
61template<
typename ViewType1,
typename ViewType2,
typename ViewType3,
62typename ViewType4,
typename ViewType5>
64 ViewType1 basisCoeffs_;
65 const ViewType2 tagToOrdinal_;
66 const ViewType3 targetEPointsRange_;
67 const ViewType4 targetAtTargetEPoints_;
68 const ViewType5 basisAtTargetEPoints_;
69 ordinal_type numVertices_;
73 ViewType4 targetAtTargetEPoints, ViewType5 basisAtTargetEPoints, ordinal_type numVertices) :
74 basisCoeffs_(basisCoeffs), tagToOrdinal_(tagToOrdinal), targetEPointsRange_(targetEPointsRange),
75 targetAtTargetEPoints_(targetAtTargetEPoints), basisAtTargetEPoints_(basisAtTargetEPoints), numVertices_(numVertices) {}
78 KOKKOS_INLINE_FUNCTION
79 operator()(
const ordinal_type ic)
const {
80 for(ordinal_type iv=0; iv<numVertices_; ++iv) {
81 ordinal_type idof = tagToOrdinal_(0, iv, 0);
82 ordinal_type pt = targetEPointsRange_(0,iv).first;
84 basisCoeffs_(ic,idof) = targetAtTargetEPoints_(ic,pt)/basisAtTargetEPoints_(ic,idof,pt,0);
90template<
typename ViewType1,
typename ViewType2,
typename ViewType3,
91typename ViewType4,
typename ViewType5,
typename ViewType6>
93 const ViewType1 basisCoeffs_;
94 const ViewType2 negPartialProj_;
95 const ViewType2 basisDofDofAtBasisEPoints_;
96 const ViewType2 basisAtBasisEPoints_;
97 const ViewType3 basisEWeights_;
98 const ViewType2 wBasisDofAtBasisEPoints_;
99 const ViewType3 targetEWeights_;
100 const ViewType2 basisAtTargetEPoints_;
101 const ViewType2 wBasisDofAtTargetEPoints_;
102 const ViewType4 computedDofs_;
103 const ViewType5 tagToOrdinal_;
104 const ViewType6 targetAtTargetEPoints_;
105 const ViewType2 targetTanAtTargetEPoints_;
106 const ViewType2 refEdgesVec_;
107 ordinal_type fieldDim_;
108 ordinal_type edgeCardinality_;
109 ordinal_type offsetBasis_;
110 ordinal_type offsetTarget_;
111 ordinal_type numVertexDofs_;
112 ordinal_type edgeDim_;
116 const ViewType2 basisAtBasisEPoints,
const ViewType3 basisEWeights,
const ViewType2 wBasisDofAtBasisEPoints,
const ViewType3 targetEWeights,
117 const ViewType2 basisAtTargetEPoints,
const ViewType2 wBasisDofAtTargetEPoints,
const ViewType4 computedDofs,
const ViewType5 tagToOrdinal,
118 const ViewType6 targetAtTargetEPoints,
const ViewType2 targetTanAtTargetEPoints,
const ViewType2 refEdgesVec,
119 ordinal_type fieldDim, ordinal_type edgeCardinality, ordinal_type offsetBasis,
120 ordinal_type offsetTarget, ordinal_type numVertexDofs, ordinal_type edgeDim, ordinal_type iedge) :
121 basisCoeffs_(basisCoeffs), negPartialProj_(negPartialProj), basisDofDofAtBasisEPoints_(basisDofDofAtBasisEPoints),
122 basisAtBasisEPoints_(basisAtBasisEPoints), basisEWeights_(basisEWeights), wBasisDofAtBasisEPoints_(wBasisDofAtBasisEPoints), targetEWeights_(targetEWeights),
123 basisAtTargetEPoints_(basisAtTargetEPoints), wBasisDofAtTargetEPoints_(wBasisDofAtTargetEPoints),
124 computedDofs_(computedDofs), tagToOrdinal_(tagToOrdinal), targetAtTargetEPoints_(targetAtTargetEPoints),
125 targetTanAtTargetEPoints_(targetTanAtTargetEPoints), refEdgesVec_(refEdgesVec),
126 fieldDim_(fieldDim), edgeCardinality_(edgeCardinality), offsetBasis_(offsetBasis),
127 offsetTarget_(offsetTarget), numVertexDofs_(numVertexDofs), edgeDim_(edgeDim), iedge_(iedge)
131 KOKKOS_INLINE_FUNCTION
132 operator()(
const ordinal_type ic)
const {
133 for(ordinal_type j=0; j <edgeCardinality_; ++j) {
134 ordinal_type jdof =tagToOrdinal_(edgeDim_, iedge_, j);
135 for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
136 for(ordinal_type d=0; d <fieldDim_; ++d)
137 basisDofDofAtBasisEPoints_(ic,j,iq) += basisAtBasisEPoints_(ic,jdof,offsetBasis_+iq,d)*refEdgesVec_(iedge_,d);
138 wBasisDofAtBasisEPoints_(ic,j,iq) = basisDofDofAtBasisEPoints_(ic,j,iq)*basisEWeights_(iq);
140 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
141 for(ordinal_type d=0; d <fieldDim_; ++d)
142 wBasisDofAtTargetEPoints_(ic,j,iq) += basisAtTargetEPoints_(ic,jdof,offsetTarget_+iq,d)*refEdgesVec_(iedge_,d)*targetEWeights_(iq);
146 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq)
147 for(ordinal_type d=0; d <fieldDim_; ++d)
148 targetTanAtTargetEPoints_(ic,iq) += targetAtTargetEPoints_(ic,offsetTarget_+iq,d)*refEdgesVec_(iedge_,d);
150 for(ordinal_type j=0; j <numVertexDofs_; ++j) {
151 ordinal_type jdof = computedDofs_(j);
152 for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq)
153 for(ordinal_type d=0; d <fieldDim_; ++d)
154 negPartialProj_(ic,iq) -= basisCoeffs_(ic,jdof)*basisAtBasisEPoints_(ic,jdof,offsetBasis_+iq,d)*refEdgesVec_(iedge_,d);
159template<
typename ViewType1,
typename ViewType2,
typename ViewType3,
160typename ViewType4,
typename ViewType5,
typename ViewType6,
typename ViewType7,
typename ViewType8>
162 const ViewType1 basisCoeffs_;
163 const ViewType2 negPartialProj_;
164 const ViewType2 faceBasisDofAtBasisEPoints_;
165 const ViewType2 basisAtBasisEPoints_;
166 const ViewType3 basisEWeights_;
167 const ViewType2 wBasisDofAtBasisEPoints_;
168 const ViewType3 targetEWeights_;
169 const ViewType2 basisAtTargetEPoints_;
170 const ViewType2 wBasisDofAtTargetEPoints_;
171 const ViewType4 computedDofs_;
172 const ViewType5 tagToOrdinal_;
173 const ViewType6 orts_;
174 const ViewType7 targetAtTargetEPoints_;
175 const ViewType2 targetDofAtTargetEPoints_;
176 const ViewType8 faceParametrization_;
177 ordinal_type fieldDim_;
178 ordinal_type faceCardinality_;
179 ordinal_type offsetBasis_;
180 ordinal_type offsetTarget_;
181 ordinal_type numVertexEdgeDofs_;
182 ordinal_type numFaces_;
183 ordinal_type faceDim_;
187 bool isHCurlBasis_, isHDivBasis_;
190 const ViewType2 basisAtBasisEPoints,
const ViewType3 basisEWeights,
const ViewType2 wBasisDofAtBasisEPoints,
const ViewType3 targetEWeights,
191 const ViewType2 basisAtTargetEPoints,
const ViewType2 wBasisDofAtTargetEPoints,
const ViewType4 computedDofs,
const ViewType5 tagToOrdinal,
192 const ViewType6 orts,
const ViewType7 targetAtTargetEPoints,
const ViewType2 targetDofAtTargetEPoints,
193 const ViewType8 faceParametrization, ordinal_type fieldDim, ordinal_type faceCardinality, ordinal_type offsetBasis,
194 ordinal_type offsetTarget, ordinal_type numVertexEdgeDofs, ordinal_type numFaces, ordinal_type faceDim,
195 ordinal_type dim, ordinal_type iface,
unsigned topoKey,
bool isHCurlBasis,
bool isHDivBasis) :
196 basisCoeffs_(basisCoeffs), negPartialProj_(negPartialProj), faceBasisDofAtBasisEPoints_(faceBasisDofAtBasisEPoints),
197 basisAtBasisEPoints_(basisAtBasisEPoints), basisEWeights_(basisEWeights), wBasisDofAtBasisEPoints_(wBasisDofAtBasisEPoints), targetEWeights_(targetEWeights),
198 basisAtTargetEPoints_(basisAtTargetEPoints), wBasisDofAtTargetEPoints_(wBasisDofAtTargetEPoints),
199 computedDofs_(computedDofs), tagToOrdinal_(tagToOrdinal), orts_(orts), targetAtTargetEPoints_(targetAtTargetEPoints),
200 targetDofAtTargetEPoints_(targetDofAtTargetEPoints), faceParametrization_(faceParametrization),
201 fieldDim_(fieldDim), faceCardinality_(faceCardinality), offsetBasis_(offsetBasis),
202 offsetTarget_(offsetTarget), numVertexEdgeDofs_(numVertexEdgeDofs), numFaces_(numFaces),
203 faceDim_(faceDim), dim_(dim), iface_(iface), topoKey_(topoKey),
204 isHCurlBasis_(isHCurlBasis), isHDivBasis_(isHDivBasis)
208 KOKKOS_INLINE_FUNCTION
209 operator()(
const ordinal_type ic)
const {
211 ordinal_type fOrt[6];
212 orts_(ic).getFaceOrientation(fOrt, numFaces_);
213 ordinal_type ort = fOrt[iface_];
216 typename ViewType3::value_type data[2*3];
217 auto tangents = ViewType3(data, 2, dim_);
219 typename ViewType3::value_type tmp[2] = {};
220 for(ordinal_type j=0; j <faceCardinality_; ++j) {
221 ordinal_type jdof = tagToOrdinal_(faceDim_, iface_, j);
222 for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
224 for(ordinal_type d=0; d <fieldDim_; ++d) {
225 tmp[0] += tangents(0,d)*basisAtBasisEPoints_(ic,jdof,offsetBasis_+iq,d);
226 tmp[1] += tangents(1,d)*basisAtBasisEPoints_(ic,jdof,offsetBasis_+iq,d);
229 for(ordinal_type d=0; d <fieldDim_; ++d) {
230 faceBasisDofAtBasisEPoints_(ic,j,iq,d) = tangents(0,d)*tmp[1]-tangents(1,d)*tmp[0];
231 wBasisDofAtBasisEPoints_(ic,j,iq,d) = faceBasisDofAtBasisEPoints_(ic,j,iq,d) * basisEWeights_(iq);
234 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
236 for(ordinal_type d=0; d <fieldDim_; ++d) {
237 tmp[0] += tangents(0,d)*basisAtTargetEPoints_(ic,jdof,offsetTarget_+iq,d);
238 tmp[1] += tangents(1,d)*basisAtTargetEPoints_(ic,jdof,offsetTarget_+iq,d);
241 for(ordinal_type d=0; d <fieldDim_; ++d)
242 wBasisDofAtTargetEPoints_(ic,j,iq,d) = (tangents(0,d)*tmp[1]-tangents(1,d)*tmp[0]) * targetEWeights_(iq);
246 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
248 for(ordinal_type d=0; d <fieldDim_; ++d) {
249 tmp[0] += tangents(0,d)*targetAtTargetEPoints_(ic,offsetTarget_+iq,d);
250 tmp[1] += tangents(1,d)*targetAtTargetEPoints_(ic,offsetTarget_+iq,d);
253 for(ordinal_type d=0; d <fieldDim_; ++d)
254 targetDofAtTargetEPoints_(ic,iq,d) = (tangents(0,d)*tmp[1]-tangents(1,d)*tmp[0]);
257 for(ordinal_type j=0; j <numVertexEdgeDofs_; ++j) {
258 ordinal_type jdof = computedDofs_(j);
259 for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
261 for(ordinal_type d=0; d <fieldDim_; ++d) {
262 tmp[0] += tangents(0,d)*basisAtBasisEPoints_(ic,jdof,offsetBasis_+iq,d);
263 tmp[1] += tangents(1,d)*basisAtBasisEPoints_(ic,jdof,offsetBasis_+iq,d);
267 for(ordinal_type d=0; d <fieldDim_; ++d)
268 negPartialProj_(ic,iq,d) -= (tangents(0,d)*tmp[1]-tangents(1,d)*tmp[0])*basisCoeffs_(ic,jdof);
272 typename ViewType3::value_type coeff[3];
274 typename ViewType3::value_type data[3*3];
275 auto tangentsAndNormal = ViewType3(data, dim_, dim_);
277 for(ordinal_type d=0; d <dim_; ++d)
278 coeff[d] = tangentsAndNormal(dim_-1,d);
283 for(ordinal_type j=0; j <faceCardinality_; ++j) {
284 ordinal_type jdof = tagToOrdinal_(faceDim_, iface_, j);
285 for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
286 faceBasisDofAtBasisEPoints_(ic,j,iq,0) = 0;
287 for(ordinal_type d=0; d <fieldDim_; ++d)
288 faceBasisDofAtBasisEPoints_(ic,j,iq,0) += coeff[d]*basisAtBasisEPoints_(ic,jdof,offsetBasis_+iq,d);
289 wBasisDofAtBasisEPoints_(ic,j,iq,0) = faceBasisDofAtBasisEPoints_(ic,j,iq,0) * basisEWeights_(iq);
291 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
292 typename ViewType2::value_type sum=0;
293 for(ordinal_type d=0; d <fieldDim_; ++d)
294 sum += coeff[d]*basisAtTargetEPoints_(ic,jdof,offsetTarget_+iq,d);
295 wBasisDofAtTargetEPoints_(ic,j,iq,0) = sum * targetEWeights_(iq);
299 for(ordinal_type d=0; d <fieldDim_; ++d)
300 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq)
301 targetDofAtTargetEPoints_(ic,iq,0) += coeff[d]*targetAtTargetEPoints_(ic,offsetTarget_+iq,d);
303 for(ordinal_type j=0; j <numVertexEdgeDofs_; ++j) {
304 ordinal_type jdof = computedDofs_(j);
305 for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq)
306 for(ordinal_type d=0; d <fieldDim_; ++d)
307 negPartialProj_(ic,iq,0) -= basisCoeffs_(ic,jdof)*coeff[d]*basisAtBasisEPoints_(ic,jdof,offsetBasis_+iq,d);
314template<
typename ViewType1,
typename ViewType2,
typename ViewType3,
315typename ViewType4,
typename ViewType5>
317 const ViewType1 basisCoeffs_;
318 const ViewType2 negPartialProj_;
319 const ViewType2 internalBasisAtBasisEPoints_;
320 const ViewType2 basisAtBasisEPoints_;
321 const ViewType3 basisEWeights_;
322 const ViewType2 wBasisAtBasisEPoints_;
323 const ViewType3 targetEWeights_;
324 const ViewType2 basisAtTargetEPoints_;
325 const ViewType2 wBasisDofAtTargetEPoints_;
326 const ViewType4 computedDofs_;
327 const ViewType5 elemDof_;
328 ordinal_type fieldDim_;
329 ordinal_type numElemDofs_;
330 ordinal_type offsetBasis_;
331 ordinal_type offsetTarget_;
332 ordinal_type numVertexEdgeFaceDofs_;
335 const ViewType2 basisAtBasisEPoints,
const ViewType3 basisEWeights,
const ViewType2 wBasisAtBasisEPoints,
const ViewType3 targetEWeights,
336 const ViewType2 basisAtTargetEPoints,
const ViewType2 wBasisDofAtTargetEPoints,
const ViewType4 computedDofs,
const ViewType5 elemDof,
337 ordinal_type fieldDim, ordinal_type numElemDofs, ordinal_type offsetBasis, ordinal_type offsetTarget, ordinal_type numVertexEdgeFaceDofs) :
338 basisCoeffs_(basisCoeffs), negPartialProj_(negPartialProj), internalBasisAtBasisEPoints_(internalBasisAtBasisEPoints),
339 basisAtBasisEPoints_(basisAtBasisEPoints), basisEWeights_(basisEWeights), wBasisAtBasisEPoints_(wBasisAtBasisEPoints), targetEWeights_(targetEWeights),
340 basisAtTargetEPoints_(basisAtTargetEPoints), wBasisDofAtTargetEPoints_(wBasisDofAtTargetEPoints),
341 computedDofs_(computedDofs), elemDof_(elemDof), fieldDim_(fieldDim), numElemDofs_(numElemDofs), offsetBasis_(offsetBasis),
342 offsetTarget_(offsetTarget), numVertexEdgeFaceDofs_(numVertexEdgeFaceDofs) {}
345 KOKKOS_INLINE_FUNCTION
346 operator()(
const ordinal_type ic)
const {
348 for(ordinal_type j=0; j <numElemDofs_; ++j) {
349 ordinal_type idof = elemDof_(j);
350 for(ordinal_type d=0; d <fieldDim_; ++d) {
351 for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
352 internalBasisAtBasisEPoints_(ic,j,iq,d) = basisAtBasisEPoints_(ic,idof,offsetBasis_+iq,d);
353 wBasisAtBasisEPoints_(ic,j,iq,d) = internalBasisAtBasisEPoints_(ic,j,iq,d) * basisEWeights_(iq);
355 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
356 wBasisDofAtTargetEPoints_(ic,j,iq,d) = basisAtTargetEPoints_(ic,idof,offsetTarget_+iq,d)* targetEWeights_(iq);
360 for(ordinal_type j=0; j < numVertexEdgeFaceDofs_; ++j) {
361 ordinal_type jdof = computedDofs_(j);
362 for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq)
363 for(ordinal_type d=0; d <fieldDim_; ++d) {
364 negPartialProj_(ic,iq,d) -= basisCoeffs_(ic,jdof)*basisAtBasisEPoints_(ic,jdof,offsetBasis_+iq,d);
371template<
typename DeviceType>
372template<
typename BasisType,
373typename ortValueType,
class ...ortProperties>
376 const Kokkos::DynRankView<ortValueType, ortProperties...> orts,
377 const BasisType* cellBasis,
379 const EvalPointsType ePointType) {
380 typedef typename BasisType::scalarType scalarType;
381 typedef Kokkos::DynRankView<scalarType,ortProperties...> ScalarViewType;
382 const auto cellTopo = cellBasis->getBaseCellTopology();
384 ordinal_type dim = cellTopo.getDimension();
385 ordinal_type numCells = ePoints.extent(0);
386 const ordinal_type edgeDim = 1;
387 const ordinal_type faceDim = 2;
389 ordinal_type numVertices = (cellBasis->getDofCount(0, 0) > 0) ? cellTopo.getVertexCount() : 0;
390 ordinal_type numEdges = (cellBasis->getDofCount(edgeDim, 0) > 0) ? cellTopo.getEdgeCount() : 0;
391 ordinal_type numFaces = (cellBasis->getDofCount(faceDim, 0) > 0) ? cellTopo.getFaceCount() : 0;
392 ordinal_type numVols = (cellBasis->getDofCount(dim, 0) > 0);
396 typename RefSubcellParametrization<DeviceType>::ConstViewType subcellParamEdge, subcellParamFace;
404 ScalarViewType workView(
"workView", numCells, projStruct->
getMaxNumEvalPoints(ePointType), dim-1);
407 for(ordinal_type iv=0; iv<numVertices; ++iv) {
408 auto vertexEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getEvalPoints(0,iv,ePointType));
410 ePointsRange(0, iv), Kokkos::ALL()), vertexEPoints);
414 for(ordinal_type ie=0; ie<numEdges; ++ie) {
415 auto edgePointsRange = ePointsRange(edgeDim, ie);
416 auto edgeEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getEvalPoints(edgeDim,ie,ePointType));
418 const auto topoKey = refTopologyKey(edgeDim,ie);
421 Kokkos::RangePolicy<ExecSpaceType, int> (0, numCells),
422 KOKKOS_LAMBDA (
const size_t ic) {
424 ordinal_type eOrt[12];
425 orts(ic).getEdgeOrientation(eOrt, numEdges);
426 ordinal_type ort = eOrt[ie];
429 edgeEPoints, subcellParamEdge, topoKey, ie, ort);
433 for(ordinal_type iface=0; iface<numFaces; ++iface) {
434 auto facePointsRange = ePointsRange(faceDim, iface);
435 auto faceEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getEvalPoints(faceDim,iface,ePointType));
437 const auto topoKey = refTopologyKey(faceDim,iface);
440 Kokkos::RangePolicy<ExecSpaceType, int> (0, numCells),
441 KOKKOS_LAMBDA (
const size_t ic) {
442 ordinal_type fOrt[6];
443 orts(ic).getFaceOrientation(fOrt, numFaces);
444 ordinal_type ort = fOrt[iface];
447 faceEPoints, subcellParamFace, topoKey, iface, ort);
453 auto pointsRange = ePointsRange(dim, 0);
454 auto cellEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getEvalPoints(dim,0,ePointType));
458 const auto topoKey = refTopologyKey(dim,0);
462 Kokkos::RangePolicy<ExecSpaceType, int> (0, numCells),
463 KOKKOS_LAMBDA (
const size_t ic) {
464 ordinal_type ort = 0;
466 orts(ic).getEdgeOrientation(&ort,1);
468 orts(ic).getFaceOrientation(&ort,1);
475template<
typename DeviceType>
476template<
typename basisCoeffsValueType,
class ...basisCoeffsProperties,
477typename funValsValueType,
class ...funValsProperties,
479typename ortValueType,
class ...ortProperties>
482 const Kokkos::DynRankView<funValsValueType,funValsProperties...> targetAtTargetEPoints,
483 const typename BasisType::ScalarViewType targetEPoints,
484 const Kokkos::DynRankView<ortValueType, ortProperties...> orts,
485 const BasisType* cellBasis,
488 typedef typename BasisType::scalarType scalarType;
489 typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
490 typedef Kokkos::pair<ordinal_type,ordinal_type> range_type;
491 const auto cellTopo = cellBasis->getBaseCellTopology();
492 ordinal_type dim = cellTopo.getDimension();
493 ordinal_type numTotalTargetEPoints(targetAtTargetEPoints.extent(1));
494 ordinal_type basisCardinality = cellBasis->getCardinality();
495 ordinal_type numCells = targetAtTargetEPoints.extent(0);
496 const ordinal_type edgeDim = 1;
497 const ordinal_type faceDim = 2;
498 const ordinal_type fieldDim = (targetAtTargetEPoints.rank()==2) ? 1 : targetAtTargetEPoints.extent(2);
500 ordinal_type numVertices = (cellBasis->getDofCount(0, 0) > 0) ? cellTopo.getVertexCount() : 0;
501 ordinal_type numEdges = (cellBasis->getDofCount(1, 0) > 0) ? cellTopo.getEdgeCount() : 0;
502 ordinal_type numFaces = (cellBasis->getDofCount(2, 0) > 0) ? cellTopo.getFaceCount() : 0;
504 ScalarViewType refEdgesVec(
"refEdgesVec", numEdges, dim);
505 ScalarViewType refFacesTangents(
"refFaceTangents", numFaces, dim, 2);
506 ScalarViewType refFacesNormal(
"refFaceNormal", numFaces, dim);
508 ordinal_type numVertexDofs = numVertices;
510 ordinal_type numEdgeDofs(0);
511 for(ordinal_type ie=0; ie<numEdges; ++ie)
512 numEdgeDofs += cellBasis->getDofCount(edgeDim,ie);
514 ordinal_type numFaceDofs(0);
515 for(ordinal_type iface=0; iface<numFaces; ++iface)
516 numFaceDofs += cellBasis->getDofCount(faceDim,iface);
518 Kokkos::View<ordinal_type*, DeviceType> computedDofs(
"computedDofs", numVertexDofs+numEdgeDofs+numFaceDofs);
525 ScalarViewType basisEPoints(
"basisEPoints",numCells,numTotalBasisEPoints, dim);
526 getL2EvaluationPoints(basisEPoints, orts, cellBasis, projStruct, EvalPointsType::BASIS);
528 auto tagToOrdinal = Kokkos::create_mirror_view_and_copy(MemSpaceType(), cellBasis->getAllDofOrdinal());
530 ScalarViewType basisAtBasisEPoints(
"basisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints, fieldDim);
531 ScalarViewType basisAtTargetEPoints(
"basisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints, fieldDim);
534 ScalarViewType nonOrientedBasisAtBasisEPoints(
"nonOrientedBasisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints);
535 ScalarViewType nonOrientedBasisAtTargetEPoints(
"nonOrientedBasisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints);
536 for(ordinal_type ic=0; ic<numCells; ++ic) {
537 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtTargetEPoints,ic,Kokkos::ALL(),Kokkos::ALL()), Kokkos::subview(targetEPoints, ic, Kokkos::ALL(), Kokkos::ALL()));
538 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtBasisEPoints,ic,Kokkos::ALL(),Kokkos::ALL()), Kokkos::subview(basisEPoints, ic, Kokkos::ALL(), Kokkos::ALL()));
541 Kokkos::ALL(),0), nonOrientedBasisAtBasisEPoints, orts, cellBasis);
543 Kokkos::ALL(), Kokkos::ALL(),0), nonOrientedBasisAtTargetEPoints, orts, cellBasis);
546 ScalarViewType nonOrientedBasisAtBasisEPoints(
"nonOrientedBasisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints,fieldDim);
547 ScalarViewType nonOrientedBasisAtTargetEPoints(
"nonOrientedBasisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints,fieldDim);
548 for(ordinal_type ic=0; ic<numCells; ++ic) {
549 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtTargetEPoints,ic,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), Kokkos::subview(targetEPoints, ic, Kokkos::ALL(), Kokkos::ALL()));
550 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtBasisEPoints,ic,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), Kokkos::subview(basisEPoints, ic, Kokkos::ALL(), Kokkos::ALL()));
558 auto hostComputedDofs = Kokkos::create_mirror_view(computedDofs);
559 ordinal_type computedDofsCount = 0;
560 for(ordinal_type iv=0; iv<numVertices; ++iv)
561 hostComputedDofs(computedDofsCount++) = cellBasis->getDofOrdinal(0, iv, 0);
563 for(ordinal_type ie=0; ie<numEdges; ++ie) {
564 ordinal_type edgeCardinality = cellBasis->getDofCount(edgeDim,ie);
565 for(ordinal_type i=0; i<edgeCardinality; ++i)
566 hostComputedDofs(computedDofsCount++) = cellBasis->getDofOrdinal(edgeDim, ie, i);
569 for(ordinal_type iface=0; iface<numFaces; ++iface) {
570 ordinal_type faceCardinality = cellBasis->getDofCount(faceDim,iface);
571 for(ordinal_type i=0; i<faceCardinality; ++i)
572 hostComputedDofs(computedDofsCount++) = cellBasis->getDofOrdinal(faceDim, iface, i);
574 Kokkos::deep_copy(computedDofs,hostComputedDofs);
577 bool isHGradBasis = (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HGRAD);
578 bool isHCurlBasis = (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HCURL);
579 bool isHDivBasis = (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HDIV);
580 ordinal_type faceDofDim = isHCurlBasis ? 3 : 1;
581 ScalarViewType edgeCoeff(
"edgeCoeff", fieldDim);
584 const Kokkos::RangePolicy<ExecSpaceType> policy(0, numCells);
588 auto targetEPointsRange = Kokkos::create_mirror_view_and_copy(MemSpaceType(), projStruct->
getTargetPointsRange());
590 decltype(targetAtTargetEPoints),
decltype(basisAtTargetEPoints)> functorType;
591 Kokkos::parallel_for(policy, functorType(basisCoeffs, tagToOrdinal, targetEPointsRange,
592 targetAtTargetEPoints, basisAtTargetEPoints, numVertices));
596 for(ordinal_type ie=0; ie<numEdges; ++ie) {
597 auto edgeVec = Kokkos::subview(refEdgesVec, ie, Kokkos::ALL());
602 }
else if(isHDivBasis) {
605 deep_copy(edgeVec, 1.0);
608 ordinal_type edgeCardinality = cellBasis->getDofCount(edgeDim,ie);
609 ordinal_type numBasisEPoints = range_size(basisEPointsRange(edgeDim, ie));
610 ordinal_type numTargetEPoints = range_size(targetEPointsRange(edgeDim, ie));
612 ScalarViewType basisDofAtBasisEPoints(
"BasisDofAtBasisEPoints",numCells,edgeCardinality, numBasisEPoints);
613 ScalarViewType tragetDofAtTargetEPoints(
"TargetDofAtTargetEPoints",numCells, numTargetEPoints);
614 ScalarViewType weightedBasisAtBasisEPoints(
"weightedTanBasisAtBasisEPoints",numCells,edgeCardinality, numBasisEPoints);
615 ScalarViewType weightedBasisAtTargetEPoints(
"weightedTanBasisAtTargetEPoints",numCells,edgeCardinality, numTargetEPoints);
616 ScalarViewType negPartialProj(
"negPartialProj", numCells, numBasisEPoints);
618 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getTargetEvalWeights(edgeDim,ie));
619 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getBasisEvalWeights(edgeDim,ie));
622 ordinal_type offsetBasis = basisEPointsRange(edgeDim, ie).first;
623 ordinal_type offsetTarget = targetEPointsRange(edgeDim, ie).first;
627 decltype(computedDofs),
decltype(tagToOrdinal),
decltype(targetAtTargetEPoints)> functorTypeEdge;
629 Kokkos::parallel_for(policy, functorTypeEdge(basisCoeffs, negPartialProj, basisDofAtBasisEPoints,
630 basisAtBasisEPoints, basisEWeights, weightedBasisAtBasisEPoints, targetEWeights,
631 basisAtTargetEPoints, weightedBasisAtTargetEPoints, computedDofs, tagToOrdinal,
632 targetAtTargetEPoints,tragetDofAtTargetEPoints, refEdgesVec, fieldDim,
633 edgeCardinality, offsetBasis, offsetTarget, numVertexDofs, edgeDim, ie));
636 ScalarViewType edgeMassMat_(
"edgeMassMat_", numCells, edgeCardinality, edgeCardinality),
637 edgeRhsMat_(
"rhsMat_", numCells, edgeCardinality);
644 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
645 ScalarViewType t_(
"t",numCells, edgeCardinality);
646 WorkArrayViewType w_(
"w",numCells, edgeCardinality);
648 auto edgeDof = Kokkos::subview(tagToOrdinal, edgeDim, ie, Kokkos::ALL());
651 edgeSystem.
solve(basisCoeffs, edgeMassMat_, edgeRhsMat_, t_, w_, edgeDof, edgeCardinality);
654 typename RefSubcellParametrization<DeviceType>::ConstViewType subcellParamFace;
658 for(ordinal_type iface=0; iface<numFaces; ++iface) {
659 const auto topoKey = refTopologyKey(faceDim,iface);
660 ordinal_type faceCardinality = cellBasis->getDofCount(faceDim,iface);
662 ordinal_type numTargetEPoints = range_size(targetEPointsRange(faceDim, iface));
663 ordinal_type numBasisEPoints = range_size(basisEPointsRange(faceDim, iface));
665 ScalarViewType faceBasisDofAtBasisEPoints(
"faceBasisDofAtBasisEPoints",numCells,faceCardinality, numBasisEPoints,faceDofDim);
666 ScalarViewType wBasisDofAtBasisEPoints(
"weightedBasisDofAtBasisEPoints",numCells,faceCardinality, numBasisEPoints,faceDofDim);
668 ScalarViewType faceBasisAtTargetEPoints(
"faceBasisDofAtTargetEPoints",numCells,faceCardinality, numTargetEPoints,faceDofDim);
669 ScalarViewType wBasisDofAtTargetEPoints(
"weightedBasisDofAtTargetEPoints",numCells,faceCardinality, numTargetEPoints,faceDofDim);
671 ScalarViewType targetDofAtTargetEPoints(
"targetDofAtTargetEPoints",numCells, numTargetEPoints,faceDofDim);
672 ScalarViewType negPartialProj(
"negComputedProjection", numCells,numBasisEPoints,faceDofDim);
674 ordinal_type offsetBasis = basisEPointsRange(faceDim, iface).first;
675 ordinal_type offsetTarget = targetEPointsRange(faceDim, iface).first;
676 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getTargetEvalWeights(faceDim,iface));
677 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getBasisEvalWeights(faceDim,iface));
681 decltype(computedDofs),
decltype(tagToOrdinal),
decltype(orts),
decltype(targetAtTargetEPoints),
decltype(subcellParamFace)> functorTypeFace;
683 Kokkos::parallel_for(policy, functorTypeFace(basisCoeffs, negPartialProj,faceBasisDofAtBasisEPoints,
684 basisAtBasisEPoints, basisEWeights, wBasisDofAtBasisEPoints, targetEWeights,
685 basisAtTargetEPoints, wBasisDofAtTargetEPoints, computedDofs, tagToOrdinal,
686 orts, targetAtTargetEPoints,targetDofAtTargetEPoints,
687 subcellParamFace, fieldDim, faceCardinality, offsetBasis,
688 offsetTarget, numVertexDofs+numEdgeDofs, numFaces, faceDim,
689 dim, iface, topoKey, isHCurlBasis, isHDivBasis));
691 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
692 ScalarViewType faceMassMat_(
"faceMassMat_", numCells, faceCardinality, faceCardinality),
693 faceRhsMat_(
"rhsMat_", numCells, faceCardinality);
699 ScalarViewType t_(
"t",numCells, faceCardinality);
700 WorkArrayViewType w_(
"w",numCells,faceCardinality);
702 auto faceDof = Kokkos::subview(tagToOrdinal, faceDim, iface, Kokkos::ALL());
705 faceSystem.
solve(basisCoeffs, faceMassMat_, faceRhsMat_, t_, w_, faceDof, faceCardinality);
708 ordinal_type numElemDofs = cellBasis->getDofCount(dim,0);
713 auto cellDofs = Kokkos::subview(tagToOrdinal, dim, 0, Kokkos::ALL());
715 range_type cellPointsRange = targetEPointsRange(dim, 0);
717 ordinal_type numTargetEPoints = range_size(targetEPointsRange(dim,0));
718 ordinal_type numBasisEPoints = range_size(basisEPointsRange(dim,0));
720 ScalarViewType internalBasisAtBasisEPoints(
"internalBasisAtBasisEPoints",numCells,numElemDofs, numBasisEPoints, fieldDim);
721 ScalarViewType negPartialProj(
"negPartialProj", numCells, numBasisEPoints, fieldDim);
723 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getTargetEvalWeights(dim,0));
724 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getBasisEvalWeights(dim,0));
725 ordinal_type offsetBasis = basisEPointsRange(dim, 0).first;
726 ordinal_type offsetTarget = targetEPointsRange(dim, 0).first;
729 ScalarViewType wBasisAtBasisEPoints(
"weightedBasisAtBasisEPoints",numCells,numElemDofs, numBasisEPoints,fieldDim);
730 ScalarViewType wBasisDofAtTargetEPoints(
"weightedBasisAtTargetEPoints",numCells,numElemDofs, numTargetEPoints,fieldDim);
732 typedef ComputeBasisCoeffsOnCells_L2<
decltype(basisCoeffs), ScalarViewType,
decltype(basisEWeights),
decltype(computedDofs),
decltype(cellDofs)> functorType;
733 Kokkos::parallel_for(policy, functorType( basisCoeffs, negPartialProj, internalBasisAtBasisEPoints,
734 basisAtBasisEPoints, basisEWeights, wBasisAtBasisEPoints, targetEWeights, basisAtTargetEPoints, wBasisDofAtTargetEPoints,
735 computedDofs, cellDofs, fieldDim, numElemDofs, offsetBasis, offsetTarget, numVertexDofs+numEdgeDofs+numFaceDofs));
737 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
738 ScalarViewType cellMassMat_(
"cellMassMat_", numCells, numElemDofs, numElemDofs),
739 cellRhsMat_(
"rhsMat_", numCells, numElemDofs);
744 Kokkos::subview(wBasisDofAtTargetEPoints,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),0));
749 ScalarViewType t_(
"t",numCells, numElemDofs);
750 WorkArrayViewType w_(
"w",numCells,numElemDofs);
752 cellSystem.
solve(basisCoeffs, cellMassMat_, cellRhsMat_, t_, w_, cellDofs, numElemDofs);
756template<
typename ViewType1,
typename ViewType2>
758 const ViewType1 basisAtBasisEPoints_;
759 const ViewType2 basisEWeights_;
760 const ViewType1 wBasisAtBasisEPoints_;
761 const ViewType2 targetEWeights_;
762 const ViewType1 basisAtTargetEPoints_;
763 const ViewType1 wBasisDofAtTargetEPoints_;
764 ordinal_type fieldDim_;
765 ordinal_type numElemDofs_;
767 MultiplyBasisByWeights(
const ViewType1 basisAtBasisEPoints,
const ViewType2 basisEWeights,
const ViewType1 wBasisAtBasisEPoints,
const ViewType2 targetEWeights,
768 const ViewType1 basisAtTargetEPoints,
const ViewType1 wBasisDofAtTargetEPoints,
769 ordinal_type fieldDim, ordinal_type numElemDofs) :
770 basisAtBasisEPoints_(basisAtBasisEPoints), basisEWeights_(basisEWeights), wBasisAtBasisEPoints_(wBasisAtBasisEPoints), targetEWeights_(targetEWeights),
771 basisAtTargetEPoints_(basisAtTargetEPoints), wBasisDofAtTargetEPoints_(wBasisDofAtTargetEPoints),
772 fieldDim_(fieldDim), numElemDofs_(numElemDofs) {}
775 KOKKOS_INLINE_FUNCTION
776 operator()(
const ordinal_type ic)
const {
778 for(ordinal_type j=0; j <numElemDofs_; ++j) {
779 for(ordinal_type d=0; d <fieldDim_; ++d) {
780 for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
781 wBasisAtBasisEPoints_(ic,j,iq,d) = basisAtBasisEPoints_(ic,j,iq,d) * basisEWeights_(iq);
783 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
784 wBasisDofAtTargetEPoints_(ic,j,iq,d) = basisAtTargetEPoints_(ic,j,iq,d)* targetEWeights_(iq);
791template<
typename DeviceType>
792template<
typename BasisType>
795 const BasisType* cellBasis,
797 const EvalPointsType ePointType) {
799 ordinal_type dim = cellBasis->getBaseCellTopology().getDimension();
800 auto cellEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getEvalPoints(dim,0,ePointType));
804template<
typename DeviceType>
805template<
typename basisCoeffsValueType,
class ...basisCoeffsProperties,
806typename funValsValueType,
class ...funValsProperties,
808typename ortValueType,
class ...ortProperties>
811 const Kokkos::DynRankView<funValsValueType,funValsProperties...> targetAtTargetEPoints,
812 const Kokkos::DynRankView<ortValueType, ortProperties...> orts,
813 const BasisType* cellBasis,
816 typedef typename BasisType::scalarType scalarType;
817 typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
818 const auto cellTopo = cellBasis->getBaseCellTopology();
820 ordinal_type dim = cellTopo.getDimension();
822 auto basisEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),
824 auto targetEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),
828 ordinal_type numTotalTargetEPoints(targetAtTargetEPoints.extent(1));
829 ordinal_type basisCardinality = cellBasis->getCardinality();
830 ordinal_type numCells = targetAtTargetEPoints.extent(0);
831 const ordinal_type fieldDim = (targetAtTargetEPoints.rank()==2) ? 1 : targetAtTargetEPoints.extent(2);
834 ScalarViewType basisAtBasisEPoints(
"basisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints, fieldDim);
835 ScalarViewType basisAtTargetEPoints(
"basisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints, fieldDim);
838 ScalarViewType nonOrientedBasisAtBasisEPoints(
"nonOrientedBasisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints);
839 ScalarViewType nonOrientedBasisAtTargetEPoints(
"nonOrientedBasisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints);
840 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL()), targetEPoints);
841 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL()), basisEPoints);
846 Kokkos::ALL(),0), nonOrientedBasisAtBasisEPoints, orts, cellBasis);
848 Kokkos::ALL(), Kokkos::ALL(),0), nonOrientedBasisAtTargetEPoints, orts, cellBasis);
851 ScalarViewType nonOrientedBasisAtBasisEPoints(
"nonOrientedBasisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints,fieldDim);
852 ScalarViewType nonOrientedBasisAtTargetEPoints(
"nonOrientedBasisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints,fieldDim);
853 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), targetEPoints);
854 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), basisEPoints);
856 RealSpaceTools<DeviceType>::clone(nonOrientedBasisAtTargetEPoints, Kokkos::subview(nonOrientedBasisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()));
863 const Kokkos::RangePolicy<ExecSpaceType> policy(0, numCells);
864 ordinal_type numElemDofs = cellBasis->getCardinality();
866 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getTargetEvalWeights(dim,0));
867 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getBasisEvalWeights(dim,0));
869 ScalarViewType wBasisAtBasisEPoints(
"weightedBasisAtBasisEPoints",numCells,numElemDofs, numTotalBasisEPoints,fieldDim);
870 ScalarViewType wBasisDofAtTargetEPoints(
"weightedBasisAtTargetEPoints",numCells,numElemDofs, numTotalTargetEPoints,fieldDim);
873 Kokkos::parallel_for(
"Multiply basis by weights", policy, functorType(basisAtBasisEPoints, basisEWeights,
874 wBasisAtBasisEPoints, targetEWeights, basisAtTargetEPoints, wBasisDofAtTargetEPoints, fieldDim, numElemDofs));
876 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
877 ScalarViewType cellMassMat_(
"cellMassMat_", numCells, numElemDofs, numElemDofs),
878 cellRhsMat_(
"rhsMat_", numCells, numElemDofs);
883 Kokkos::subview(wBasisDofAtTargetEPoints,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),0));
887 Kokkos::DynRankView<ordinal_type,Kokkos::HostSpace> hCellDofs(
"cellDoFs", numElemDofs);
888 for(ordinal_type i=0; i<numElemDofs; ++i) hCellDofs(i) = i;
889 auto cellDofs = Kokkos::create_mirror_view_and_copy(MemSpaceType(),hCellDofs);
891 ScalarViewType t_(
"t",numCells, numElemDofs);
892 WorkArrayViewType w_(
"w",numCells,numElemDofs);
894 cellSystem.
solve(basisCoeffs, cellMassMat_, cellRhsMat_, t_, w_, cellDofs, numElemDofs);
897template<
typename DeviceType>
898template<
typename basisViewType,
typename targetViewType,
typename BasisType>
901 const targetViewType targetAtTargetEPoints,
902 const BasisType* cellBasis,
905 typedef typename BasisType::scalarType scalarType;
906 typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
907 const auto cellTopo = cellBasis->getBaseCellTopology();
909 ordinal_type dim = cellTopo.getDimension();
911 auto basisEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),
913 auto targetEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),
916 ordinal_type numTotalTargetEPoints(targetAtTargetEPoints.extent(1));
917 ordinal_type basisCardinality = cellBasis->getCardinality();
918 ordinal_type numCells = targetAtTargetEPoints.extent(0);
919 const ordinal_type fieldDim = (targetAtTargetEPoints.rank()==2) ? 1 : targetAtTargetEPoints.extent(2);
922 ScalarViewType basisAtBasisEPoints(
"basisAtBasisEPoints",1,basisCardinality, numTotalBasisEPoints, fieldDim);
923 ScalarViewType basisAtTargetEPoints(
"basisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints, fieldDim);
926 cellBasis->getValues(Kokkos::subview(basisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),0), targetEPoints);
927 cellBasis->getValues(Kokkos::subview(basisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL(),0), basisEPoints);
932 cellBasis->getValues(Kokkos::subview(basisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), targetEPoints);
933 cellBasis->getValues(Kokkos::subview(basisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), basisEPoints);
939 const Kokkos::RangePolicy<ExecSpaceType> policy(0, numCells);
940 ordinal_type numElemDofs = cellBasis->getCardinality();
942 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getTargetEvalWeights(dim,0));
943 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->
getBasisEvalWeights(dim,0));
945 ScalarViewType wBasisAtBasisEPoints(
"weightedBasisAtBasisEPoints",1,numElemDofs, numTotalBasisEPoints,fieldDim);
946 ScalarViewType wBasisDofAtTargetEPoints(
"weightedBasisAtTargetEPoints",numCells,numElemDofs, numTotalTargetEPoints,fieldDim);
947 Kokkos::DynRankView<ordinal_type, DeviceType> cellDofs(
"cellDoFs", numElemDofs);
949 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,numElemDofs),
950 KOKKOS_LAMBDA (
const int &j) {
951 for(ordinal_type d=0; d <fieldDim; ++d) {
952 for(ordinal_type iq=0; iq < numTotalBasisEPoints; ++iq)
953 wBasisAtBasisEPoints(0,j,iq,d) = basisAtBasisEPoints(0,j,iq,d) * basisEWeights(iq);
954 for(ordinal_type iq=0; iq <numTotalTargetEPoints; ++iq) {
955 wBasisDofAtTargetEPoints(0,j,iq,d) = basisAtTargetEPoints(0,j,iq,d)* targetEWeights(iq);
962 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
963 ScalarViewType cellMassMat_(
"cellMassMat_", 1, numElemDofs, numElemDofs),
964 cellRhsMat_(
"rhsMat_", numCells, numElemDofs);
969 Kokkos::subview(wBasisDofAtTargetEPoints,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),0));
973 ScalarViewType t_(
"t",1, numElemDofs);
974 WorkArrayViewType w_(
"w",numCells,numElemDofs);
976 cellSystem.
solve(basisCoeffs, cellMassMat_, cellRhsMat_, t_, w_, cellDofs, numElemDofs);
Header file for the abstract base class Intrepid2::DefaultCubatureFactory.
An helper class to compute the evaluation points and weights needed for performing projections.
ordinal_type getNumBasisEvalPoints()
Returns number of basis evaluation points.
const range_tag getBasisPointsRange() const
Returns the range tag of the basis evaluation points subcells.
const range_tag getTargetPointsRange() const
Returns the range of the target function evaluation points on subcells.
view_type getEvalPoints(const ordinal_type subCellDim, const ordinal_type subCellId, EvalPointsType type) const
Returns the basis/target evaluation points on a subcell.
const key_tag getTopologyKey() const
Returns the key tag for subcells.
view_type getBasisEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the basis evaluation weights on a subcell.
ordinal_type getMaxNumEvalPoints(const EvalPointsType type) const
Returns the maximum number of evaluation points across all the subcells.
view_type getTargetEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the function evaluation weights on a subcell.
const range_tag getPointsRange(const EvalPointsType type) const
Returns the range tag of the basis/target evaluation points in subcells.
static ConstViewType get(const ordinal_type subcellDim, const unsigned parentCellKey)
Returns a Kokkos view with the coefficients of the parametrization maps for the edges or faces of a r...