Intrepid2
Intrepid2_ProjectionToolsDefL2.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Intrepid2 Package
5 // Copyright (2007) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Kyungjoo Kim (kyukim@sandia.gov), or
38 // Mauro Perego (mperego@sandia.gov)
39 //
40 // ************************************************************************
41 // @HEADER
42 
49 #ifndef __INTREPID2_PROJECTIONTOOLSDEFL2_HPP__
50 #define __INTREPID2_PROJECTIONTOOLSDEFL2_HPP__
51 
52 #include "Intrepid2_DefaultCubatureFactory.hpp"
53 #include "Intrepid2_ArrayTools.hpp"
54 #include "Intrepid2_FunctionSpaceTools.hpp"
55 
56 
57 namespace Intrepid2 {
58 
59 namespace FunctorsProjectionTools {
60 
61 template<typename ViewType1, typename ViewType2, typename ViewType3,
62 typename ViewType4>
63 struct ComputeBasisCoeffsOnVertices_L2 {
64  ViewType1 basisCoeffs_;
65  const ViewType2 tagToOrdinal_;
66  const ViewType3 targetEPointsRange_;
67  const ViewType4 targetAtTargetEPoints_;
68  const ViewType1 basisAtTargetEPoints_;
69  ordinal_type numVertices_;
70 
71 
72  ComputeBasisCoeffsOnVertices_L2(ViewType1 basisCoeffs, ViewType2 tagToOrdinal, ViewType3 targetEPointsRange,
73  ViewType4 targetAtTargetEPoints, ViewType1 basisAtTargetEPoints, ordinal_type numVertices) :
74  basisCoeffs_(basisCoeffs), tagToOrdinal_(tagToOrdinal), targetEPointsRange_(targetEPointsRange),
75  targetAtTargetEPoints_(targetAtTargetEPoints), basisAtTargetEPoints_(basisAtTargetEPoints), numVertices_(numVertices) {}
76 
77  void
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;
83  //the value of the basis at the vertex might be different than 1; HGrad basis, so the function is scalar
84  basisCoeffs_(ic,idof) = targetAtTargetEPoints_(ic,pt)/basisAtTargetEPoints_(idof,pt,0);
85  }
86  }
87 };
88 
89 
90 template<typename ViewType1, typename ViewType2, typename ViewType3,
91 typename ViewType4, typename ViewType5>
92 struct ComputeBasisCoeffsOnEdges_L2 {
93  const ViewType1 basisCoeffs_;
94  const ViewType1 negPartialProj_;
95  const ViewType1 basisDofDofAtBasisEPoints_;
96  const ViewType1 basisAtBasisEPoints_;
97  const ViewType2 basisEWeights_;
98  const ViewType1 wBasisDofAtBasisEPoints_;
99  const ViewType2 targetEWeights_;
100  const ViewType1 basisAtTargetEPoints_;
101  const ViewType1 wBasisDofAtTargetEPoints_;
102  const ViewType3 computedDofs_;
103  const ViewType4 tagToOrdinal_;
104  const ViewType5 targetAtTargetEPoints_;
105  const ViewType1 targetTanAtTargetEPoints_;
106  const ViewType1 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_;
113  ordinal_type iedge_;
114 
115  ComputeBasisCoeffsOnEdges_L2(const ViewType1 basisCoeffs, ViewType1 negPartialProj, const ViewType1 basisDofDofAtBasisEPoints,
116  const ViewType1 basisAtBasisEPoints, const ViewType2 basisEWeights, const ViewType1 wBasisDofAtBasisEPoints, const ViewType2 targetEWeights,
117  const ViewType1 basisAtTargetEPoints, const ViewType1 wBasisDofAtTargetEPoints, const ViewType3 computedDofs, const ViewType4 tagToOrdinal,
118  const ViewType5 targetAtTargetEPoints, const ViewType1 targetTanAtTargetEPoints, const ViewType1 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)
128  {}
129 
130  void
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  typename ViewType1::value_type tmp =0;
137  for(ordinal_type d=0; d <fieldDim_; ++d)
138  tmp += basisAtBasisEPoints_(jdof,offsetBasis_+iq,d)*refEdgesVec_(iedge_,d);
139  basisDofDofAtBasisEPoints_(0,j,iq) = tmp;
140  wBasisDofAtBasisEPoints_(ic,j,iq) = tmp*basisEWeights_(iq);
141  }
142  for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
143  for(ordinal_type d=0; d <fieldDim_; ++d)
144  wBasisDofAtTargetEPoints_(ic,j,iq) += basisAtTargetEPoints_(jdof,offsetTarget_+iq,d)*refEdgesVec_(iedge_,d)*targetEWeights_(iq);
145  }
146  }
147 
148  for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq)
149  for(ordinal_type d=0; d <fieldDim_; ++d)
150  targetTanAtTargetEPoints_(ic,iq) += targetAtTargetEPoints_(ic,offsetTarget_+iq,d)*refEdgesVec_(iedge_,d);
151 
152  for(ordinal_type j=0; j <numVertexDofs_; ++j) {
153  ordinal_type jdof = computedDofs_(j);
154  for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq)
155  for(ordinal_type d=0; d <fieldDim_; ++d)
156  negPartialProj_(ic,iq) -= basisCoeffs_(ic,jdof)*basisAtBasisEPoints_(jdof,offsetBasis_+iq,d)*refEdgesVec_(iedge_,d);
157  }
158  }
159 };
160 
161 template<typename ViewType1, typename ViewType2, typename ViewType3,
162 typename ViewType4, typename ViewType5>
163 struct ComputeBasisCoeffsOnFaces_L2 {
164  const ViewType1 basisCoeffs_;
165  const ViewType1 negPartialProj_;
166  const ViewType1 faceBasisDofAtBasisEPoints_;
167  const ViewType1 basisAtBasisEPoints_;
168  const ViewType2 basisEWeights_;
169  const ViewType1 wBasisDofAtBasisEPoints_;
170  const ViewType2 targetEWeights_;
171  const ViewType1 basisAtTargetEPoints_;
172  const ViewType1 wBasisDofAtTargetEPoints_;
173  const ViewType3 computedDofs_;
174  const ViewType4 tagToOrdinal_;
175  const ViewType5 targetAtTargetEPoints_;
176  const ViewType1 targetDofAtTargetEPoints_;
177  const ViewType1 refSideNormal_;
178  ordinal_type fieldDim_;
179  ordinal_type faceCardinality_;
180  ordinal_type offsetBasis_;
181  ordinal_type offsetTarget_;
182  ordinal_type numVertexEdgeDofs_;
183  ordinal_type numFaces_;
184  ordinal_type faceDim_;
185  ordinal_type dim_;
186  ordinal_type iface_;
187  bool isHCurlBasis_, isHDivBasis_;
188 
189  ComputeBasisCoeffsOnFaces_L2(const ViewType1 basisCoeffs, ViewType1 negPartialProj, const ViewType1 faceBasisDofAtBasisEPoints,
190  const ViewType1 basisAtBasisEPoints, const ViewType2 basisEWeights, const ViewType1 wBasisDofAtBasisEPoints, const ViewType2 targetEWeights,
191  const ViewType1 basisAtTargetEPoints, const ViewType1 wBasisDofAtTargetEPoints, const ViewType3 computedDofs, const ViewType4 tagToOrdinal,
192  const ViewType5 targetAtTargetEPoints, const ViewType1 targetDofAtTargetEPoints,
193  const ViewType1 refSideNormal, 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, 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), targetAtTargetEPoints_(targetAtTargetEPoints),
200  targetDofAtTargetEPoints_(targetDofAtTargetEPoints), refSideNormal_(refSideNormal),
201  fieldDim_(fieldDim), faceCardinality_(faceCardinality), offsetBasis_(offsetBasis),
202  offsetTarget_(offsetTarget), numVertexEdgeDofs_(numVertexEdgeDofs), numFaces_(numFaces),
203  faceDim_(faceDim), dim_(dim), iface_(iface),
204  isHCurlBasis_(isHCurlBasis), isHDivBasis_(isHDivBasis)
205  {}
206 
207  void
208  KOKKOS_INLINE_FUNCTION
209  operator()(const ordinal_type ic) const {
210 
211  if(isHCurlBasis_) {
212  //normal
213  typename ViewType1::value_type n[3] = {refSideNormal_(0), refSideNormal_(1), refSideNormal_(2)};
214 
215  for(ordinal_type j=0; j <faceCardinality_; ++j) {
216  ordinal_type jdof = tagToOrdinal_(faceDim_, iface_, j);
217  for(ordinal_type d=0; d <fieldDim_; ++d) {
218  ordinal_type dp1 = (d+1) % dim_;
219  ordinal_type dp2 = (d+2) % dim_;
220  for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
221 
222  faceBasisDofAtBasisEPoints_(0,j,iq,d) = basisAtBasisEPoints_(jdof,offsetBasis_+iq,dp1)*n[dp2] - basisAtBasisEPoints_(jdof,offsetBasis_+iq,dp2)*n[dp1];
223  // basis \times n
224  wBasisDofAtBasisEPoints_(ic,j,iq,d) = faceBasisDofAtBasisEPoints_(0,j,iq,d) * basisEWeights_(iq);
225  }
226  for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
227  // basis \times n
228  wBasisDofAtTargetEPoints_(ic,j,iq,d) = (basisAtTargetEPoints_(jdof,offsetTarget_+iq,dp1)*n[dp2] - basisAtTargetEPoints_(jdof,offsetTarget_+iq,dp2)*n[dp1]) * targetEWeights_(iq);
229  }
230  }
231  }
232 
233  for(ordinal_type d=0; d <fieldDim_; ++d) {
234  ordinal_type dp1 = (d+1) % dim_;
235  ordinal_type dp2 = (d+2) % dim_;
236  for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
237  // target \times n
238  targetDofAtTargetEPoints_(ic,iq,d) = targetAtTargetEPoints_(ic,offsetTarget_+iq,dp1)*n[dp2] - targetAtTargetEPoints_(ic,offsetTarget_+iq,dp2)*n[dp1];
239  }
240  }
241 
242  for(ordinal_type j=0; j <numVertexEdgeDofs_; ++j) {
243  ordinal_type jdof = computedDofs_(j);
244  for(ordinal_type d=0; d <fieldDim_; ++d) {
245  ordinal_type dp1 = (d+1) % dim_;
246  ordinal_type dp2 = (d+2) % dim_;
247  for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
248  // basis \times n
249  negPartialProj_(ic,iq,d) -= (basisAtBasisEPoints_(jdof,offsetBasis_+iq,dp1)*n[dp2] - basisAtBasisEPoints_(jdof,offsetBasis_+iq,dp2)*n[dp1])*basisCoeffs_(ic,jdof);
250  }
251  }
252  }
253  } else { // isHDivBasis_ || isHGradBasis_
254  typename ViewType1::value_type coeff[3] = {1,0,0}; //only need first component for HGrad basis
255  if (isHDivBasis_) {
256  for (ordinal_type d=0; d<3; d++)
257  coeff[d] = refSideNormal_(d);
258  }
259 
260  for(ordinal_type j=0; j <faceCardinality_; ++j) {
261  ordinal_type jdof = tagToOrdinal_(faceDim_, iface_, j);
262  for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
263  typename ViewType1::value_type tmp =0;
264  for(ordinal_type d=0; d <fieldDim_; ++d)
265  tmp += coeff[d]*basisAtBasisEPoints_(jdof,offsetBasis_+iq,d);
266  faceBasisDofAtBasisEPoints_(0,j,iq,0) = tmp;
267  wBasisDofAtBasisEPoints_(ic,j,iq,0) = tmp * basisEWeights_(iq);
268  }
269  for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
270  typename ViewType2::value_type sum=0;
271  for(ordinal_type d=0; d <fieldDim_; ++d)
272  sum += coeff[d]*basisAtTargetEPoints_(jdof,offsetTarget_+iq,d);
273  wBasisDofAtTargetEPoints_(ic,j,iq,0) = sum * targetEWeights_(iq);
274  }
275  }
276 
277  for(ordinal_type d=0; d <fieldDim_; ++d)
278  for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq)
279  targetDofAtTargetEPoints_(ic,iq,0) += coeff[d]*targetAtTargetEPoints_(ic,offsetTarget_+iq,d);
280 
281  for(ordinal_type j=0; j <numVertexEdgeDofs_; ++j) {
282  ordinal_type jdof = computedDofs_(j);
283  for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq)
284  for(ordinal_type d=0; d <fieldDim_; ++d)
285  negPartialProj_(ic,iq,0) -= basisCoeffs_(ic,jdof)*coeff[d]*basisAtBasisEPoints_(jdof,offsetBasis_+iq,d);
286  }
287  }
288  }
289 };
290 
291 
292 template<typename ViewType1, typename ViewType2, typename ViewType3, typename ViewType4>
293 struct ComputeBasisCoeffsOnCells_L2 {
294  const ViewType1 basisCoeffs_;
295  const ViewType1 negPartialProj_;
296  const ViewType1 internalBasisAtBasisEPoints_;
297  const ViewType1 basisAtBasisEPoints_;
298  const ViewType2 basisEWeights_;
299  const ViewType1 wBasisAtBasisEPoints_;
300  const ViewType2 targetEWeights_;
301  const ViewType1 basisAtTargetEPoints_;
302  const ViewType1 wBasisDofAtTargetEPoints_;
303  const ViewType3 computedDofs_;
304  const ViewType4 elemDof_;
305  ordinal_type fieldDim_;
306  ordinal_type numElemDofs_;
307  ordinal_type offsetBasis_;
308  ordinal_type offsetTarget_;
309  ordinal_type numVertexEdgeFaceDofs_;
310 
311  ComputeBasisCoeffsOnCells_L2(const ViewType1 basisCoeffs, ViewType1 negPartialProj, const ViewType1 internalBasisAtBasisEPoints,
312  const ViewType1 basisAtBasisEPoints, const ViewType2 basisEWeights, const ViewType1 wBasisAtBasisEPoints, const ViewType2 targetEWeights,
313  const ViewType1 basisAtTargetEPoints, const ViewType1 wBasisDofAtTargetEPoints, const ViewType3 computedDofs, const ViewType4 elemDof,
314  ordinal_type fieldDim, ordinal_type numElemDofs, ordinal_type offsetBasis, ordinal_type offsetTarget, ordinal_type numVertexEdgeFaceDofs) :
315  basisCoeffs_(basisCoeffs), negPartialProj_(negPartialProj), internalBasisAtBasisEPoints_(internalBasisAtBasisEPoints),
316  basisAtBasisEPoints_(basisAtBasisEPoints), basisEWeights_(basisEWeights), wBasisAtBasisEPoints_(wBasisAtBasisEPoints), targetEWeights_(targetEWeights),
317  basisAtTargetEPoints_(basisAtTargetEPoints), wBasisDofAtTargetEPoints_(wBasisDofAtTargetEPoints),
318  computedDofs_(computedDofs), elemDof_(elemDof), fieldDim_(fieldDim), numElemDofs_(numElemDofs), offsetBasis_(offsetBasis),
319  offsetTarget_(offsetTarget), numVertexEdgeFaceDofs_(numVertexEdgeFaceDofs) {}
320 
321  void
322  KOKKOS_INLINE_FUNCTION
323  operator()(const ordinal_type ic) const {
324 
325  for(ordinal_type j=0; j <numElemDofs_; ++j) {
326  ordinal_type idof = elemDof_(j);
327  for(ordinal_type d=0; d <fieldDim_; ++d) {
328  for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
329  internalBasisAtBasisEPoints_(0,j,iq,d) = basisAtBasisEPoints_(idof,offsetBasis_+iq,d);
330  wBasisAtBasisEPoints_(ic,j,iq,d) = internalBasisAtBasisEPoints_(0,j,iq,d) * basisEWeights_(iq);
331  }
332  for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
333  wBasisDofAtTargetEPoints_(ic,j,iq,d) = basisAtTargetEPoints_(idof,offsetTarget_+iq,d)* targetEWeights_(iq);
334  }
335  }
336  }
337  for(ordinal_type j=0; j < numVertexEdgeFaceDofs_; ++j) {
338  ordinal_type jdof = computedDofs_(j);
339  for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq)
340  for(ordinal_type d=0; d <fieldDim_; ++d) {
341  negPartialProj_(ic,iq,d) -= basisCoeffs_(ic,jdof)*basisAtBasisEPoints_(jdof,offsetBasis_+iq,d);
342  }
343  }
344  }
345 };
346 
347 template<typename ViewType1, typename ViewType2>
348 struct MultiplyBasisByWeights {
349  const ViewType1 basisAtBasisEPoints_;
350  const ViewType2 basisEWeights_;
351  const ViewType1 wBasisAtBasisEPoints_;
352  const ViewType2 targetEWeights_;
353  const ViewType1 basisAtTargetEPoints_;
354  const ViewType1 wBasisDofAtTargetEPoints_;
355  ordinal_type fieldDim_;
356  ordinal_type numElemDofs_;
357 
358  MultiplyBasisByWeights(const ViewType1 basisAtBasisEPoints, const ViewType2 basisEWeights, const ViewType1 wBasisAtBasisEPoints, const ViewType2 targetEWeights,
359  const ViewType1 basisAtTargetEPoints, const ViewType1 wBasisDofAtTargetEPoints,
360  ordinal_type fieldDim, ordinal_type numElemDofs) :
361  basisAtBasisEPoints_(basisAtBasisEPoints), basisEWeights_(basisEWeights), wBasisAtBasisEPoints_(wBasisAtBasisEPoints), targetEWeights_(targetEWeights),
362  basisAtTargetEPoints_(basisAtTargetEPoints), wBasisDofAtTargetEPoints_(wBasisDofAtTargetEPoints),
363  fieldDim_(fieldDim), numElemDofs_(numElemDofs) {}
364 
365  void
366  KOKKOS_INLINE_FUNCTION
367  operator()(const ordinal_type ic) const {
368 
369  for(ordinal_type j=0; j <numElemDofs_; ++j) {
370  for(ordinal_type d=0; d <fieldDim_; ++d) {
371  for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
372  wBasisAtBasisEPoints_(ic,j,iq,d) = basisAtBasisEPoints_(j,iq,d) * basisEWeights_(iq);
373  }
374  for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
375  wBasisDofAtTargetEPoints_(ic,j,iq,d) = basisAtTargetEPoints_(j,iq,d)* targetEWeights_(iq);
376  }
377  }
378  }
379  }
380 };
381 
382 } // FunctorsProjectionTools namespace
383 
384 
385 template<typename DeviceType>
386 template<typename basisCoeffsValueType, class ...basisCoeffsProperties,
387 typename funValsValueType, class ...funValsProperties,
388 typename BasisType,
389 typename ortValueType,class ...ortProperties>
390 void
391 ProjectionTools<DeviceType>::getL2BasisCoeffs(Kokkos::DynRankView<basisCoeffsValueType,basisCoeffsProperties...> basisCoeffs,
392  const Kokkos::DynRankView<funValsValueType,funValsProperties...> targetAtTargetEPoints,
393  const Kokkos::DynRankView<ortValueType, ortProperties...> orts,
394  const BasisType* cellBasis,
395  ProjectionStruct<DeviceType, typename BasisType::scalarType> * projStruct){
396 
397  typedef typename BasisType::scalarType scalarType;
398  typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
399  typedef Kokkos::pair<ordinal_type,ordinal_type> range_type;
400  const auto cellTopo = cellBasis->getBaseCellTopology();
401  ordinal_type dim = cellTopo.getDimension();
402  ordinal_type basisCardinality = cellBasis->getCardinality();
403  ordinal_type numCells = targetAtTargetEPoints.extent(0);
404  const ordinal_type edgeDim = 1;
405  const ordinal_type faceDim = 2;
406  const ordinal_type fieldDim = (targetAtTargetEPoints.rank()==2) ? 1 : targetAtTargetEPoints.extent(2);
407 
408  ordinal_type numVertices = (cellBasis->getDofCount(0, 0) > 0) ? cellTopo.getVertexCount() : 0;
409  ordinal_type numEdges = (cellBasis->getDofCount(1, 0) > 0) ? cellTopo.getEdgeCount() : 0;
410  ordinal_type numFaces = (cellBasis->getDofCount(2, 0) > 0) ? cellTopo.getFaceCount() : 0;
411 
412  ScalarViewType refEdgesVec("refEdgesVec", numEdges, dim);
413  ScalarViewType refFacesTangents("refFaceTangents", numFaces, dim, 2);
414  ScalarViewType refFacesNormal("refFaceNormal", numFaces, dim);
415 
416  ordinal_type numVertexDofs = numVertices;
417 
418  ordinal_type numEdgeDofs(0);
419  for(ordinal_type ie=0; ie<numEdges; ++ie)
420  numEdgeDofs += cellBasis->getDofCount(edgeDim,ie);
421 
422  ordinal_type numFaceDofs(0);
423  for(ordinal_type iface=0; iface<numFaces; ++iface)
424  numFaceDofs += cellBasis->getDofCount(faceDim,iface);
425 
426  Kokkos::View<ordinal_type*, DeviceType> computedDofs("computedDofs", numVertexDofs+numEdgeDofs+numFaceDofs);
427 
428  auto basisEPointsRange = projStruct->getBasisPointsRange();
429 
430  ordinal_type numTotalBasisEPoints = projStruct->getNumBasisEvalPoints();
431  auto basisEPoints = projStruct->getAllEvalPoints(EvalPointsType::BASIS);
432 
433  ordinal_type numTotalTargetEPoints = projStruct->getNumTargetEvalPoints();
434  auto targetEPoints = projStruct->getAllEvalPoints(EvalPointsType::TARGET);
435 
436  auto tagToOrdinal = Kokkos::create_mirror_view_and_copy(MemSpaceType(), cellBasis->getAllDofOrdinal());
437 
438  ScalarViewType basisAtBasisEPoints("basisAtBasisEPoints",basisCardinality, numTotalBasisEPoints, fieldDim);
439  ScalarViewType basisAtTargetEPoints("basisAtTargetEPoints",basisCardinality, numTotalTargetEPoints, fieldDim);
440  {
441  if(fieldDim == 1) {
442  cellBasis->getValues(Kokkos::subview(basisAtTargetEPoints, Kokkos::ALL(), Kokkos::ALL(), 0), targetEPoints);
443  cellBasis->getValues(Kokkos::subview(basisAtBasisEPoints, Kokkos::ALL(), Kokkos::ALL(), 0), basisEPoints);
444  }
445  else {
446  cellBasis->getValues(basisAtTargetEPoints, targetEPoints);
447  cellBasis->getValues(basisAtBasisEPoints, basisEPoints);
448  }
449  }
450 
451  {
452  auto hostComputedDofs = Kokkos::create_mirror_view(computedDofs);
453  ordinal_type computedDofsCount = 0;
454  for(ordinal_type iv=0; iv<numVertices; ++iv)
455  hostComputedDofs(computedDofsCount++) = cellBasis->getDofOrdinal(0, iv, 0);
456 
457  for(ordinal_type ie=0; ie<numEdges; ++ie) {
458  ordinal_type edgeCardinality = cellBasis->getDofCount(edgeDim,ie);
459  for(ordinal_type i=0; i<edgeCardinality; ++i)
460  hostComputedDofs(computedDofsCount++) = cellBasis->getDofOrdinal(edgeDim, ie, i);
461  }
462 
463  for(ordinal_type iface=0; iface<numFaces; ++iface) {
464  ordinal_type faceCardinality = cellBasis->getDofCount(faceDim,iface);
465  for(ordinal_type i=0; i<faceCardinality; ++i)
466  hostComputedDofs(computedDofsCount++) = cellBasis->getDofOrdinal(faceDim, iface, i);
467  }
468  Kokkos::deep_copy(computedDofs,hostComputedDofs);
469  }
470 
471  bool isHGradBasis = (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HGRAD);
472  bool isHCurlBasis = (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HCURL);
473  bool isHDivBasis = (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HDIV);
474  ordinal_type faceDofDim = isHCurlBasis ? 3 : 1;
475  ScalarViewType edgeCoeff("edgeCoeff", fieldDim);
476 
477 
478  const Kokkos::RangePolicy<ExecSpaceType> policy(0, numCells);
479  ScalarViewType refBasisCoeffs("refBasisCoeffs", basisCoeffs.extent(0), basisCoeffs.extent(1));
480 
481  if(isHGradBasis) {
482 
483  auto targetEPointsRange = Kokkos::create_mirror_view_and_copy(MemSpaceType(), projStruct->getTargetPointsRange());
484  using functorType = FunctorsProjectionTools::ComputeBasisCoeffsOnVertices_L2<ScalarViewType, decltype(tagToOrdinal), decltype(targetEPointsRange),
485  decltype(targetAtTargetEPoints)>;
486  Kokkos::parallel_for(policy, functorType(refBasisCoeffs, tagToOrdinal, targetEPointsRange,
487  targetAtTargetEPoints, basisAtTargetEPoints, numVertices));
488  }
489 
490  auto targetEPointsRange = projStruct->getTargetPointsRange();
491  for(ordinal_type ie=0; ie<numEdges; ++ie) {
492  auto edgeVec = Kokkos::subview(refEdgesVec, ie, Kokkos::ALL());
493  //auto edgeVecHost = Kokkos::create_mirror_view(edgeVec);
494 
495  if(isHCurlBasis) {
496  CellTools<DeviceType>::getReferenceEdgeTangent(edgeVec, ie, cellTopo);
497  } else if(isHDivBasis) {
498  CellTools<DeviceType>::getReferenceSideNormal(edgeVec, ie, cellTopo);
499  } else {
500  deep_copy(edgeVec, 1.0);
501  }
502 
503  ordinal_type edgeCardinality = cellBasis->getDofCount(edgeDim,ie);
504  ordinal_type numBasisEPoints = range_size(basisEPointsRange(edgeDim, ie));
505  ordinal_type numTargetEPoints = range_size(targetEPointsRange(edgeDim, ie));
506 
507  ScalarViewType basisDofAtBasisEPoints("BasisDofAtBasisEPoints",1,edgeCardinality, numBasisEPoints);
508  ScalarViewType tragetDofAtTargetEPoints("TargetDofAtTargetEPoints",numCells, numTargetEPoints);
509  ScalarViewType weightedBasisAtBasisEPoints("weightedTanBasisAtBasisEPoints",numCells,edgeCardinality, numBasisEPoints);
510  ScalarViewType weightedBasisAtTargetEPoints("weightedTanBasisAtTargetEPoints",numCells,edgeCardinality, numTargetEPoints);
511  ScalarViewType negPartialProj("negPartialProj", numCells, numBasisEPoints);
512 
513  auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(edgeDim,ie));
514  auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(edgeDim,ie));
515 
516  //Note: we are not considering the jacobian of the orientation map since it is simply a scalar term for the integrals and it does not affect the projection
517  ordinal_type offsetBasis = basisEPointsRange(edgeDim, ie).first;
518  ordinal_type offsetTarget = targetEPointsRange(edgeDim, ie).first;
519 
520 
521  using functorTypeEdge = FunctorsProjectionTools::ComputeBasisCoeffsOnEdges_L2<ScalarViewType, decltype(basisEWeights),
522  decltype(computedDofs), decltype(tagToOrdinal), decltype(targetAtTargetEPoints)>;
523 
524  Kokkos::parallel_for(policy, functorTypeEdge(refBasisCoeffs, negPartialProj, basisDofAtBasisEPoints,
525  basisAtBasisEPoints, basisEWeights, weightedBasisAtBasisEPoints, targetEWeights,
526  basisAtTargetEPoints, weightedBasisAtTargetEPoints, computedDofs, tagToOrdinal,
527  targetAtTargetEPoints,tragetDofAtTargetEPoints, refEdgesVec, fieldDim,
528  edgeCardinality, offsetBasis, offsetTarget, numVertexDofs, edgeDim, ie));
529 
530 
531  ScalarViewType edgeMassMat_("edgeMassMat_", 1, edgeCardinality, edgeCardinality),
532  edgeRhsMat_("rhsMat_", numCells, edgeCardinality);
533 
534  FunctionSpaceTools<DeviceType >::integrate(edgeMassMat_, basisDofAtBasisEPoints, Kokkos::subview(weightedBasisAtBasisEPoints, std::make_pair(0,1), Kokkos::ALL(),Kokkos::ALL()) );
535  FunctionSpaceTools<DeviceType >::integrate(edgeRhsMat_, tragetDofAtTargetEPoints, weightedBasisAtTargetEPoints);
536  FunctionSpaceTools<DeviceType >::integrate(edgeRhsMat_, negPartialProj, weightedBasisAtBasisEPoints,true);
537 
538 
539  typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
540  ScalarViewType t_("t",numCells, edgeCardinality);
541  WorkArrayViewType w_("w",numCells, edgeCardinality);
542 
543  auto edgeDof = Kokkos::subview(tagToOrdinal, edgeDim, ie, Kokkos::ALL());
544 
545  ElemSystem edgeSystem("edgeSystem", true);
546  edgeSystem.solve(refBasisCoeffs, edgeMassMat_, edgeRhsMat_, t_, w_, edgeDof, edgeCardinality);
547  }
548 
549  for(ordinal_type iface=0; iface<numFaces; ++iface) {
550  ordinal_type faceCardinality = cellBasis->getDofCount(faceDim,iface);
551 
552  ordinal_type numTargetEPoints = range_size(targetEPointsRange(faceDim, iface));
553  ordinal_type numBasisEPoints = range_size(basisEPointsRange(faceDim, iface));
554 
555  ScalarViewType faceBasisDofAtBasisEPoints("faceBasisDofAtBasisEPoints",1,faceCardinality, numBasisEPoints,faceDofDim);
556  ScalarViewType wBasisDofAtBasisEPoints("weightedBasisDofAtBasisEPoints",numCells,faceCardinality, numBasisEPoints,faceDofDim);
557 
558  ScalarViewType faceBasisAtTargetEPoints("faceBasisDofAtTargetEPoints",numCells,faceCardinality, numTargetEPoints,faceDofDim);
559  ScalarViewType wBasisDofAtTargetEPoints("weightedBasisDofAtTargetEPoints",numCells,faceCardinality, numTargetEPoints,faceDofDim);
560 
561  ScalarViewType targetDofAtTargetEPoints("targetDofAtTargetEPoints",numCells, numTargetEPoints,faceDofDim);
562  ScalarViewType negPartialProj("negComputedProjection", numCells,numBasisEPoints,faceDofDim);
563 
564  ordinal_type offsetBasis = basisEPointsRange(faceDim, iface).first;
565  ordinal_type offsetTarget = targetEPointsRange(faceDim, iface).first;
566  auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(faceDim,iface));
567  auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(faceDim,iface));
568 
569  ScalarViewType refSideNormal("refSideNormal", dim);
570  CellTools<DeviceType>::getReferenceSideNormal(refSideNormal, iface, cellTopo);
571 
572  using functorTypeFace = FunctorsProjectionTools::ComputeBasisCoeffsOnFaces_L2<ScalarViewType, decltype(basisEWeights),
573  decltype(computedDofs), decltype(tagToOrdinal), decltype(targetAtTargetEPoints)>;
574 
575  Kokkos::parallel_for(policy, functorTypeFace(refBasisCoeffs, negPartialProj,faceBasisDofAtBasisEPoints,
576  basisAtBasisEPoints, basisEWeights, wBasisDofAtBasisEPoints, targetEWeights,
577  basisAtTargetEPoints, wBasisDofAtTargetEPoints, computedDofs, tagToOrdinal,
578  targetAtTargetEPoints,targetDofAtTargetEPoints,
579  refSideNormal, fieldDim, faceCardinality, offsetBasis,
580  offsetTarget, numVertexDofs+numEdgeDofs, numFaces, faceDim,
581  dim, iface, isHCurlBasis, isHDivBasis));
582 
583  typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
584  ScalarViewType faceMassMat_("faceMassMat_", 1, faceCardinality, faceCardinality),
585  faceRhsMat_("rhsMat_", numCells, faceCardinality);
586 
587  FunctionSpaceTools<DeviceType >::integrate(faceMassMat_, faceBasisDofAtBasisEPoints, Kokkos::subview(wBasisDofAtBasisEPoints, std::make_pair(0,1), Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL()) );
588  FunctionSpaceTools<DeviceType >::integrate(faceRhsMat_, targetDofAtTargetEPoints, wBasisDofAtTargetEPoints);
589  FunctionSpaceTools<DeviceType >::integrate(faceRhsMat_, negPartialProj, wBasisDofAtBasisEPoints,true);
590 
591  ScalarViewType t_("t",numCells, faceCardinality);
592  WorkArrayViewType w_("w",numCells,faceCardinality);
593 
594  auto faceDof = Kokkos::subview(tagToOrdinal, faceDim, iface, Kokkos::ALL());
595 
596  ElemSystem faceSystem("faceSystem", true);
597  faceSystem.solve(refBasisCoeffs, faceMassMat_, faceRhsMat_, t_, w_, faceDof, faceCardinality);
598  }
599 
600  ordinal_type numElemDofs = cellBasis->getDofCount(dim,0);
601 
602 
603  if(numElemDofs>0) {
604 
605  auto cellDofs = Kokkos::subview(tagToOrdinal, dim, 0, Kokkos::ALL());
606 
607  range_type cellPointsRange = targetEPointsRange(dim, 0);
608 
609  ordinal_type numTargetEPoints = range_size(targetEPointsRange(dim,0));
610  ordinal_type numBasisEPoints = range_size(basisEPointsRange(dim,0));
611 
612  ScalarViewType internalBasisAtBasisEPoints("internalBasisAtBasisEPoints",1,numElemDofs, numBasisEPoints, fieldDim);
613  ScalarViewType negPartialProj("negPartialProj", numCells, numBasisEPoints, fieldDim);
614 
615  auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(dim,0));
616  auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(dim,0));
617  ordinal_type offsetBasis = basisEPointsRange(dim, 0).first;
618  ordinal_type offsetTarget = targetEPointsRange(dim, 0).first;
619 
620  ScalarViewType wBasisAtBasisEPoints("weightedBasisAtBasisEPoints",numCells,numElemDofs, numBasisEPoints,fieldDim);
621  ScalarViewType wBasisDofAtTargetEPoints("weightedBasisAtTargetEPoints",numCells,numElemDofs, numTargetEPoints,fieldDim);
622 
623  using functorType = FunctorsProjectionTools::ComputeBasisCoeffsOnCells_L2<ScalarViewType, decltype(basisEWeights), decltype(computedDofs), decltype(cellDofs)>;
624  Kokkos::parallel_for(policy, functorType( refBasisCoeffs, negPartialProj, internalBasisAtBasisEPoints,
625  basisAtBasisEPoints, basisEWeights, wBasisAtBasisEPoints, targetEWeights, basisAtTargetEPoints, wBasisDofAtTargetEPoints,
626  computedDofs, cellDofs, fieldDim, numElemDofs, offsetBasis, offsetTarget, numVertexDofs+numEdgeDofs+numFaceDofs));
627 
628  typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
629  ScalarViewType cellMassMat_("cellMassMat_", 1, numElemDofs, numElemDofs),
630  cellRhsMat_("rhsMat_", numCells, numElemDofs);
631 
632  FunctionSpaceTools<DeviceType >::integrate(cellMassMat_, internalBasisAtBasisEPoints, Kokkos::subview(wBasisAtBasisEPoints, std::make_pair(0,1), Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL()) );
633  if(fieldDim==1)
634  FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, Kokkos::subview(targetAtTargetEPoints,Kokkos::ALL(),cellPointsRange,Kokkos::ALL()),
635  Kokkos::subview(wBasisDofAtTargetEPoints,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),0));
636  else
637  FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, Kokkos::subview(targetAtTargetEPoints,Kokkos::ALL(),cellPointsRange,Kokkos::ALL()), wBasisDofAtTargetEPoints);
638  FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, negPartialProj, wBasisAtBasisEPoints, true);
639 
640  ScalarViewType t_("t",numCells, numElemDofs);
641  WorkArrayViewType w_("w",numCells,numElemDofs);
642  ElemSystem cellSystem("cellSystem", true);
643  cellSystem.solve(refBasisCoeffs, cellMassMat_, cellRhsMat_, t_, w_, cellDofs, numElemDofs);
644  }
645 
646  OrientationTools<DeviceType>::modifyBasisByOrientationInverse(basisCoeffs, refBasisCoeffs, orts, cellBasis, true);
647 }
648 
649 
650 template<typename DeviceType>
651 template<typename basisCoeffsValueType, class ...basisCoeffsProperties,
652 typename funValsValueType, class ...funValsProperties,
653 typename BasisType,
654 typename ortValueType,class ...ortProperties>
655 void
656 ProjectionTools<DeviceType>::getL2DGBasisCoeffs(Kokkos::DynRankView<basisCoeffsValueType,basisCoeffsProperties...> basisCoeffs,
657  const Kokkos::DynRankView<funValsValueType,funValsProperties...> targetAtTargetEPoints,
658  const Kokkos::DynRankView<ortValueType, ortProperties...> orts,
659  const BasisType* cellBasis,
660  ProjectionStruct<DeviceType, typename BasisType::scalarType> * projStruct){
661 
662  Kokkos::DynRankView<typename BasisType::scalarType,DeviceType> refBasisCoeffs("refBasisCoeffs", basisCoeffs.extent(0), basisCoeffs.extent(1));
663  getL2DGBasisCoeffs(refBasisCoeffs, targetAtTargetEPoints, cellBasis, projStruct);
664 
665  OrientationTools<DeviceType>::modifyBasisByOrientationInverse(basisCoeffs, refBasisCoeffs, orts, cellBasis, true);
666 }
667 
668 
669 template<typename DeviceType>
670 template<typename basisViewType, typename targetViewType, typename BasisType>
671 void
672 ProjectionTools<DeviceType>::getL2DGBasisCoeffs(basisViewType basisCoeffs,
673  const targetViewType targetAtTargetEPoints,
674  const BasisType* cellBasis,
675  ProjectionStruct<DeviceType, typename BasisType::scalarType> * projStruct){
676 
677  typedef typename BasisType::scalarType scalarType;
678  typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
679  const auto cellTopo = cellBasis->getBaseCellTopology();
680 
681  ordinal_type dim = cellTopo.getDimension();
682 
683  auto basisEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),
684  projStruct->getEvalPoints(dim,0,EvalPointsType::BASIS));
685  auto targetEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),
686  projStruct->getEvalPoints(dim,0,EvalPointsType::TARGET));
687 
688  ordinal_type numTotalTargetEPoints(targetAtTargetEPoints.extent(1));
689  ordinal_type basisCardinality = cellBasis->getCardinality();
690  ordinal_type numCells = targetAtTargetEPoints.extent(0);
691  const ordinal_type fieldDim = (targetAtTargetEPoints.rank()==2) ? 1 : targetAtTargetEPoints.extent(2);
692 
693  ordinal_type numTotalBasisEPoints = projStruct->getNumBasisEvalPoints();
694  ScalarViewType basisAtBasisEPoints("basisAtBasisEPoints",1,basisCardinality, numTotalBasisEPoints, fieldDim);
695  ScalarViewType basisAtTargetEPoints("basisAtTargetEPoints",basisCardinality, numTotalTargetEPoints, fieldDim);
696  {
697  if(fieldDim == 1) {
698  cellBasis->getValues(Kokkos::subview(basisAtTargetEPoints,Kokkos::ALL(),Kokkos::ALL(),0), targetEPoints);
699  cellBasis->getValues(Kokkos::subview(basisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL(),0), basisEPoints);
700  }
701  else {
702  cellBasis->getValues(basisAtTargetEPoints, targetEPoints);
703  cellBasis->getValues(Kokkos::subview(basisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), basisEPoints);
704  }
705  }
706 
707  const Kokkos::RangePolicy<ExecSpaceType> policy(0, numCells);
708  ordinal_type numElemDofs = cellBasis->getCardinality();
709 
710  auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(dim,0));
711  auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(dim,0));
712 
713  ScalarViewType wBasisAtBasisEPoints("weightedBasisAtBasisEPoints",1,numElemDofs, numTotalBasisEPoints,fieldDim);
714  ScalarViewType wBasisDofAtTargetEPoints("weightedBasisAtTargetEPoints",numCells,numElemDofs, numTotalTargetEPoints,fieldDim);
715  Kokkos::DynRankView<ordinal_type, DeviceType> cellDofs("cellDoFs", numElemDofs);
716 
717  Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,numElemDofs),
718  KOKKOS_LAMBDA (const int &j) {
719  for(ordinal_type d=0; d <fieldDim; ++d) {
720  for(ordinal_type iq=0; iq < numTotalBasisEPoints; ++iq)
721  wBasisAtBasisEPoints(0,j,iq,d) = basisAtBasisEPoints(0,j,iq,d) * basisEWeights(iq);
722  for(ordinal_type iq=0; iq <numTotalTargetEPoints; ++iq) {
723  wBasisDofAtTargetEPoints(0,j,iq,d) = basisAtTargetEPoints(j,iq,d)* targetEWeights(iq);
724  }
725  }
726  cellDofs(j) = j;
727  });
728  RealSpaceTools<DeviceType>::clone(wBasisDofAtTargetEPoints, Kokkos::subview(wBasisDofAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()));
729 
730  typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
731  ScalarViewType cellMassMat_("cellMassMat_", 1, numElemDofs, numElemDofs),
732  cellRhsMat_("rhsMat_", numCells, numElemDofs);
733 
734  FunctionSpaceTools<DeviceType >::integrate(cellMassMat_, basisAtBasisEPoints, wBasisAtBasisEPoints);
735  if(fieldDim==1)
736  FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, targetAtTargetEPoints,
737  Kokkos::subview(wBasisDofAtTargetEPoints,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),0));
738  else
739  FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, targetAtTargetEPoints, wBasisDofAtTargetEPoints);
740 
741  ScalarViewType t_("t",1, numElemDofs);
742  WorkArrayViewType w_("w",numCells,numElemDofs);
743  ElemSystem cellSystem("cellSystem", true);
744  cellSystem.solve(basisCoeffs, cellMassMat_, cellRhsMat_, t_, w_, cellDofs, numElemDofs);
745 }
746 
747 } //Intrepid2 namespace
748 
749 #endif
750 
Intrepid2_DefaultCubatureFactory.hpp
Header file for the abstract base class Intrepid2::DefaultCubatureFactory.
Intrepid2::FunctorsProjectionTools::ComputeBasisCoeffsOnEdges_L2
Definition: Intrepid2_ProjectionToolsDefL2.hpp:92
Intrepid2::ProjectionStruct::getBasisPointsRange
const range_tag getBasisPointsRange() const
Returns the range tag of the basis evaluation points subcells.
Definition: Intrepid2_ProjectionStruct.hpp:377
Intrepid2::ProjectionStruct::getEvalPoints
host_view_type getEvalPoints(const ordinal_type subCellDim, const ordinal_type subCellId, EvalPointsType type) const
Returns the basis/target evaluation points on a subcell.
Definition: Intrepid2_ProjectionStruct.hpp:242
Intrepid2::ProjectionStruct::getNumBasisEvalPoints
ordinal_type getNumBasisEvalPoints()
Returns number of basis evaluation points.
Definition: Intrepid2_ProjectionStruct.hpp:112
Intrepid2::ProjectionTools::ElemSystem::solve
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...
Definition: Intrepid2_ProjectionTools.hpp:516
Intrepid2::RealSpaceTools::clone
static void clone(Kokkos::DynRankView< outputValueType, outputProperties...> output, const Kokkos::DynRankView< inputValueType, inputProperties...> input)
Clone input array.
Definition: Intrepid2_RealSpaceToolsDef.hpp:334
Intrepid2::FunctorsProjectionTools::ComputeBasisCoeffsOnVertices_L2
Definition: Intrepid2_ProjectionToolsDefL2.hpp:63
Intrepid2::ProjectionTools::getL2DGBasisCoeffs
static void getL2DGBasisCoeffs(Kokkos::DynRankView< basisCoeffsValueType, basisCoeffsProperties...> basisCoeffs, const Kokkos::DynRankView< funValsValueType, funValsProperties...> targetAtEvalPoints, const Kokkos::DynRankView< ortValueType, ortProperties...> cellOrientations, const BasisType *cellBasis, ProjectionStruct< DeviceType, typename BasisType::scalarType > *projStruct)
Computes evaluation points for local L2 projection for broken HGRAD HCURL HDIV and HVOL spaces...
Definition: Intrepid2_ProjectionToolsDefL2.hpp:656
Intrepid2::FunctorsProjectionTools::ComputeBasisCoeffsOnCells_L2
Definition: Intrepid2_ProjectionToolsDefL2.hpp:293
Intrepid2::CellTools::getReferenceEdgeTangent
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.
Definition: Intrepid2_CellToolsDefNodeInfo.hpp:244
Intrepid2::FunctorsProjectionTools::ComputeBasisCoeffsOnFaces_L2
Definition: Intrepid2_ProjectionToolsDefL2.hpp:163
Intrepid2_FunctionSpaceTools.hpp
Header file for the Intrepid2::FunctionSpaceTools class.
Intrepid2::CellTools::getReferenceSideNormal
static void getReferenceSideNormal(RefSideNormalViewType refSideNormal, const ordinal_type sideOrd, const shards::CellTopology parentCell)
Computes constant normal vectors to sides of 2D or 3D reference cells.
Definition: Intrepid2_CellToolsDefNodeInfo.hpp:327
Intrepid2::ProjectionStruct::getAllEvalPoints
view_type getAllEvalPoints(EvalPointsType type=TARGET) const
Returns the basis/target evaluation points in the cell.
Definition: Intrepid2_ProjectionStruct.hpp:261
Intrepid2::ProjectionStruct::getBasisEvalWeights
host_view_type getBasisEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the basis evaluation weights on a subcell.
Definition: Intrepid2_ProjectionStruct.hpp:320
Intrepid2::ProjectionStruct::getNumTargetEvalPoints
ordinal_type getNumTargetEvalPoints()
Returns number of points where to evaluate the target function.
Definition: Intrepid2_ProjectionStruct.hpp:124
Intrepid2_ArrayTools.hpp
Header file for Intrepid2::ArrayTools class providing utilities for array operations.
Intrepid2::FunctionSpaceTools::integrate
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...
Definition: Intrepid2_FunctionSpaceToolsDef.hpp:426
Intrepid2::ProjectionStruct::getTargetPointsRange
const range_tag getTargetPointsRange() const
Returns the range of the target function evaluation points on subcells.
Definition: Intrepid2_ProjectionStruct.hpp:421
Intrepid2::ProjectionStruct::getTargetEvalWeights
host_view_type getTargetEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the function evaluation weights on a subcell.
Definition: Intrepid2_ProjectionStruct.hpp:352
Intrepid2::ProjectionTools::getL2BasisCoeffs
static void getL2BasisCoeffs(Kokkos::DynRankView< basisCoeffsValueType, basisCoeffsProperties...> basisCoeffs, const Kokkos::DynRankView< funValsValueType, funValsProperties...> targetAtEvalPoints, const Kokkos::DynRankView< ortValueType, ortProperties...> cellOrientations, const BasisType *cellBasis, ProjectionStruct< DeviceType, typename BasisType::scalarType > *projStruct)
Computes the basis coefficients of the L2 projection of the target function.
Definition: Intrepid2_ProjectionToolsDefL2.hpp:391
Intrepid2::ProjectionTools::ElemSystem
Class to solve a square system A x = b on each cell A is expected to be saddle a point (KKT) matrix o...
Definition: Intrepid2_ProjectionTools.hpp:472
Intrepid2::ProjectionStruct
An helper class to compute the evaluation points and weights needed for performing projections...
Definition: Intrepid2_ProjectionStruct.hpp:88
Intrepid2::OrientationTools::modifyBasisByOrientationInverse
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...
Definition: Intrepid2_OrientationToolsDefModifyBasis.hpp:383
Generated by   doxygen 1.8.5