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
48
49#ifndef __INTREPID2_PROJECTIONTOOLSDEFL2_HPP__
50#define __INTREPID2_PROJECTIONTOOLSDEFL2_HPP__
51
55
56
57namespace Intrepid2 {
58namespace Experimental {
59
60
61template<typename ViewType1, typename ViewType2, typename ViewType3,
62typename ViewType4, typename ViewType5>
63struct ComputeBasisCoeffsOnVertices_L2 {
64 ViewType1 basisCoeffs_;
65 const ViewType2 tagToOrdinal_;
66 const ViewType3 targetEPointsRange_;
67 const ViewType4 targetAtTargetEPoints_;
68 const ViewType5 basisAtTargetEPoints_;
69 ordinal_type numVertices_;
70
71
72 ComputeBasisCoeffsOnVertices_L2(ViewType1 basisCoeffs, ViewType2 tagToOrdinal, ViewType3 targetEPointsRange,
73 ViewType4 targetAtTargetEPoints, ViewType5 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_(ic,idof,pt,0);
85 }
86 }
87};
88
89
90template<typename ViewType1, typename ViewType2, typename ViewType3,
91typename ViewType4, typename ViewType5, typename ViewType6>
92struct ComputeBasisCoeffsOnEdges_L2 {
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_;
113 ordinal_type iedge_;
114
115 ComputeBasisCoeffsOnEdges_L2(const ViewType1 basisCoeffs, ViewType2 negPartialProj, const ViewType2 basisDofDofAtBasisEPoints,
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)
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 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);
139 }
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);
143 }
144 }
145
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);
149
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);
155 }
156 }
157};
158
159template<typename ViewType1, typename ViewType2, typename ViewType3,
160typename ViewType4, typename ViewType5, typename ViewType6, typename ViewType7, typename ViewType8>
161struct ComputeBasisCoeffsOnFaces_L2 {
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_;
184 ordinal_type dim_;
185 ordinal_type iface_;
186 unsigned topoKey_;
187 bool isHCurlBasis_, isHDivBasis_;
188
189 ComputeBasisCoeffsOnFaces_L2(const ViewType1 basisCoeffs, ViewType2 negPartialProj, const ViewType2 faceBasisDofAtBasisEPoints,
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)
205 {}
206
207 void
208 KOKKOS_INLINE_FUNCTION
209 operator()(const ordinal_type ic) const {
210
211 ordinal_type fOrt[6];
212 orts_(ic).getFaceOrientation(fOrt, numFaces_);
213 ordinal_type ort = fOrt[iface_];
214
215 if(isHCurlBasis_) {
216 typename ViewType3::value_type data[2*3];
217 auto tangents = ViewType3(data, 2, dim_);
218 Impl::OrientationTools::getRefSubcellTangents(tangents,faceParametrization_,topoKey_,iface_,ort);
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) {
223 tmp[0] = tmp[1] = 0;
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);
227 }
228 // u \times n = t_0 (u \dot t_1) - t1 (u \dot t_0)
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);
232 }
233 }
234 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
235 tmp[0] = tmp[1] = 0;
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);
239 }
240 // u \times n = t_0 (u \dot t_1) - t1 (u \dot t_0)
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);
243 }
244 }
245
246 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
247 tmp[0] = tmp[1] = 0;
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);
251 }
252 // u \times n = t_0 (u \dot t_1) - t1 (u \dot t_0)
253 for(ordinal_type d=0; d <fieldDim_; ++d)
254 targetDofAtTargetEPoints_(ic,iq,d) = (tangents(0,d)*tmp[1]-tangents(1,d)*tmp[0]);
255 }
256
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) {
260 tmp[0] = tmp[1] = 0;
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);
264 }
265
266 // u \times n = t_0 (u \dot t_1) - t1 (u \dot t_0)
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);
269 }
270 }
271 } else {
272 typename ViewType3::value_type coeff[3];
273 if (isHDivBasis_) {
274 typename ViewType3::value_type data[3*3];
275 auto tangentsAndNormal = ViewType3(data, dim_, dim_);
276 Impl::OrientationTools::getRefSideTangentsAndNormal(tangentsAndNormal, faceParametrization_,topoKey_, iface_, ort);
277 for(ordinal_type d=0; d <dim_; ++d)
278 coeff[d] = tangentsAndNormal(dim_-1,d);
279 }
280 else //isHGradBasis
281 coeff[0] = 1;
282
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);
290 }
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);
296 }
297 }
298
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);
302
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);
308 }
309 }
310 }
311};
312
313
314template<typename ViewType1, typename ViewType2, typename ViewType3,
315typename ViewType4, typename ViewType5>
316struct ComputeBasisCoeffsOnCells_L2 {
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_;
333
334 ComputeBasisCoeffsOnCells_L2(const ViewType1 basisCoeffs, ViewType2 negPartialProj, const ViewType2 internalBasisAtBasisEPoints,
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) {}
343
344 void
345 KOKKOS_INLINE_FUNCTION
346 operator()(const ordinal_type ic) const {
347
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);
354 }
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);
357 }
358 }
359 }
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);
365 }
366 }
367 }
368};
369
370
371template<typename DeviceType>
372template<typename BasisType,
373typename ortValueType, class ...ortProperties>
374void
375ProjectionTools<DeviceType>::getL2EvaluationPoints(typename BasisType::ScalarViewType ePoints,
376 const Kokkos::DynRankView<ortValueType, ortProperties...> orts,
377 const BasisType* cellBasis,
379 const EvalPointsType ePointType) {
380 const auto cellTopo = cellBasis->getBaseCellTopology();
381 //const auto cellTopoKey = cellBasis->getBaseCellTopology().getKey();
382 ordinal_type dim = cellTopo.getDimension();
383 ordinal_type numCells = ePoints.extent(0);
384 const ordinal_type edgeDim = 1;
385 const ordinal_type faceDim = 2;
386
387 ordinal_type numVertices = (cellBasis->getDofCount(0, 0) > 0) ? cellTopo.getVertexCount() : 0;
388 ordinal_type numEdges = (cellBasis->getDofCount(edgeDim, 0) > 0) ? cellTopo.getEdgeCount() : 0;
389 ordinal_type numFaces = (cellBasis->getDofCount(faceDim, 0) > 0) ? cellTopo.getFaceCount() : 0;
390 ordinal_type numVols = (cellBasis->getDofCount(dim, 0) > 0);
391
392 auto ePointsRange = projStruct->getPointsRange(ePointType);
393
394 typename RefSubcellParametrization<DeviceType>::ConstViewType subcellParamEdge, subcellParamFace;
395 if(numEdges>0)
396 subcellParamEdge = RefSubcellParametrization<DeviceType>::get(edgeDim, cellTopo.getKey());
397 if(numFaces>0)
398 subcellParamFace = RefSubcellParametrization<DeviceType>::get(faceDim, cellTopo.getKey());
399
400 auto refTopologyKey = projStruct->getTopologyKey();
401
402 if(numVertices>0) {
403 for(ordinal_type iv=0; iv<numVertices; ++iv) {
404 auto vertexEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getEvalPoints(0,iv,ePointType));
405 RealSpaceTools<DeviceType>::clone(Kokkos::subview(ePoints, Kokkos::ALL(),
406 ePointsRange(0, iv), Kokkos::ALL()), vertexEPoints);
407 }
408 }
409
410 for(ordinal_type ie=0; ie<numEdges; ++ie) {
411 auto edgePointsRange = ePointsRange(edgeDim, ie);
412 auto edgeEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getEvalPoints(edgeDim,ie,ePointType));
413
414 const auto topoKey = refTopologyKey(edgeDim,ie);
415 Kokkos::parallel_for
416 ("Evaluate Points",
417 Kokkos::RangePolicy<ExecSpaceType, int> (0, numCells),
418 KOKKOS_LAMBDA (const size_t ic) {
419
420 ordinal_type eOrt[12];
421 orts(ic).getEdgeOrientation(eOrt, numEdges);
422 ordinal_type ort = eOrt[ie];
423
424 Impl::OrientationTools::mapSubcellCoordsToRefCell(Kokkos::subview(ePoints,ic,edgePointsRange,Kokkos::ALL()),
425 edgeEPoints, subcellParamEdge, topoKey, ie, ort);
426 });
427 }
428
429 for(ordinal_type iface=0; iface<numFaces; ++iface) {
430 auto facePointsRange = ePointsRange(faceDim, iface);
431 auto faceEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getEvalPoints(faceDim,iface,ePointType));
432
433 const auto topoKey = refTopologyKey(faceDim,iface);
434 Kokkos::parallel_for
435 ("Evaluate Points",
436 Kokkos::RangePolicy<ExecSpaceType, int> (0, numCells),
437 KOKKOS_LAMBDA (const size_t ic) {
438 ordinal_type fOrt[6];
439 orts(ic).getFaceOrientation(fOrt, numFaces);
440 ordinal_type ort = fOrt[iface];
441
442 Impl::OrientationTools::mapSubcellCoordsToRefCell(Kokkos::subview(ePoints, ic, facePointsRange, Kokkos::ALL()),
443 faceEPoints, subcellParamFace, topoKey, iface, ort);
444 });
445 }
446
447
448 if(numVols > 0) {
449 auto pointsRange = ePointsRange(dim, 0);
450 auto cellEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getEvalPoints(dim,0,ePointType));
451 if(dim == 3)
452 RealSpaceTools<DeviceType>::clone(Kokkos::subview(ePoints, Kokkos::ALL(), pointsRange, Kokkos::ALL()), cellEPoints);
453 else { //if the cell is a side also internal points need orientation
454 const auto topoKey = refTopologyKey(dim,0);
455 RealSpaceTools<DeviceType>::clone(Kokkos::subview(ePoints, Kokkos::ALL(), pointsRange, Kokkos::ALL()), cellEPoints);
456 Kokkos::parallel_for
457 ("Evaluate Points",
458 Kokkos::RangePolicy<ExecSpaceType, int> (0, numCells),
459 KOKKOS_LAMBDA (const size_t ic) {
460 ordinal_type ort = 0;
461 if(dim == 1)
462 orts(ic).getEdgeOrientation(&ort,1);
463 else if (dim == 2)
464 orts(ic).getFaceOrientation(&ort,1);
465 Impl::OrientationTools::mapToModifiedReference(Kokkos::subview(ePoints, ic, pointsRange, Kokkos::ALL()), cellEPoints,topoKey,ort);
466 });
467 }
468 }
469}
470
471template<typename DeviceType>
472template<typename basisCoeffsValueType, class ...basisCoeffsProperties,
473typename funValsValueType, class ...funValsProperties,
474typename BasisType,
475typename ortValueType,class ...ortProperties>
476void
477ProjectionTools<DeviceType>::getL2BasisCoeffs(Kokkos::DynRankView<basisCoeffsValueType,basisCoeffsProperties...> basisCoeffs,
478 const Kokkos::DynRankView<funValsValueType,funValsProperties...> targetAtTargetEPoints,
479 const typename BasisType::ScalarViewType targetEPoints,
480 const Kokkos::DynRankView<ortValueType, ortProperties...> orts,
481 const BasisType* cellBasis,
483
484 typedef typename BasisType::scalarType scalarType;
485 typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
486 typedef Kokkos::pair<ordinal_type,ordinal_type> range_type;
487 const auto cellTopo = cellBasis->getBaseCellTopology();
488 ordinal_type dim = cellTopo.getDimension();
489 ordinal_type numTotalTargetEPoints(targetAtTargetEPoints.extent(1));
490 ordinal_type basisCardinality = cellBasis->getCardinality();
491 ordinal_type numCells = targetAtTargetEPoints.extent(0);
492 const ordinal_type edgeDim = 1;
493 const ordinal_type faceDim = 2;
494 const ordinal_type fieldDim = (targetAtTargetEPoints.rank()==2) ? 1 : targetAtTargetEPoints.extent(2);
495
496 ordinal_type numVertices = (cellBasis->getDofCount(0, 0) > 0) ? cellTopo.getVertexCount() : 0;
497 ordinal_type numEdges = (cellBasis->getDofCount(1, 0) > 0) ? cellTopo.getEdgeCount() : 0;
498 ordinal_type numFaces = (cellBasis->getDofCount(2, 0) > 0) ? cellTopo.getFaceCount() : 0;
499
500 ScalarViewType refEdgesVec("refEdgesVec", numEdges, dim);
501 ScalarViewType refFacesTangents("refFaceTangents", numFaces, dim, 2);
502 ScalarViewType refFacesNormal("refFaceNormal", numFaces, dim);
503
504 ordinal_type numVertexDofs = numVertices;
505
506 ordinal_type numEdgeDofs(0);
507 for(ordinal_type ie=0; ie<numEdges; ++ie)
508 numEdgeDofs += cellBasis->getDofCount(edgeDim,ie);
509
510 ordinal_type numFaceDofs(0);
511 for(ordinal_type iface=0; iface<numFaces; ++iface)
512 numFaceDofs += cellBasis->getDofCount(faceDim,iface);
513
514 Kokkos::View<ordinal_type*, DeviceType> computedDofs("computedDofs", numVertexDofs+numEdgeDofs+numFaceDofs);
515
516 auto basisEPointsRange = projStruct->getBasisPointsRange();
517
518 auto refTopologyKey = projStruct->getTopologyKey();
519
520 ordinal_type numTotalBasisEPoints = projStruct->getNumBasisEvalPoints();
521 ScalarViewType basisEPoints("basisEPoints",numCells,numTotalBasisEPoints, dim);
522 getL2EvaluationPoints(basisEPoints, orts, cellBasis, projStruct, EvalPointsType::BASIS);
523
524 auto tagToOrdinal = Kokkos::create_mirror_view_and_copy(MemSpaceType(), cellBasis->getAllDofOrdinal());
525
526 ScalarViewType basisAtBasisEPoints("basisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints, fieldDim);
527 ScalarViewType basisAtTargetEPoints("basisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints, fieldDim);
528 {
529 if(fieldDim == 1) {
530 ScalarViewType nonOrientedBasisAtBasisEPoints("nonOrientedBasisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints);
531 ScalarViewType nonOrientedBasisAtTargetEPoints("nonOrientedBasisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints);
532 for(ordinal_type ic=0; ic<numCells; ++ic) {
533 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtTargetEPoints,ic,Kokkos::ALL(),Kokkos::ALL()), Kokkos::subview(targetEPoints, ic, Kokkos::ALL(), Kokkos::ALL()));
534 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtBasisEPoints,ic,Kokkos::ALL(),Kokkos::ALL()), Kokkos::subview(basisEPoints, ic, Kokkos::ALL(), Kokkos::ALL()));
535 }
536 OrientationTools<DeviceType>::modifyBasisByOrientation(Kokkos::subview(basisAtBasisEPoints, Kokkos::ALL(), Kokkos::ALL(),
537 Kokkos::ALL(),0), nonOrientedBasisAtBasisEPoints, orts, cellBasis);
538 OrientationTools<DeviceType>::modifyBasisByOrientation(Kokkos::subview(basisAtTargetEPoints, Kokkos::ALL(),
539 Kokkos::ALL(), Kokkos::ALL(),0), nonOrientedBasisAtTargetEPoints, orts, cellBasis);
540 }
541 else {
542 ScalarViewType nonOrientedBasisAtBasisEPoints("nonOrientedBasisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints,fieldDim);
543 ScalarViewType nonOrientedBasisAtTargetEPoints("nonOrientedBasisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints,fieldDim);
544 for(ordinal_type ic=0; ic<numCells; ++ic) {
545 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtTargetEPoints,ic,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), Kokkos::subview(targetEPoints, ic, Kokkos::ALL(), Kokkos::ALL()));
546 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtBasisEPoints,ic,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), Kokkos::subview(basisEPoints, ic, Kokkos::ALL(), Kokkos::ALL()));
547 }
548 OrientationTools<DeviceType>::modifyBasisByOrientation(basisAtBasisEPoints, nonOrientedBasisAtBasisEPoints, orts, cellBasis);
549 OrientationTools<DeviceType>::modifyBasisByOrientation(basisAtTargetEPoints, nonOrientedBasisAtTargetEPoints, orts, cellBasis);
550 }
551 }
552
553 {
554 auto hostComputedDofs = Kokkos::create_mirror_view(computedDofs);
555 ordinal_type computedDofsCount = 0;
556 for(ordinal_type iv=0; iv<numVertices; ++iv)
557 hostComputedDofs(computedDofsCount++) = cellBasis->getDofOrdinal(0, iv, 0);
558
559 for(ordinal_type ie=0; ie<numEdges; ++ie) {
560 ordinal_type edgeCardinality = cellBasis->getDofCount(edgeDim,ie);
561 for(ordinal_type i=0; i<edgeCardinality; ++i)
562 hostComputedDofs(computedDofsCount++) = cellBasis->getDofOrdinal(edgeDim, ie, i);
563 }
564
565 for(ordinal_type iface=0; iface<numFaces; ++iface) {
566 ordinal_type faceCardinality = cellBasis->getDofCount(faceDim,iface);
567 for(ordinal_type i=0; i<faceCardinality; ++i)
568 hostComputedDofs(computedDofsCount++) = cellBasis->getDofOrdinal(faceDim, iface, i);
569 }
570 Kokkos::deep_copy(computedDofs,hostComputedDofs);
571 }
572
573 bool isHGradBasis = (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HGRAD);
574 bool isHCurlBasis = (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HCURL);
575 bool isHDivBasis = (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HDIV);
576 ordinal_type faceDofDim = isHCurlBasis ? 3 : 1;
577 ScalarViewType edgeCoeff("edgeCoeff", fieldDim);
578
579
580 const Kokkos::RangePolicy<ExecSpaceType> policy(0, numCells);
581
582 if(isHGradBasis) {
583
584 auto targetEPointsRange = Kokkos::create_mirror_view_and_copy(MemSpaceType(), projStruct->getTargetPointsRange());
585 typedef ComputeBasisCoeffsOnVertices_L2<decltype(basisCoeffs), decltype(tagToOrdinal), decltype(targetEPointsRange),
586 decltype(targetAtTargetEPoints), decltype(basisAtTargetEPoints)> functorType;
587 Kokkos::parallel_for(policy, functorType(basisCoeffs, tagToOrdinal, targetEPointsRange,
588 targetAtTargetEPoints, basisAtTargetEPoints, numVertices));
589 }
590
591 auto targetEPointsRange = projStruct->getTargetPointsRange();
592 for(ordinal_type ie=0; ie<numEdges; ++ie) {
593 auto edgeVec = Kokkos::subview(refEdgesVec, ie, Kokkos::ALL());
594 //auto edgeVecHost = Kokkos::create_mirror_view(edgeVec);
595
596 if(isHCurlBasis) {
598 } else if(isHDivBasis) {
600 } else {
601 deep_copy(edgeVec, 1.0);
602 }
603
604 ordinal_type edgeCardinality = cellBasis->getDofCount(edgeDim,ie);
605 ordinal_type numBasisEPoints = range_size(basisEPointsRange(edgeDim, ie));
606 ordinal_type numTargetEPoints = range_size(targetEPointsRange(edgeDim, ie));
607
608 ScalarViewType basisDofAtBasisEPoints("BasisDofAtBasisEPoints",numCells,edgeCardinality, numBasisEPoints);
609 ScalarViewType tragetDofAtTargetEPoints("TargetDofAtTargetEPoints",numCells, numTargetEPoints);
610 ScalarViewType weightedBasisAtBasisEPoints("weightedTanBasisAtBasisEPoints",numCells,edgeCardinality, numBasisEPoints);
611 ScalarViewType weightedBasisAtTargetEPoints("weightedTanBasisAtTargetEPoints",numCells,edgeCardinality, numTargetEPoints);
612 ScalarViewType negPartialProj("negPartialProj", numCells, numBasisEPoints);
613
614 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(edgeDim,ie));
615 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(edgeDim,ie));
616
617 //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
618 ordinal_type offsetBasis = basisEPointsRange(edgeDim, ie).first;
619 ordinal_type offsetTarget = targetEPointsRange(edgeDim, ie).first;
620
621
622 typedef ComputeBasisCoeffsOnEdges_L2<decltype(basisCoeffs), ScalarViewType, decltype(basisEWeights),
623 decltype(computedDofs), decltype(tagToOrdinal), decltype(targetAtTargetEPoints)> functorTypeEdge;
624
625 Kokkos::parallel_for(policy, functorTypeEdge(basisCoeffs, negPartialProj, basisDofAtBasisEPoints,
626 basisAtBasisEPoints, basisEWeights, weightedBasisAtBasisEPoints, targetEWeights,
627 basisAtTargetEPoints, weightedBasisAtTargetEPoints, computedDofs, tagToOrdinal,
628 targetAtTargetEPoints,tragetDofAtTargetEPoints, refEdgesVec, fieldDim,
629 edgeCardinality, offsetBasis, offsetTarget, numVertexDofs, edgeDim, ie));
630
631
632 ScalarViewType edgeMassMat_("edgeMassMat_", numCells, edgeCardinality, edgeCardinality),
633 edgeRhsMat_("rhsMat_", numCells, edgeCardinality);
634
635 FunctionSpaceTools<DeviceType >::integrate(edgeMassMat_, basisDofAtBasisEPoints, weightedBasisAtBasisEPoints);
636 FunctionSpaceTools<DeviceType >::integrate(edgeRhsMat_, tragetDofAtTargetEPoints, weightedBasisAtTargetEPoints);
637 FunctionSpaceTools<DeviceType >::integrate(edgeRhsMat_, negPartialProj, weightedBasisAtBasisEPoints,true);
638
639
640 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
641 ScalarViewType t_("t",numCells, edgeCardinality);
642 WorkArrayViewType w_("w",numCells, edgeCardinality);
643
644 auto edgeDof = Kokkos::subview(tagToOrdinal, edgeDim, ie, Kokkos::ALL());
645
646 ElemSystem edgeSystem("edgeSystem", false);
647 edgeSystem.solve(basisCoeffs, edgeMassMat_, edgeRhsMat_, t_, w_, edgeDof, edgeCardinality);
648 }
649
650 typename RefSubcellParametrization<DeviceType>::ConstViewType subcellParamFace;
651 if(numFaces>0)
652 subcellParamFace = RefSubcellParametrization<DeviceType>::get(faceDim, cellBasis->getBaseCellTopology().getKey());
653
654 for(ordinal_type iface=0; iface<numFaces; ++iface) {
655 const auto topoKey = refTopologyKey(faceDim,iface);
656 ordinal_type faceCardinality = cellBasis->getDofCount(faceDim,iface);
657
658 ordinal_type numTargetEPoints = range_size(targetEPointsRange(faceDim, iface));
659 ordinal_type numBasisEPoints = range_size(basisEPointsRange(faceDim, iface));
660
661 ScalarViewType faceBasisDofAtBasisEPoints("faceBasisDofAtBasisEPoints",numCells,faceCardinality, numBasisEPoints,faceDofDim);
662 ScalarViewType wBasisDofAtBasisEPoints("weightedBasisDofAtBasisEPoints",numCells,faceCardinality, numBasisEPoints,faceDofDim);
663
664 ScalarViewType faceBasisAtTargetEPoints("faceBasisDofAtTargetEPoints",numCells,faceCardinality, numTargetEPoints,faceDofDim);
665 ScalarViewType wBasisDofAtTargetEPoints("weightedBasisDofAtTargetEPoints",numCells,faceCardinality, numTargetEPoints,faceDofDim);
666
667 ScalarViewType targetDofAtTargetEPoints("targetDofAtTargetEPoints",numCells, numTargetEPoints,faceDofDim);
668 ScalarViewType negPartialProj("negComputedProjection", numCells,numBasisEPoints,faceDofDim);
669
670 ordinal_type offsetBasis = basisEPointsRange(faceDim, iface).first;
671 ordinal_type offsetTarget = targetEPointsRange(faceDim, iface).first;
672 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(faceDim,iface));
673 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(faceDim,iface));
674
675
676 typedef ComputeBasisCoeffsOnFaces_L2<decltype(basisCoeffs), ScalarViewType, decltype(basisEWeights),
677 decltype(computedDofs), decltype(tagToOrdinal), decltype(orts), decltype(targetAtTargetEPoints), decltype(subcellParamFace)> functorTypeFace;
678
679 Kokkos::parallel_for(policy, functorTypeFace(basisCoeffs, negPartialProj,faceBasisDofAtBasisEPoints,
680 basisAtBasisEPoints, basisEWeights, wBasisDofAtBasisEPoints, targetEWeights,
681 basisAtTargetEPoints, wBasisDofAtTargetEPoints, computedDofs, tagToOrdinal,
682 orts, targetAtTargetEPoints,targetDofAtTargetEPoints,
683 subcellParamFace, fieldDim, faceCardinality, offsetBasis,
684 offsetTarget, numVertexDofs+numEdgeDofs, numFaces, faceDim,
685 dim, iface, topoKey, isHCurlBasis, isHDivBasis));
686
687 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
688 ScalarViewType faceMassMat_("faceMassMat_", numCells, faceCardinality, faceCardinality),
689 faceRhsMat_("rhsMat_", numCells, faceCardinality);
690
691 FunctionSpaceTools<DeviceType >::integrate(faceMassMat_, faceBasisDofAtBasisEPoints, wBasisDofAtBasisEPoints);
692 FunctionSpaceTools<DeviceType >::integrate(faceRhsMat_, targetDofAtTargetEPoints, wBasisDofAtTargetEPoints);
693 FunctionSpaceTools<DeviceType >::integrate(faceRhsMat_, negPartialProj, wBasisDofAtBasisEPoints,true);
694
695 ScalarViewType t_("t",numCells, faceCardinality);
696 WorkArrayViewType w_("w",numCells,faceCardinality);
697
698 auto faceDof = Kokkos::subview(tagToOrdinal, faceDim, iface, Kokkos::ALL());
699
700 ElemSystem faceSystem("faceSystem", false);
701 faceSystem.solve(basisCoeffs, faceMassMat_, faceRhsMat_, t_, w_, faceDof, faceCardinality);
702 }
703
704 ordinal_type numElemDofs = cellBasis->getDofCount(dim,0);
705
706
707 if(numElemDofs>0) {
708
709 auto cellDofs = Kokkos::subview(tagToOrdinal, dim, 0, Kokkos::ALL());
710
711 range_type cellPointsRange = targetEPointsRange(dim, 0);
712
713 ordinal_type numTargetEPoints = range_size(targetEPointsRange(dim,0));
714 ordinal_type numBasisEPoints = range_size(basisEPointsRange(dim,0));
715
716 ScalarViewType internalBasisAtBasisEPoints("internalBasisAtBasisEPoints",numCells,numElemDofs, numBasisEPoints, fieldDim);
717 ScalarViewType negPartialProj("negPartialProj", numCells, numBasisEPoints, fieldDim);
718
719 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(dim,0));
720 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(dim,0));
721 ordinal_type offsetBasis = basisEPointsRange(dim, 0).first;
722 ordinal_type offsetTarget = targetEPointsRange(dim, 0).first;
723
724
725 ScalarViewType wBasisAtBasisEPoints("weightedBasisAtBasisEPoints",numCells,numElemDofs, numBasisEPoints,fieldDim);
726 ScalarViewType wBasisDofAtTargetEPoints("weightedBasisAtTargetEPoints",numCells,numElemDofs, numTargetEPoints,fieldDim);
727
728 typedef ComputeBasisCoeffsOnCells_L2<decltype(basisCoeffs), ScalarViewType, decltype(basisEWeights), decltype(computedDofs), decltype(cellDofs)> functorType;
729 Kokkos::parallel_for(policy, functorType( basisCoeffs, negPartialProj, internalBasisAtBasisEPoints,
730 basisAtBasisEPoints, basisEWeights, wBasisAtBasisEPoints, targetEWeights, basisAtTargetEPoints, wBasisDofAtTargetEPoints,
731 computedDofs, cellDofs, fieldDim, numElemDofs, offsetBasis, offsetTarget, numVertexDofs+numEdgeDofs+numFaceDofs));
732
733 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
734 ScalarViewType cellMassMat_("cellMassMat_", numCells, numElemDofs, numElemDofs),
735 cellRhsMat_("rhsMat_", numCells, numElemDofs);
736
737 FunctionSpaceTools<DeviceType >::integrate(cellMassMat_, internalBasisAtBasisEPoints, wBasisAtBasisEPoints);
738 if(fieldDim==1)
739 FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, Kokkos::subview(targetAtTargetEPoints,Kokkos::ALL(),cellPointsRange,Kokkos::ALL()),
740 Kokkos::subview(wBasisDofAtTargetEPoints,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),0));
741 else
742 FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, Kokkos::subview(targetAtTargetEPoints,Kokkos::ALL(),cellPointsRange,Kokkos::ALL()), wBasisDofAtTargetEPoints);
743 FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, negPartialProj, wBasisAtBasisEPoints, true);
744
745 ScalarViewType t_("t",numCells, numElemDofs);
746 WorkArrayViewType w_("w",numCells,numElemDofs);
747 ElemSystem cellSystem("cellSystem", false);
748 cellSystem.solve(basisCoeffs, cellMassMat_, cellRhsMat_, t_, w_, cellDofs, numElemDofs);
749 }
750}
751
752template<typename ViewType1, typename ViewType2>
753struct MultiplyBasisByWeights {
754 const ViewType1 basisAtBasisEPoints_;
755 const ViewType2 basisEWeights_;
756 const ViewType1 wBasisAtBasisEPoints_;
757 const ViewType2 targetEWeights_;
758 const ViewType1 basisAtTargetEPoints_;
759 const ViewType1 wBasisDofAtTargetEPoints_;
760 ordinal_type fieldDim_;
761 ordinal_type numElemDofs_;
762
763 MultiplyBasisByWeights(const ViewType1 basisAtBasisEPoints, const ViewType2 basisEWeights, const ViewType1 wBasisAtBasisEPoints, const ViewType2 targetEWeights,
764 const ViewType1 basisAtTargetEPoints, const ViewType1 wBasisDofAtTargetEPoints,
765 ordinal_type fieldDim, ordinal_type numElemDofs) :
766 basisAtBasisEPoints_(basisAtBasisEPoints), basisEWeights_(basisEWeights), wBasisAtBasisEPoints_(wBasisAtBasisEPoints), targetEWeights_(targetEWeights),
767 basisAtTargetEPoints_(basisAtTargetEPoints), wBasisDofAtTargetEPoints_(wBasisDofAtTargetEPoints),
768 fieldDim_(fieldDim), numElemDofs_(numElemDofs) {}
769
770 void
771 KOKKOS_INLINE_FUNCTION
772 operator()(const ordinal_type ic) const {
773
774 for(ordinal_type j=0; j <numElemDofs_; ++j) {
775 for(ordinal_type d=0; d <fieldDim_; ++d) {
776 for(ordinal_type iq=0; iq <ordinal_type(basisEWeights_.extent(0)); ++iq) {
777 wBasisAtBasisEPoints_(ic,j,iq,d) = basisAtBasisEPoints_(ic,j,iq,d) * basisEWeights_(iq);
778 }
779 for(ordinal_type iq=0; iq <ordinal_type(targetEWeights_.extent(0)); ++iq) {
780 wBasisDofAtTargetEPoints_(ic,j,iq,d) = basisAtTargetEPoints_(ic,j,iq,d)* targetEWeights_(iq);
781 }
782 }
783 }
784 }
785};
786
787template<typename DeviceType>
788template<typename BasisType>
789void
790ProjectionTools<DeviceType>::getL2DGEvaluationPoints(typename BasisType::ScalarViewType ePoints,
791 const BasisType* cellBasis,
793 const EvalPointsType ePointType) {
794
795 ordinal_type dim = cellBasis->getBaseCellTopology().getDimension();
796 auto cellEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getEvalPoints(dim,0,ePointType));
797 RealSpaceTools<DeviceType>::clone(ePoints, cellEPoints);
798}
799
800template<typename DeviceType>
801template<typename basisCoeffsValueType, class ...basisCoeffsProperties,
802typename funValsValueType, class ...funValsProperties,
803typename BasisType,
804typename ortValueType,class ...ortProperties>
805void
806ProjectionTools<DeviceType>::getL2DGBasisCoeffs(Kokkos::DynRankView<basisCoeffsValueType,basisCoeffsProperties...> basisCoeffs,
807 const Kokkos::DynRankView<funValsValueType,funValsProperties...> targetAtTargetEPoints,
808 const Kokkos::DynRankView<ortValueType, ortProperties...> orts,
809 const BasisType* cellBasis,
811
812 typedef typename BasisType::scalarType scalarType;
813 typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
814 const auto cellTopo = cellBasis->getBaseCellTopology();
815
816 ordinal_type dim = cellTopo.getDimension();
817
818 auto basisEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),
819 projStruct->getEvalPoints(dim,0,EvalPointsType::BASIS));
820 auto targetEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),
821 projStruct->getEvalPoints(dim,0,EvalPointsType::TARGET));
822
823
824 ordinal_type numTotalTargetEPoints(targetAtTargetEPoints.extent(1));
825 ordinal_type basisCardinality = cellBasis->getCardinality();
826 ordinal_type numCells = targetAtTargetEPoints.extent(0);
827 const ordinal_type fieldDim = (targetAtTargetEPoints.rank()==2) ? 1 : targetAtTargetEPoints.extent(2);
828
829 ordinal_type numTotalBasisEPoints = projStruct->getNumBasisEvalPoints();
830 ScalarViewType basisAtBasisEPoints("basisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints, fieldDim);
831 ScalarViewType basisAtTargetEPoints("basisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints, fieldDim);
832 {
833 if(fieldDim == 1) {
834 ScalarViewType nonOrientedBasisAtBasisEPoints("nonOrientedBasisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints);
835 ScalarViewType nonOrientedBasisAtTargetEPoints("nonOrientedBasisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints);
836 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL()), targetEPoints);
837 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL()), basisEPoints);
838
839 RealSpaceTools<DeviceType>::clone(nonOrientedBasisAtTargetEPoints, Kokkos::subview(nonOrientedBasisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL()));
840 RealSpaceTools<DeviceType>::clone(nonOrientedBasisAtBasisEPoints, Kokkos::subview(nonOrientedBasisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL()));
841 OrientationTools<DeviceType>::modifyBasisByOrientation(Kokkos::subview(basisAtBasisEPoints, Kokkos::ALL(), Kokkos::ALL(),
842 Kokkos::ALL(),0), nonOrientedBasisAtBasisEPoints, orts, cellBasis);
843 OrientationTools<DeviceType>::modifyBasisByOrientation(Kokkos::subview(basisAtTargetEPoints, Kokkos::ALL(),
844 Kokkos::ALL(), Kokkos::ALL(),0), nonOrientedBasisAtTargetEPoints, orts, cellBasis);
845 }
846 else {
847 ScalarViewType nonOrientedBasisAtBasisEPoints("nonOrientedBasisAtBasisEPoints",numCells,basisCardinality, numTotalBasisEPoints,fieldDim);
848 ScalarViewType nonOrientedBasisAtTargetEPoints("nonOrientedBasisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints,fieldDim);
849 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), targetEPoints);
850 cellBasis->getValues(Kokkos::subview(nonOrientedBasisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), basisEPoints);
851
852 RealSpaceTools<DeviceType>::clone(nonOrientedBasisAtTargetEPoints, Kokkos::subview(nonOrientedBasisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()));
853 RealSpaceTools<DeviceType>::clone(nonOrientedBasisAtBasisEPoints, Kokkos::subview(nonOrientedBasisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()));
854 OrientationTools<DeviceType>::modifyBasisByOrientation(basisAtBasisEPoints, nonOrientedBasisAtBasisEPoints, orts, cellBasis);
855 OrientationTools<DeviceType>::modifyBasisByOrientation(basisAtTargetEPoints, nonOrientedBasisAtTargetEPoints, orts, cellBasis);
856 }
857 }
858
859 const Kokkos::RangePolicy<ExecSpaceType> policy(0, numCells);
860 ordinal_type numElemDofs = cellBasis->getCardinality();
861
862 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(dim,0));
863 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(dim,0));
864
865 ScalarViewType wBasisAtBasisEPoints("weightedBasisAtBasisEPoints",numCells,numElemDofs, numTotalBasisEPoints,fieldDim);
866 ScalarViewType wBasisDofAtTargetEPoints("weightedBasisAtTargetEPoints",numCells,numElemDofs, numTotalTargetEPoints,fieldDim);
867
868 typedef MultiplyBasisByWeights<decltype(basisAtBasisEPoints), decltype(basisEWeights)> functorType;
869 Kokkos::parallel_for( "Multiply basis by weights", policy, functorType(basisAtBasisEPoints, basisEWeights,
870 wBasisAtBasisEPoints, targetEWeights, basisAtTargetEPoints, wBasisDofAtTargetEPoints, fieldDim, numElemDofs));// )){
871
872 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
873 ScalarViewType cellMassMat_("cellMassMat_", numCells, numElemDofs, numElemDofs),
874 cellRhsMat_("rhsMat_", numCells, numElemDofs);
875
876 FunctionSpaceTools<DeviceType >::integrate(cellMassMat_, basisAtBasisEPoints, wBasisAtBasisEPoints);
877 if(fieldDim==1)
878 FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, targetAtTargetEPoints,
879 Kokkos::subview(wBasisDofAtTargetEPoints,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),0));
880 else
881 FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, targetAtTargetEPoints, wBasisDofAtTargetEPoints);
882
883 Kokkos::DynRankView<ordinal_type,Kokkos::HostSpace> hCellDofs("cellDoFs", numElemDofs);
884 for(ordinal_type i=0; i<numElemDofs; ++i) hCellDofs(i) = i;
885 auto cellDofs = Kokkos::create_mirror_view_and_copy(MemSpaceType(),hCellDofs);
886
887 ScalarViewType t_("t",numCells, numElemDofs);
888 WorkArrayViewType w_("w",numCells,numElemDofs);
889 ElemSystem cellSystem("cellSystem", false);
890 cellSystem.solve(basisCoeffs, cellMassMat_, cellRhsMat_, t_, w_, cellDofs, numElemDofs);
891}
892
893template<typename DeviceType>
894template<typename basisViewType, typename targetViewType, typename BasisType>
895void
897 const targetViewType targetAtTargetEPoints,
898 const BasisType* cellBasis,
900
901 typedef typename BasisType::scalarType scalarType;
902 typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
903 const auto cellTopo = cellBasis->getBaseCellTopology();
904
905 ordinal_type dim = cellTopo.getDimension();
906
907 auto basisEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),
908 projStruct->getEvalPoints(dim,0,EvalPointsType::BASIS));
909 auto targetEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),
910 projStruct->getEvalPoints(dim,0,EvalPointsType::TARGET));
911
912 ordinal_type numTotalTargetEPoints(targetAtTargetEPoints.extent(1));
913 ordinal_type basisCardinality = cellBasis->getCardinality();
914 ordinal_type numCells = targetAtTargetEPoints.extent(0);
915 const ordinal_type fieldDim = (targetAtTargetEPoints.rank()==2) ? 1 : targetAtTargetEPoints.extent(2);
916
917 ordinal_type numTotalBasisEPoints = projStruct->getNumBasisEvalPoints();
918 ScalarViewType basisAtBasisEPoints("basisAtBasisEPoints",1,basisCardinality, numTotalBasisEPoints, fieldDim);
919 ScalarViewType basisAtTargetEPoints("basisAtTargetEPoints",numCells,basisCardinality, numTotalTargetEPoints, fieldDim);
920 {
921 if(fieldDim == 1) {
922 cellBasis->getValues(Kokkos::subview(basisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),0), targetEPoints);
923 cellBasis->getValues(Kokkos::subview(basisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL(),0), basisEPoints);
924
925 RealSpaceTools<DeviceType>::clone(basisAtTargetEPoints, Kokkos::subview(basisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()));
926 }
927 else {
928 cellBasis->getValues(Kokkos::subview(basisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), targetEPoints);
929 cellBasis->getValues(Kokkos::subview(basisAtBasisEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()), basisEPoints);
930
931 RealSpaceTools<DeviceType>::clone(basisAtTargetEPoints, Kokkos::subview(basisAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()));
932 }
933 }
934
935 const Kokkos::RangePolicy<ExecSpaceType> policy(0, numCells);
936 ordinal_type numElemDofs = cellBasis->getCardinality();
937
938 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(dim,0));
939 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(dim,0));
940
941 ScalarViewType wBasisAtBasisEPoints("weightedBasisAtBasisEPoints",1,numElemDofs, numTotalBasisEPoints,fieldDim);
942 ScalarViewType wBasisDofAtTargetEPoints("weightedBasisAtTargetEPoints",numCells,numElemDofs, numTotalTargetEPoints,fieldDim);
943 Kokkos::DynRankView<ordinal_type, DeviceType> cellDofs("cellDoFs", numElemDofs);
944
945 Kokkos::parallel_for(Kokkos::RangePolicy<ExecSpaceType>(0,numElemDofs),
946 KOKKOS_LAMBDA (const int &j) {
947 for(ordinal_type d=0; d <fieldDim; ++d) {
948 for(ordinal_type iq=0; iq < numTotalBasisEPoints; ++iq)
949 wBasisAtBasisEPoints(0,j,iq,d) = basisAtBasisEPoints(0,j,iq,d) * basisEWeights(iq);
950 for(ordinal_type iq=0; iq <numTotalTargetEPoints; ++iq) {
951 wBasisDofAtTargetEPoints(0,j,iq,d) = basisAtTargetEPoints(0,j,iq,d)* targetEWeights(iq);
952 }
953 }
954 cellDofs(j) = j;
955 });
956 RealSpaceTools<DeviceType>::clone(wBasisDofAtTargetEPoints, Kokkos::subview(wBasisDofAtTargetEPoints,0,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL()));
957
958 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
959 ScalarViewType cellMassMat_("cellMassMat_", 1, numElemDofs, numElemDofs),
960 cellRhsMat_("rhsMat_", numCells, numElemDofs);
961
962 FunctionSpaceTools<DeviceType >::integrate(cellMassMat_, basisAtBasisEPoints, wBasisAtBasisEPoints);
963 if(fieldDim==1)
964 FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, targetAtTargetEPoints,
965 Kokkos::subview(wBasisDofAtTargetEPoints,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),0));
966 else
967 FunctionSpaceTools<DeviceType >::integrate(cellRhsMat_, targetAtTargetEPoints, wBasisDofAtTargetEPoints);
968
969 ScalarViewType t_("t",1, numElemDofs);
970 WorkArrayViewType w_("w",numCells,numElemDofs);
971 ElemSystem cellSystem("cellSystem", true);
972 cellSystem.solve(basisCoeffs, cellMassMat_, cellRhsMat_, t_, w_, cellDofs, numElemDofs);
973}
974
975}
976}
977
978#endif
979
Header file for Intrepid2::ArrayTools class providing utilities for array operations.
Header file for the abstract base class Intrepid2::DefaultCubatureFactory.
Header file for the Intrepid2::FunctionSpaceTools class.
shards::CellTopology getBaseCellTopology() const
Returns the base cell topology for which the basis is defined. See Shards documentation https://trili...
static void getReferenceEdgeTangent(Kokkos::DynRankView< refEdgeTangentValueType, refEdgeTangentProperties... > refEdgeTangent, const ordinal_type edgeOrd, const shards::CellTopology parentCell)
Computes constant tangent vectors to edges of 2D or 3D reference cells.
static void getReferenceSideNormal(Kokkos::DynRankView< refSideNormalValueType, refSideNormalProperties... > refSideNormal, const ordinal_type sideOrd, const shards::CellTopology parentCell)
Computes constant normal vectors to sides of 2D or 3D reference cells.
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.
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 void getL2BasisCoeffs(Kokkos::DynRankView< basisCoeffsValueType, basisCoeffsProperties... > basisCoeffs, const Kokkos::DynRankView< funValsValueType, funValsProperties... > targetAtEvalPoints, const typename BasisType::ScalarViewType evaluationPoints, 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.
static void getL2DGEvaluationPoints(typename BasisType::ScalarViewType evaluationPoints, const BasisType *cellBasis, ProjectionStruct< DeviceType, typename BasisType::scalarType > *projStruct, const EvalPointsType evalPointType=EvalPointsType::TARGET)
Computes evaluation points for local L2 projection for broken HGRAD HCURL HDIV and HVOL spaces.
static void getL2EvaluationPoints(typename BasisType::ScalarViewType evaluationPoints, const Kokkos::DynRankView< ortValueType, ortProperties... > cellOrientations, const BasisType *cellBasis, ProjectionStruct< DeviceType, typename BasisType::scalarType > *projStruct, const EvalPointsType evalPointType=EvalPointsType::TARGET)
Computes evaluation points for L2 projection.
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.
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...
static KOKKOS_INLINE_FUNCTION void getRefSubcellTangents(TanViewType tangents, const ParamViewType subCellParametrization, const unsigned subcellTopoKey, const ordinal_type subCellOrd, const ordinal_type ort)
Computes the (oriented) subCell tangents.
static void mapToModifiedReference(outPointViewType outPoints, const refPointViewType refPoints, const shards::CellTopology cellTopo, const ordinal_type cellOrt=0)
Computes modified parameterization maps of 1- and 2-subcells with orientation.
static KOKKOS_INLINE_FUNCTION void getRefSideTangentsAndNormal(TanNormViewType tangentsAndNormal, const ParamViewType subCellParametrization, const unsigned subcellTopoKey, const ordinal_type subCellOrd, const ordinal_type ort)
Computes the (oriented) tangents and normal of the side subCell The normals are only defined for side...
static KOKKOS_INLINE_FUNCTION void mapSubcellCoordsToRefCell(coordsViewType cellCoords, const subcellCoordsViewType subCellCoords, const ParamViewType subcellParametrization, const unsigned subcellTopoKey, const ordinal_type subCellOrd, const ordinal_type ort)
Maps points defined on the subCell manifold into the parent Cell.
static void modifyBasisByOrientation(Kokkos::DynRankView< outputValueType, outputProperties... > output, const Kokkos::DynRankView< inputValueType, inputProperties... > input, const OrientationViewType orts, const BasisType *basis)
Modify basis due to orientation.
static void clone(Kokkos::DynRankView< outputValueType, outputProperties... > output, const Kokkos::DynRankView< inputValueType, inputProperties... > input)
Clone input array.
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...
Class to solve a square system A x = b on each cell A is expected to be saddle a point (KKT) matrix o...
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...