Intrepid2
Intrepid2_ProjectionToolsDefHGRAD.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Intrepid2 Package
4 //
5 // Copyright 2007 NTESS and the Intrepid2 contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
16 #ifndef __INTREPID2_PROJECTIONTOOLSDEFHGRAD_HPP__
17 #define __INTREPID2_PROJECTIONTOOLSDEFHGRAD_HPP__
18 
20 #include "Intrepid2_ArrayTools.hpp"
22 
23 #include "Intrepid2_TestUtils.hpp"
24 
25 namespace Intrepid2 {
26 
27 namespace FunctorsProjectionTools {
28 
29 template<typename ViewType1, typename ViewType2, typename ViewType3,
30 typename ViewType4, typename ViewType5>
32  ViewType1 basisCoeffs_;
33  const ViewType2 tagToOrdinal_;
34  const ViewType3 targetEPointsRange_;
35  const ViewType4 targetAtTargetEPoints_;
36  const ViewType5 basisAtTargetEPoints_;
37  ordinal_type numVertices_;
38 
39 
40  ComputeBasisCoeffsOnVertices_HGRAD(ViewType1 basisCoeffs, ViewType2 tagToOrdinal, ViewType3 targetEPointsRange,
41  ViewType4 targetAtTargetEPoints, ViewType5 basisAtTargetEPoints, ordinal_type numVertices) :
42  basisCoeffs_(basisCoeffs), tagToOrdinal_(tagToOrdinal), targetEPointsRange_(targetEPointsRange),
43  targetAtTargetEPoints_(targetAtTargetEPoints), basisAtTargetEPoints_(basisAtTargetEPoints), numVertices_(numVertices) {}
44 
45  void
46  KOKKOS_INLINE_FUNCTION
47  operator()(const ordinal_type ic) const {
48 
49 
50  for(ordinal_type iv=0; iv<numVertices_; ++iv) {
51  ordinal_type idof = tagToOrdinal_(0, iv, 0);
52  ordinal_type pt = targetEPointsRange_(0,iv).first;
53  //the value of the basis at the vertex might be different than 1; HGrad basis, so the function is scalar
54  basisCoeffs_(ic,idof) = targetAtTargetEPoints_(ic,pt)/basisAtTargetEPoints_(idof,pt,0);
55  }
56  }
57 };
58 
59 template<typename ViewType1, typename ViewType2, typename ViewType3,
60 typename ViewType4, typename ViewType5>
62  const ViewType1 basisCoeffs_;
63  const ViewType1 negPartialProjGrad_;
64  const ViewType1 basisTanAtEPoints_;
65  const ViewType1 basisGradAtBasisGradEPoints_;
66  const ViewType2 basisGradEWeights_;
67  const ViewType1 wBasisAtBasisGradEPoints_;
68  const ViewType2 targetGradEWeights_;
69  const ViewType1 basisGradAtTargetGradEPoints_;
70  const ViewType1 wBasisAtTargetGradEPoints_;
71  const ViewType3 computedDofs_;
72  const ViewType4 tagToOrdinal_;
73  const ViewType1 targetGradTanAtTargetGradEPoints_;
74  const ViewType5 targetGradAtTargetGradEPoints_;
75  const ViewType1 refEdgeTan_;
76  ordinal_type edgeCardinality_;
77  ordinal_type offsetBasis_;
78  ordinal_type offsetTarget_;
79  ordinal_type numVertexDofs_;
80  ordinal_type edgeDim_;
81  ordinal_type dim_;
82  ordinal_type iedge_;
83 
84  ComputeBasisCoeffsOnEdges_HGRAD(const ViewType1 basisCoeffs, ViewType1 negPartialProjGrad, const ViewType1 basisTanAtEPoints,
85  const ViewType1 basisGradAtBasisGradEPoints, const ViewType2 basisGradEWeights, const ViewType1 wBasisAtBasisGradEPoints, const ViewType2 targetGradEWeights,
86  const ViewType1 basisGradAtTargetGradEPoints, const ViewType1 wBasisAtTargetGradEPoints, const ViewType3 computedDofs, const ViewType4 tagToOrdinal,
87  const ViewType1 targetGradTanAtTargetGradEPoints, const ViewType5 targetGradAtTargetGradEPoints, const ViewType1 refEdgeTan,
88  ordinal_type edgeCardinality, ordinal_type offsetBasis,
89  ordinal_type offsetTarget, ordinal_type numVertexDofs, ordinal_type edgeDim, ordinal_type dim, ordinal_type iedge) :
90  basisCoeffs_(basisCoeffs), negPartialProjGrad_(negPartialProjGrad), basisTanAtEPoints_(basisTanAtEPoints),
91  basisGradAtBasisGradEPoints_(basisGradAtBasisGradEPoints), basisGradEWeights_(basisGradEWeights), wBasisAtBasisGradEPoints_(wBasisAtBasisGradEPoints), targetGradEWeights_(targetGradEWeights),
92  basisGradAtTargetGradEPoints_(basisGradAtTargetGradEPoints), wBasisAtTargetGradEPoints_(wBasisAtTargetGradEPoints),
93  computedDofs_(computedDofs), tagToOrdinal_(tagToOrdinal), targetGradTanAtTargetGradEPoints_(targetGradTanAtTargetGradEPoints),
94  targetGradAtTargetGradEPoints_(targetGradAtTargetGradEPoints), refEdgeTan_(refEdgeTan),
95  edgeCardinality_(edgeCardinality), offsetBasis_(offsetBasis),
96  offsetTarget_(offsetTarget), numVertexDofs_(numVertexDofs), edgeDim_(edgeDim), dim_(dim), iedge_(iedge)
97  {}
98 
99  void
100  KOKKOS_INLINE_FUNCTION
101  operator()(const ordinal_type ic) const {
102 
103  ordinal_type numBasisGradEPoints = basisGradEWeights_.extent(0);
104  ordinal_type numTargetGradEPoints = targetGradEWeights_.extent(0);
105 
106  //this loop could be moved outside of cell loop
107  for(ordinal_type j=0; j <edgeCardinality_; ++j) {
108  ordinal_type jdof = tagToOrdinal_(edgeDim_, iedge_, j);
109 
110  for(ordinal_type iq=0; iq <numBasisGradEPoints; ++iq) {
111  typename ViewType1::value_type tmp =0;
112  for(ordinal_type d=0; d <dim_; ++d)
113  tmp += basisGradAtBasisGradEPoints_(jdof,offsetBasis_+iq,d)*refEdgeTan_(d);
114  basisTanAtEPoints_(0,j,iq) = tmp;
115  wBasisAtBasisGradEPoints_(ic,j,iq) = tmp*basisGradEWeights_(iq);
116  }
117  for(ordinal_type iq=0; iq <numTargetGradEPoints; ++iq) {
118  for(ordinal_type d=0; d <dim_; ++d)
119  wBasisAtTargetGradEPoints_(ic,j,iq) += basisGradAtTargetGradEPoints_(jdof,offsetTarget_+iq,d)*refEdgeTan_(d)*targetGradEWeights_(iq);
120  }
121  }
122 
123  for(ordinal_type iq=0; iq <numTargetGradEPoints; ++iq)
124  for(ordinal_type d=0; d <dim_; ++d)
125  targetGradTanAtTargetGradEPoints_(ic,iq) += targetGradAtTargetGradEPoints_(ic,offsetTarget_+iq,d)*refEdgeTan_(d);
126 
127  for(ordinal_type j=0; j <numVertexDofs_; ++j) {
128  ordinal_type jdof = computedDofs_(j);
129  for(ordinal_type iq=0; iq <numBasisGradEPoints; ++iq)
130  for(ordinal_type d=0; d <dim_; ++d)
131  negPartialProjGrad_(ic,iq) -= basisCoeffs_(ic,jdof)*basisGradAtBasisGradEPoints_(jdof,offsetBasis_+iq,d)*refEdgeTan_(d);
132  }
133  }
134 };
135 
136 template<typename ViewType1, typename ViewType2, typename ViewType3,
137 typename ViewType4, typename ViewType5>
139  const ViewType1 basisCoeffs_;
140  const ViewType1 negPartialProjGrad_;
141  const ViewType1 faceBasisGradAtGradEPoints_;
142  const ViewType1 basisGradAtBasisGradEPoints_;
143  const ViewType2 basisGradEWeights_;
144  const ViewType1 wBasisGradAtGradEPoints_;
145  const ViewType2 targetGradEWeights_;
146  const ViewType1 basisGradAtTargetGradEPoints_;
147  const ViewType1 wBasisGradBasisAtTargetGradEPoints_;
148  const ViewType3 computedDofs_;
149  const ViewType4 tagToOrdinal_;
150  const ViewType1 targetGradTanAtTargetGradEPoints_;
151  const ViewType5 targetGradAtTargetGradEPoints_;
152  const ViewType1 refFaceNormal_;
153  ordinal_type faceCardinality_;
154  ordinal_type offsetBasisGrad_;
155  ordinal_type offsetTargetGrad_;
156  ordinal_type numVertexEdgeDofs_;
157  ordinal_type numFaces_;
158  ordinal_type faceDim_;
159  ordinal_type dim_;
160  ordinal_type iface_;
161 
162  ComputeBasisCoeffsOnFaces_HGRAD(const ViewType1 basisCoeffs, ViewType1 negPartialProjGrad, const ViewType1 faceBasisGradAtGradEPoints,
163  const ViewType1 basisGradAtBasisGradEPoints, const ViewType2 basisGradEWeights, const ViewType1 wBasisGradAtGradEPoints, const ViewType2 targetGradEWeights,
164  const ViewType1 basisGradAtTargetGradEPoints, const ViewType1 wBasisGradBasisAtTargetGradEPoints, const ViewType3 computedDofs, const ViewType4 tagToOrdinal,
165  const ViewType1 targetGradTanAtTargetGradEPoints, const ViewType5 targetGradAtTargetGradEPoints,
166  const ViewType1 refFaceNormal, ordinal_type faceCardinality, ordinal_type offsetBasisGrad,
167  ordinal_type offsetTargetGrad, ordinal_type numVertexEdgeDofs, ordinal_type numFaces, ordinal_type faceDim,
168  ordinal_type dim, ordinal_type iface) :
169  basisCoeffs_(basisCoeffs), negPartialProjGrad_(negPartialProjGrad), faceBasisGradAtGradEPoints_(faceBasisGradAtGradEPoints),
170  basisGradAtBasisGradEPoints_(basisGradAtBasisGradEPoints), basisGradEWeights_(basisGradEWeights), wBasisGradAtGradEPoints_(wBasisGradAtGradEPoints),
171  targetGradEWeights_(targetGradEWeights),
172  basisGradAtTargetGradEPoints_(basisGradAtTargetGradEPoints), wBasisGradBasisAtTargetGradEPoints_(wBasisGradBasisAtTargetGradEPoints),
173  computedDofs_(computedDofs), tagToOrdinal_(tagToOrdinal), targetGradTanAtTargetGradEPoints_(targetGradTanAtTargetGradEPoints),
174  targetGradAtTargetGradEPoints_(targetGradAtTargetGradEPoints), refFaceNormal_(refFaceNormal),
175  faceCardinality_(faceCardinality), offsetBasisGrad_(offsetBasisGrad),
176  offsetTargetGrad_(offsetTargetGrad), numVertexEdgeDofs_(numVertexEdgeDofs), numFaces_(numFaces),
177  faceDim_(faceDim), dim_(dim), iface_(iface)
178  {}
179 
180  void
181  KOKKOS_INLINE_FUNCTION
182  operator()(const ordinal_type ic) const {
183 
184  ordinal_type numBasisGradEPoints = basisGradEWeights_.extent(0);
185  ordinal_type numTargetGradEPoints = targetGradEWeights_.extent(0);
186 
187  //normal
188  typename ViewType2::value_type n[3] = {refFaceNormal_(0), refFaceNormal_(1), refFaceNormal_(2)};
189 
190  for(ordinal_type j=0; j <faceCardinality_; ++j) {
191  ordinal_type jdof = tagToOrdinal_(faceDim_, iface_, j);
192  for(ordinal_type d=0; d <dim_; ++d) {
193  ordinal_type dp1 = (d+1) % dim_;
194  ordinal_type dp2 = (d+2) % dim_;
195  for(ordinal_type iq=0; iq <numBasisGradEPoints; ++iq) {
196  faceBasisGradAtGradEPoints_(0,j,iq,d) = basisGradAtBasisGradEPoints_(jdof,offsetBasisGrad_+iq,dp1)*n[dp2] - basisGradAtBasisGradEPoints_(jdof,offsetBasisGrad_+iq,dp2)*n[dp1];
197  wBasisGradAtGradEPoints_(ic,j,iq,d) = faceBasisGradAtGradEPoints_(0,j,iq,d) * basisGradEWeights_(iq);
198  }
199  for(ordinal_type iq=0; iq <numTargetGradEPoints; ++iq) {
200  wBasisGradBasisAtTargetGradEPoints_(ic,j,iq,d) = (basisGradAtTargetGradEPoints_(jdof,offsetTargetGrad_+iq,dp1)*n[dp2] - basisGradAtTargetGradEPoints_(jdof,offsetTargetGrad_+iq,dp2)*n[dp1]) * targetGradEWeights_(iq);
201  }
202  }
203  }
204 
205  for(ordinal_type d=0; d <dim_; ++d) {
206  ordinal_type dp1 = (d+1) % dim_;
207  ordinal_type dp2 = (d+2) % dim_;
208  for(ordinal_type iq=0; iq <numTargetGradEPoints; ++iq)
209  targetGradTanAtTargetGradEPoints_(ic,iq,d) = targetGradAtTargetGradEPoints_(ic,offsetTargetGrad_+iq,dp1)*n[dp2] - targetGradAtTargetGradEPoints_(ic,offsetTargetGrad_+iq,dp2)*n[dp1];
210  }
211 
212  for(ordinal_type j=0; j <numVertexEdgeDofs_; ++j) {
213  ordinal_type jdof = computedDofs_(j);
214  for(ordinal_type d=0; d <dim_; ++d) {
215  ordinal_type dp1 = (d+1) % dim_;
216  ordinal_type dp2 = (d+2) % dim_;
217  for(ordinal_type iq=0; iq <numBasisGradEPoints; ++iq)
218  negPartialProjGrad_(ic,iq,d) -= (basisGradAtBasisGradEPoints_(jdof,offsetBasisGrad_+iq,dp1)*n[dp2] - basisGradAtBasisGradEPoints_(jdof,offsetBasisGrad_+iq,dp2)*n[dp1]) * basisCoeffs_(ic,jdof);
219  }
220  }
221  }
222 };
223 
224 
225 template<typename ViewType1, typename ViewType2, typename ViewType3,
226 typename ViewType4>
228  const ViewType1 basisCoeffs_;
229  const ViewType1 negPartialProjGrad_;
230  const ViewType1 cellBasisGradAtGradEPoints_;
231  const ViewType1 basisGradAtBasisGradEPoints_;
232  const ViewType2 basisGradEWeights_;
233  const ViewType1 wBasisGradAtGradEPoints_;
234  const ViewType2 targetGradEWeights_;
235  const ViewType1 basisGradAtTargetGradEPoints_;
236  const ViewType1 wBasisGradBasisAtTargetGradEPoints_;
237  const ViewType3 computedDofs_;
238  const ViewType4 elemDof_;
239  ordinal_type dim_;
240  ordinal_type numElemDofs_;
241  ordinal_type offsetBasisGrad_;
242  ordinal_type offsetTargetGrad_;
243  ordinal_type numVertexEdgeFaceDofs_;
244 
245  ComputeBasisCoeffsOnCells_HGRAD(const ViewType1 basisCoeffs, ViewType1 negPartialProjGrad, const ViewType1 cellBasisGradAtGradEPoints,
246  const ViewType1 basisGradAtBasisGradEPoints, const ViewType2 basisGradEWeights, const ViewType1 wBasisGradAtGradEPoints, const ViewType2 targetGradEWeights,
247  const ViewType1 basisGradAtTargetGradEPoints, const ViewType1 wBasisGradBasisAtTargetGradEPoints, const ViewType3 computedDofs, const ViewType4 elemDof,
248  ordinal_type dim, ordinal_type numElemDofs, ordinal_type offsetBasisGrad, ordinal_type offsetTargetGrad, ordinal_type numVertexEdgeFaceDofs) :
249  basisCoeffs_(basisCoeffs), negPartialProjGrad_(negPartialProjGrad), cellBasisGradAtGradEPoints_(cellBasisGradAtGradEPoints),
250  basisGradAtBasisGradEPoints_(basisGradAtBasisGradEPoints), basisGradEWeights_(basisGradEWeights), wBasisGradAtGradEPoints_(wBasisGradAtGradEPoints), targetGradEWeights_(targetGradEWeights),
251  basisGradAtTargetGradEPoints_(basisGradAtTargetGradEPoints), wBasisGradBasisAtTargetGradEPoints_(wBasisGradBasisAtTargetGradEPoints),
252  computedDofs_(computedDofs), elemDof_(elemDof), dim_(dim), numElemDofs_(numElemDofs), offsetBasisGrad_(offsetBasisGrad),
253  offsetTargetGrad_(offsetTargetGrad), numVertexEdgeFaceDofs_(numVertexEdgeFaceDofs) {}
254 
255  void
256  KOKKOS_INLINE_FUNCTION
257  operator()(const ordinal_type ic) const {
258  ordinal_type numBasisGradEPoints = basisGradEWeights_.extent(0);
259  ordinal_type numTargetGradEPoints = targetGradEWeights_.extent(0);
260  for(ordinal_type j=0; j <numElemDofs_; ++j) {
261  ordinal_type idof = elemDof_(j);
262  for(ordinal_type d=0; d <dim_; ++d) {
263  for(ordinal_type iq=0; iq <numBasisGradEPoints; ++iq) {
264  cellBasisGradAtGradEPoints_(0,j,iq,d) = basisGradAtBasisGradEPoints_(idof,offsetBasisGrad_+iq,d);
265  wBasisGradAtGradEPoints_(ic,j,iq,d) = cellBasisGradAtGradEPoints_(0,j,iq,d) * basisGradEWeights_(iq);
266  }
267  for(ordinal_type iq=0; iq <numTargetGradEPoints; ++iq) {
268  wBasisGradBasisAtTargetGradEPoints_(ic,j,iq,d )= basisGradAtTargetGradEPoints_(idof,offsetTargetGrad_+iq,d) * targetGradEWeights_(iq);
269  }
270  }
271  }
272  for(ordinal_type j=0; j <numVertexEdgeFaceDofs_; ++j) {
273  ordinal_type jdof = computedDofs_(j);
274  for(ordinal_type iq=0; iq <numBasisGradEPoints; ++iq)
275  for(ordinal_type d=0; d <dim_; ++d) {
276  negPartialProjGrad_(ic,iq,d) -= basisCoeffs_(ic,jdof)*basisGradAtBasisGradEPoints_(jdof,offsetBasisGrad_+iq,d);
277  }
278  }
279  }
280 };
281 
282 } // FunctorsProjectionTools namespace
283 
284 
285 template<typename DeviceType>
286 template<class BasisCoeffsViewType, class TargetValueViewType, class TargetGradViewType, class BasisType, class OrientationViewType>
287 void
288 ProjectionTools<DeviceType>::getHGradBasisCoeffs(BasisCoeffsViewType basisCoeffs,
289  const TargetValueViewType targetAtTargetEPoints,
290  const TargetGradViewType targetGradAtTargetGradEPoints,
291  const OrientationViewType orts,
292  const BasisType* cellBasis,
294 {
295  using funValsValueType = typename TargetValueViewType::value_type;
296  static_assert(std::is_same<funValsValueType,typename TargetGradViewType::value_type>::value,
297  "targetGradAtTargetGradEPoints and targetAtTargetEPoints must agree on their value type" );
298 
299  typedef typename BasisType::scalarType scalarType;
300  typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
301  typedef Kokkos::pair<ordinal_type,ordinal_type> range_type;
302  const auto cellTopo = cellBasis->getBaseCellTopology();
303  ordinal_type dim = cellTopo.getDimension();
304  ordinal_type basisCardinality = cellBasis->getCardinality();
305  ordinal_type numCells = targetAtTargetEPoints.extent(0);
306  const ordinal_type edgeDim = 1;
307  const ordinal_type faceDim = 2;
308 
309  ordinal_type numVertices = (cellBasis->getDofCount(0, 0) > 0) ? cellTopo.getVertexCount() : 0;
310  ordinal_type numEdges = (cellBasis->getDofCount(1, 0) > 0) ? cellTopo.getEdgeCount() : 0;
311  ordinal_type numFaces = (cellBasis->getDofCount(2, 0) > 0) ? cellTopo.getFaceCount() : 0;
312 
313  ScalarViewType refEdgeTan("refTan", dim);
314  ScalarViewType refFaceNormal("refNormal", dim);
315 
316  ordinal_type numVertexDofs = numVertices;
317 
318  ordinal_type numEdgeDofs(0);
319  for(ordinal_type ie=0; ie<numEdges; ++ie)
320  numEdgeDofs += cellBasis->getDofCount(edgeDim,ie);
321 
322  ordinal_type numFaceDofs(0);
323  for(ordinal_type iface=0; iface<numFaces; ++iface)
324  numFaceDofs += cellBasis->getDofCount(faceDim,iface);
325 
326  Kokkos::View<ordinal_type*, DeviceType> computedDofs("computedDofs",numVertexDofs+numEdgeDofs+numFaceDofs);
327 
328  ordinal_type numTotalBasisGradEPoints = projStruct->getNumBasisDerivEvalPoints();
329  auto basisGradEPoints = projStruct->getAllDerivEvalPoints(EvalPointsType::BASIS);
330 
331  ordinal_type numTotalTargetEPoints = projStruct->getNumTargetEvalPoints(), numTotalTargetGradEPoints = projStruct->getNumTargetDerivEvalPoints();
332  auto targetEPoints = projStruct->getAllEvalPoints(EvalPointsType::TARGET);
333  auto targetGradEPoints = projStruct->getAllDerivEvalPoints(EvalPointsType::TARGET);
334 
335 
336  ScalarViewType basisAtTargetEPoints("basisAtTargetEPoints",basisCardinality, numTotalTargetEPoints);
337  cellBasis->getValues(basisAtTargetEPoints, targetEPoints);
338 
339  ScalarViewType basisGradAtBasisGradEPoints;
340  ScalarViewType basisGradAtTargetGradEPoints;
341  if(numTotalBasisGradEPoints>0) {
342  basisGradAtBasisGradEPoints = ScalarViewType ("basisGradAtBasisGradEPoints",basisCardinality, numTotalBasisGradEPoints, dim);
343  basisGradAtTargetGradEPoints = ScalarViewType("basisGradAtTargetGradEPoints",basisCardinality, numTotalTargetGradEPoints, dim);
344 
345  cellBasis->getValues(basisGradAtBasisGradEPoints, basisGradEPoints, OPERATOR_GRAD);
346  cellBasis->getValues(basisGradAtTargetGradEPoints, targetGradEPoints, OPERATOR_GRAD);
347  }
348 
349  auto targetEPointsRange = Kokkos::create_mirror_view_and_copy(MemSpaceType(), projStruct->getTargetPointsRange());
350  auto targetGradEPointsRange = projStruct->getTargetDerivPointsRange();
351  auto basisGradEPointsRange = projStruct->getBasisDerivPointsRange();
352 
353  auto tagToOrdinal = Kokkos::create_mirror_view_and_copy(MemSpaceType(), cellBasis->getAllDofOrdinal());
354 
355  typename RefSubcellParametrization<DeviceType>::ConstViewType subcellParamFace;
356  if(numFaces>0)
357  subcellParamFace = RefSubcellParametrization<DeviceType>::get(faceDim, cellBasis->getBaseCellTopology().getKey());
358 
359 
360  Kokkos::parallel_for("Compute Dofs ", Kokkos::RangePolicy<ExecSpaceType, int> (0, numVertices),
361  KOKKOS_LAMBDA (const size_t iv) {
362  computedDofs(iv) = tagToOrdinal(0, iv, 0);
363  });
364  ordinal_type computedDofsCount = numVertices;
365 
366  ScalarViewType refBasisCoeffs("refBasisCoeffs", basisCoeffs.extent(0), basisCoeffs.extent(1));
367 
368  const Kokkos::RangePolicy<ExecSpaceType> policy(0, numCells);
369  using functorType = FunctorsProjectionTools::ComputeBasisCoeffsOnVertices_HGRAD<ScalarViewType, decltype(tagToOrdinal), decltype(targetEPointsRange),
370  decltype(targetAtTargetEPoints), decltype(basisAtTargetEPoints)>;
371  Kokkos::parallel_for(policy, functorType(refBasisCoeffs, tagToOrdinal, targetEPointsRange,
372  targetAtTargetEPoints, basisAtTargetEPoints, numVertices));
373 
374  for(ordinal_type ie=0; ie<numEdges; ++ie) {
375 
376  ordinal_type edgeCardinality = cellBasis->getDofCount(edgeDim,ie);
377  ordinal_type offsetBasis = basisGradEPointsRange(edgeDim, ie).first;
378  ordinal_type offsetTarget = targetGradEPointsRange(edgeDim, ie).first;
379  ordinal_type numBasisGradEPoints = range_size(basisGradEPointsRange(edgeDim, ie));
380  ordinal_type numTargetGradEPoints = range_size(targetGradEPointsRange(edgeDim, ie));
381  auto basisGradEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisDerivEvalWeights(edgeDim,ie));
382  auto targetGradEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetDerivEvalWeights(edgeDim,ie));
383 
384  CellTools<DeviceType>::getReferenceEdgeTangent(refEdgeTan, ie, cellTopo);
385 
386  ScalarViewType basisTanAtEPoints("basisTanAtEPoints",1,edgeCardinality, numBasisGradEPoints);
387  ScalarViewType targetGradTanAtTargetGradEPoints("tanBasisAtTargetGradEPoints",numCells, numTargetGradEPoints);
388  ScalarViewType wBasisAtBasisGradEPoints("wTanBasisAtBasisGradEPoints",numCells,edgeCardinality, numBasisGradEPoints);
389  ScalarViewType wBasisAtTargetGradEPoints("wTanBasisAtTargetGradEPoints",numCells,edgeCardinality, numTargetGradEPoints);
390  ScalarViewType negPartialProjGrad("negPartialProjGrad", numCells, numBasisGradEPoints);
391 
392  using functorTypeEdge = FunctorsProjectionTools::ComputeBasisCoeffsOnEdges_HGRAD<ScalarViewType, decltype(basisGradEWeights),
393  decltype(computedDofs), decltype(tagToOrdinal), decltype(targetGradAtTargetGradEPoints)>;
394 
395  Kokkos::parallel_for(policy, functorTypeEdge(refBasisCoeffs, negPartialProjGrad, basisTanAtEPoints,
396  basisGradAtBasisGradEPoints, basisGradEWeights, wBasisAtBasisGradEPoints, targetGradEWeights,
397  basisGradAtTargetGradEPoints, wBasisAtTargetGradEPoints, computedDofs, tagToOrdinal,
398  targetGradTanAtTargetGradEPoints, targetGradAtTargetGradEPoints, refEdgeTan,
399  edgeCardinality, offsetBasis,
400  offsetTarget, numVertexDofs, edgeDim, dim, ie));
401 
402  ScalarViewType edgeMassMat_("edgeMassMat_", 1, edgeCardinality, edgeCardinality),
403  edgeRhsMat_("rhsMat_", numCells, edgeCardinality);
404 
405  FunctionSpaceTools<DeviceType >::integrate(edgeMassMat_, basisTanAtEPoints, Kokkos::subview(wBasisAtBasisGradEPoints, std::make_pair(0,1), Kokkos::ALL(), Kokkos::ALL()) );
406  FunctionSpaceTools<DeviceType >::integrate(edgeRhsMat_, targetGradTanAtTargetGradEPoints, wBasisAtTargetGradEPoints);
407  FunctionSpaceTools<DeviceType >::integrate(edgeRhsMat_, negPartialProjGrad, wBasisAtBasisGradEPoints,true);
408 
409  typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
410  ScalarViewType t_("t",numCells, edgeCardinality);
411  WorkArrayViewType w_("w",numCells, edgeCardinality);
412 
413  auto edgeDofs = Kokkos::subview(tagToOrdinal, edgeDim, ie, Kokkos::ALL());
414 
415  ElemSystem edgeSystem("edgeSystem", true);
416  edgeSystem.solve(refBasisCoeffs, edgeMassMat_, edgeRhsMat_, t_, w_, edgeDofs, edgeCardinality);
417 
418  Kokkos::parallel_for("Compute Dofs ", Kokkos::RangePolicy<ExecSpaceType, int> (0, edgeCardinality),
419  KOKKOS_LAMBDA (const size_t i) {
420  computedDofs(computedDofsCount+i) = tagToOrdinal(edgeDim, ie, i);
421  });
422  computedDofsCount += edgeCardinality;
423 
424  }
425 
426  for(ordinal_type iface=0; iface<numFaces; ++iface) {
427 
428  ordinal_type faceCardinality = cellBasis->getDofCount(faceDim,iface);
429 
430  ordinal_type numGradEPoints = range_size(basisGradEPointsRange(faceDim, iface));
431  ordinal_type numTargetGradEPoints = range_size(targetGradEPointsRange(faceDim, iface));
432  ordinal_type offsetBasisGrad = basisGradEPointsRange(faceDim, iface).first;
433  ordinal_type offsetTargetGrad = targetGradEPointsRange(faceDim, iface).first;
434  auto basisGradEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisDerivEvalWeights(faceDim,iface));
435  auto targetGradEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetDerivEvalWeights(faceDim,iface));
436 
437  CellTools<DeviceType>::getReferenceFaceNormal(refFaceNormal, iface, cellTopo);
438 
439  ScalarViewType faceBasisGradAtGradEPoints("faceBasisGradAtGradEPoints",1,faceCardinality, numGradEPoints,dim);
440  ScalarViewType wBasisGradAtGradEPoints("wNormalBasisGradAtGradEPoints",numCells,faceCardinality, numGradEPoints,dim);
441  ScalarViewType wBasisGradBasisAtTargetGradEPoints("wNormalBasisGradAtTargetGradEPoints",numCells,faceCardinality, numTargetGradEPoints,dim);
442  ScalarViewType targetGradTanAtTargetGradEPoints("targetGradTanAtTargetGradEPoints",numCells, numTargetGradEPoints,dim);
443  ScalarViewType negPartialProjGrad("mNormalComputedProjection", numCells,numGradEPoints,dim);
444 
445  using functorTypeFace_HGRAD = FunctorsProjectionTools::ComputeBasisCoeffsOnFaces_HGRAD<ScalarViewType, decltype(basisGradEWeights),
446  decltype(computedDofs), decltype(tagToOrdinal), decltype(targetGradAtTargetGradEPoints)>;
447 
448  Kokkos::parallel_for(policy, functorTypeFace_HGRAD(refBasisCoeffs, negPartialProjGrad, faceBasisGradAtGradEPoints,
449  basisGradAtBasisGradEPoints, basisGradEWeights, wBasisGradAtGradEPoints, targetGradEWeights,
450  basisGradAtTargetGradEPoints,wBasisGradBasisAtTargetGradEPoints, computedDofs, tagToOrdinal,
451  targetGradTanAtTargetGradEPoints,targetGradAtTargetGradEPoints,
452  refFaceNormal, faceCardinality, offsetBasisGrad,
453  offsetTargetGrad, numVertexDofs+numEdgeDofs, numFaces, faceDim,
454  dim, iface));
455 
456  ScalarViewType faceMassMat_("faceMassMat_", 1, faceCardinality, faceCardinality),
457  faceRhsMat_("rhsMat_", numCells, faceCardinality);
458 
459  FunctionSpaceTools<DeviceType >::integrate(faceMassMat_, faceBasisGradAtGradEPoints, Kokkos::subview(wBasisGradAtGradEPoints, std::make_pair(0,1), Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL()) );
460  FunctionSpaceTools<DeviceType >::integrate(faceRhsMat_, targetGradTanAtTargetGradEPoints, wBasisGradBasisAtTargetGradEPoints);
461  FunctionSpaceTools<DeviceType >::integrate(faceRhsMat_, negPartialProjGrad, wBasisGradAtGradEPoints,true);
462 
463  typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
464  ScalarViewType t_("t",numCells, faceCardinality);
465  WorkArrayViewType w_("w",numCells, faceCardinality);
466 
467  auto faceDofs = Kokkos::subview(tagToOrdinal, faceDim, iface, Kokkos::ALL());
468 
469  ElemSystem faceSystem("faceSystem", true);
470  faceSystem.solve(refBasisCoeffs, faceMassMat_, faceRhsMat_, t_, w_, faceDofs, faceCardinality);
471 
472  Kokkos::parallel_for("Compute Face Dofs ", Kokkos::RangePolicy<ExecSpaceType, int> (0, faceCardinality),
473  KOKKOS_LAMBDA (const size_t i) {
474  computedDofs(computedDofsCount+i) = tagToOrdinal(faceDim, iface, i);
475  });
476  computedDofsCount += faceCardinality;
477  }
478 
479  ordinal_type numElemDofs = cellBasis->getDofCount(dim,0);
480  if(numElemDofs>0) {
481 
482  range_type cellTargetGradEPointsRange = targetGradEPointsRange(dim, 0);
483  ordinal_type numTargetGradEPoints = range_size(cellTargetGradEPointsRange);
484  ordinal_type numGradEPoints = range_size(basisGradEPointsRange(dim, 0));
485  ordinal_type offsetBasisGrad = basisGradEPointsRange(dim, 0).first;
486  ordinal_type offsetTargetGrad = cellTargetGradEPointsRange.first;
487  auto targetGradEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetDerivEvalWeights(dim,0));
488  auto basisGradEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisDerivEvalWeights(dim,0));
489 
490  ScalarViewType cellBasisGradAtGradEPoints("internalBasisGradAtEPoints",1,numElemDofs, numGradEPoints, dim);
491  ScalarViewType negPartialProjGrad("negPartialProjGrad", numCells, numGradEPoints, dim);
492  ScalarViewType wBasisGradAtGradEPoints("wBasisGradAtGradEPoints",numCells,numElemDofs, numGradEPoints,dim);
493  ScalarViewType wBasisGradBasisAtTargetGradEPoints("wBasisGradAtTargetGradEPoints",numCells,numElemDofs, numTargetGradEPoints,dim);
494 
495  auto elemDof = Kokkos::subview(tagToOrdinal, dim, 0, Kokkos::ALL());
497  Kokkos::parallel_for(policy, functorTypeCell_HGRAD(refBasisCoeffs, negPartialProjGrad, cellBasisGradAtGradEPoints,
498  basisGradAtBasisGradEPoints, basisGradEWeights, wBasisGradAtGradEPoints, targetGradEWeights,
499  basisGradAtTargetGradEPoints, wBasisGradBasisAtTargetGradEPoints, computedDofs, elemDof,
500  dim, numElemDofs, offsetBasisGrad, offsetTargetGrad, numVertexDofs+numEdgeDofs+numFaceDofs));
501  ScalarViewType cellMassMat_("cellMassMat_", 1, numElemDofs, numElemDofs),
502  cellRhsMat_("rhsMat_", numCells, numElemDofs);
503 
504  FunctionSpaceTools<DeviceType >::integrate(cellMassMat_, cellBasisGradAtGradEPoints, Kokkos::subview(wBasisGradAtGradEPoints, std::make_pair(0,1), Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL()) );
505  FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, Kokkos::subview(targetGradAtTargetGradEPoints,Kokkos::ALL(),cellTargetGradEPointsRange,Kokkos::ALL()), wBasisGradBasisAtTargetGradEPoints);
506  FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, negPartialProjGrad, wBasisGradAtGradEPoints, true);
507 
508  typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
509  ScalarViewType t_("t",numCells, numElemDofs);
510  WorkArrayViewType w_("w",numCells, numElemDofs);
511 
512  auto cellDofs = Kokkos::subview(tagToOrdinal, dim, 0, Kokkos::ALL());
513  ElemSystem cellSystem("cellSystem", true);
514  cellSystem.solve(refBasisCoeffs, cellMassMat_, cellRhsMat_, t_, w_, cellDofs, numElemDofs);
515  }
516 
517  OrientationTools<DeviceType>::modifyBasisByOrientationInverse(basisCoeffs, refBasisCoeffs, orts, cellBasis, true);
518 
519 }
520 
521 } // Intrepid2 namespace
522 
523 #endif
524 
static void getHGradBasisCoeffs(BasisCoeffsViewType basisCoeffs, const TargetValueViewType targetAtEvalPoints, const TargetGradViewType targetGradAtGradEvalPoints, const OrientationViewType cellOrientations, const BasisType *cellBasis, ProjectionStruct< DeviceType, typename BasisType::scalarType > *projStruct)
Computes the basis coefficients of the HGrad projection of the target function.
Header file for the abstract base class Intrepid2::DefaultCubatureFactory.
ordinal_type getNumBasisDerivEvalPoints()
Returns number of evaluation points for basis derivatives.
const range_tag getTargetDerivPointsRange() const
Returns the range tag of the target function derivative evaluation points on 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...
host_view_type getTargetDerivEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the function derivatives evaluation weights on a subcell.
void solve(ViewType1 basisCoeffs, ViewType2 elemMat, ViewType2 elemRhs, ViewType2 tau, ViewType3 w, const ViewType4 elemDof, ordinal_type n, ordinal_type m=0)
Solve the system and returns the basis coefficients solve the system either using Kokkos Kernel QR or...
static void getReferenceEdgeTangent(RefEdgeTangentViewType refEdgeTangent, const ordinal_type edgeOrd, const shards::CellTopology parentCell)
Computes constant tangent vectors to edges of 2D or 3D reference cells.
Header file for the Intrepid2::FunctionSpaceTools class.
view_type getAllDerivEvalPoints(EvalPointsType type=TARGET) const
Returns all the evaluation points for basis/target derivatives in the cell.
view_type getAllEvalPoints(EvalPointsType type=TARGET) const
Returns the basis/target evaluation points in the cell.
static void getReferenceFaceNormal(RefFaceNormalViewType refFaceNormal, const ordinal_type faceOrd, const shards::CellTopology parentCell)
Computes constant normal vectors to faces of 3D reference cell.
ordinal_type getNumTargetEvalPoints()
Returns number of points where to evaluate the target function.
const range_tag getBasisDerivPointsRange() const
Returns the range tag of the derivative evaluation points on subcell.
Header file for Intrepid2::ArrayTools class providing utilities for array operations.
static void integrate(Kokkos::DynRankView< outputValueValueType, outputValueProperties...> outputValues, const Kokkos::DynRankView< leftValueValueType, leftValueProperties...> leftValues, const Kokkos::DynRankView< rightValueValueType, rightValueProperties...> rightValues, const bool sumInto=false)
Contracts leftValues and rightValues arrays on the point and possibly space dimensions and stores the...
const range_tag getTargetPointsRange() const
Returns the range of the target function evaluation points on subcells.
host_view_type getBasisDerivEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the basis derivatives evaluation weights on a subcell.
Utility methods for Intrepid2 unit tests.
ordinal_type getNumTargetDerivEvalPoints()
Returns number of points where to evaluate the derivatives of the target function.
Class to solve a square system A x = b on each cell A is expected to be saddle a point (KKT) matrix o...
An helper class to compute the evaluation points and weights needed for performing projections...
static void modifyBasisByOrientationInverse(Kokkos::DynRankView< outputValueType, outputProperties...> output, const Kokkos::DynRankView< inputValueType, inputProperties...> input, const OrientationViewType orts, const BasisType *basis, const bool transpose=false)
Modify basis due to orientation, applying the inverse of the operator applied in modifyBasisByOrienta...