Intrepid2
Intrepid2_ProjectionToolsDefHCURL.hpp
Go to the documentation of this file.
1// @HEADER
2// *****************************************************************************
3// Intrepid2 Package
4//
5// Copyright 2007 NTESS and the Intrepid2 contributors.
6// SPDX-License-Identifier: BSD-3-Clause
7// *****************************************************************************
8// @HEADER
9
16#ifndef __INTREPID2_PROJECTIONTOOLSDEFHCURL_HPP__
17#define __INTREPID2_PROJECTIONTOOLSDEFHCURL_HPP__
18
23
24
25namespace Intrepid2 {
26
27namespace FunctorsProjectionTools {
28
29
30template<typename ViewType1, typename ViewType2, typename ViewType3>
32 ViewType1 basisTanAtBasisEPoints_;
33 ViewType1 weightedTanBasisAtBasisEPoints_;
34 const ViewType1 basisAtBasisEPoints_;
35 const ViewType2 basisEWeights_;
36 const ViewType1 refEdgeTangent_;
37 const ViewType3 tagToOrdinal_;
38 ordinal_type edgeDim_;
39 ordinal_type iedge_;
40 ordinal_type offsetBasis_;
41
42 ComputeWBasisEdge_HCurl(ViewType1 basisTanAtBasisEPoints, ViewType1 weightedTanBasisAtBasisEPoints,
43 ViewType1 basisAtBasisEPoints, ViewType2 basisEWeights, ViewType1 refEdgeTangent, ViewType3 tagToOrdinal,
44 ordinal_type edgeDim, ordinal_type iedge, ordinal_type offsetBasis) :
45 basisTanAtBasisEPoints_(basisTanAtBasisEPoints), weightedTanBasisAtBasisEPoints_(weightedTanBasisAtBasisEPoints), basisAtBasisEPoints_(basisAtBasisEPoints),
46 basisEWeights_(basisEWeights), refEdgeTangent_(refEdgeTangent), tagToOrdinal_(tagToOrdinal),
47 edgeDim_(edgeDim), iedge_(iedge), offsetBasis_(offsetBasis) {}
48
49 void
50 KOKKOS_INLINE_FUNCTION
51 operator()(const ordinal_type j, const ordinal_type iq) const {
52 ordinal_type jdof = tagToOrdinal_(edgeDim_, iedge_, j);
53 for(ordinal_type d=0; d <ordinal_type(refEdgeTangent_.extent(0)); ++d)
54 basisTanAtBasisEPoints_(0,j,iq) += refEdgeTangent_(d)*basisAtBasisEPoints_(jdof,offsetBasis_+iq,d);
55 weightedTanBasisAtBasisEPoints_(0,j,iq) = basisTanAtBasisEPoints_(0,j,iq)*basisEWeights_(iq);
56 }
57};
58
59template<typename ViewType1, typename ViewType2, typename ViewType3, typename ViewType4>
61 const ViewType2 targetEWeights_;
62 const ViewType1 basisAtTargetEPoints_;
63 const ViewType1 wTanBasisAtTargetEPoints_;
64 const ViewType3 tagToOrdinal_;
65 const ViewType4 targetAtTargetEPoints_;
66 const ViewType1 targetTanAtTargetEPoints_;
67 const ViewType1 refEdgeTangent_;
68 ordinal_type edgeCardinality_;
69 ordinal_type offsetTarget_;
70 ordinal_type edgeDim_;
71 ordinal_type dim_;
72 ordinal_type iedge_;
73
75 const ViewType2 targetEWeights,
76 const ViewType1 basisAtTargetEPoints, const ViewType1 wTanBasisAtTargetEPoints, const ViewType3 tagToOrdinal,
77 const ViewType4 targetAtTargetEPoints, const ViewType1 targetTanAtTargetEPoints,
78 const ViewType1 refEdgeTangent, ordinal_type edgeCardinality,
79 ordinal_type offsetTarget, ordinal_type edgeDim,
80 ordinal_type dim, ordinal_type iedge) :
81 targetEWeights_(targetEWeights),
82 basisAtTargetEPoints_(basisAtTargetEPoints), wTanBasisAtTargetEPoints_(wTanBasisAtTargetEPoints),
83 tagToOrdinal_(tagToOrdinal), targetAtTargetEPoints_(targetAtTargetEPoints),
84 targetTanAtTargetEPoints_(targetTanAtTargetEPoints),
85 refEdgeTangent_(refEdgeTangent), edgeCardinality_(edgeCardinality),
86 offsetTarget_(offsetTarget), edgeDim_(edgeDim), dim_(dim), iedge_(iedge)
87 {}
88
89 void
90 KOKKOS_INLINE_FUNCTION
91 operator()(const ordinal_type ic) const {
92
93 ordinal_type numTargetEPoints = targetEWeights_.extent(0);
94 typename ViewType1::value_type tmp = 0;
95 for(ordinal_type j=0; j <edgeCardinality_; ++j) {
96 ordinal_type jdof = tagToOrdinal_(edgeDim_, iedge_, j);
97 for(ordinal_type iq=0; iq <numTargetEPoints; ++iq) {
98 tmp = 0;
99 for(ordinal_type d=0; d <dim_; ++d)
100 tmp += refEdgeTangent_(d)*basisAtTargetEPoints_(jdof,offsetTarget_+iq,d);
101 wTanBasisAtTargetEPoints_(ic,j,iq) = tmp*targetEWeights_(iq);
102 }
103 }
104 for(ordinal_type iq=0; iq <numTargetEPoints; ++iq)
105 for(ordinal_type d=0; d <dim_; ++d)
106 targetTanAtTargetEPoints_(ic,iq) += refEdgeTangent_(d)*targetAtTargetEPoints_(ic,offsetTarget_+iq,d);
107 }
108};
109
110
111template<typename ViewType1, typename ViewType2, typename ViewType3, typename ViewType4,
112typename ViewType5>
114 const ViewType1 basisCoeffs_;
115 const ViewType1 negPartialProjTan_;
116 const ViewType1 negPartialProjCurlNormal_;
117 const ViewType1 hgradBasisGradAtBasisEPoints_;
118 const ViewType1 wHgradBasisGradAtBasisEPoints_;
119 const ViewType1 basisCurlAtBasisCurlEPoints_;
120 const ViewType1 basisCurlNormalAtBasisCurlEPoints_;
121 const ViewType1 basisAtBasisEPoints_;
122 const ViewType1 normalTargetCurlAtTargetEPoints_;
123 const ViewType1 basisTanAtBasisEPoints_;
124 const ViewType1 hgradBasisGradAtTargetEPoints_;
125 const ViewType1 wHgradBasisGradAtTargetEPoints_;
126 const ViewType1 wNormalBasisCurlAtBasisCurlEPoints_;
127 const ViewType1 basisCurlAtTargetCurlEPoints_;
128 const ViewType1 wNormalBasisCurlBasisAtTargetCurlEPoints_;
129 const ViewType2 targetAtTargetEPoints_;
130 const ViewType1 targetTanAtTargetEPoints_;
131 const ViewType2 targetCurlAtTargetCurlEPoints_;
132 const ViewType3 basisEWeights_;
133 const ViewType3 targetEWeights_;
134 const ViewType3 basisCurlEWeights_;
135 const ViewType3 targetCurlEWeights_;
136 const ViewType4 tagToOrdinal_;
137 const ViewType4 hGradTagToOrdinal_;
138 const ViewType5 computedDofs_;
139 const ViewType1 refFaceNormal_;
140 ordinal_type offsetBasis_;
141 ordinal_type offsetBasisCurl_;
142 ordinal_type offsetTarget_;
143 ordinal_type offsetTargetCurl_;
144 ordinal_type iface_;
145 ordinal_type hgradCardinality_;
146 ordinal_type numFaces_;
147 ordinal_type numFaceDofs_;
148 ordinal_type numEdgeDofs_;
149 ordinal_type faceDim_;
150 ordinal_type dim_;
151
152
153
154
155 ComputeBasisCoeffsOnFaces_HCurl(const ViewType1 basisCoeffs,
156 const ViewType1 negPartialProjTan, const ViewType1 negPartialProjCurlNormal,
157 const ViewType1 hgradBasisGradAtBasisEPoints, const ViewType1 wHgradBasisGradAtBasisEPoints,
158 const ViewType1 basisCurlAtBasisCurlEPoints, const ViewType1 basisCurlNormalAtBasisCurlEPoints,
159 const ViewType1 basisAtBasisEPoints,
160 const ViewType1 normalTargetCurlAtTargetEPoints,
161 const ViewType1 basisTanAtBasisEPoints,
162 const ViewType1 hgradBasisGradAtTargetEPoints, const ViewType1 wHgradBasisGradAtTargetEPoints,
163 const ViewType1 wNormalBasisCurlAtBasisCurlEPoints, const ViewType1 basisCurlAtTargetCurlEPoints,
164 const ViewType1 wNormalBasisCurlBasisAtTargetCurlEPoints, const ViewType2 targetAtTargetEPoints,
165 const ViewType1 targetTanAtTargetEPoints, const ViewType2 targetCurlAtTargetCurlEPoints,
166 const ViewType3 basisEWeights, const ViewType3 targetEWeights,
167 const ViewType3 basisCurlEWeights, const ViewType3 targetCurlEWeights, const ViewType4 tagToOrdinal,
168 const ViewType4 hGradTagToOrdinal, const ViewType5 computedDofs,
169 const ViewType1 refFaceNormal , ordinal_type offsetBasis,
170 ordinal_type offsetBasisCurl, ordinal_type offsetTarget,
171 ordinal_type offsetTargetCurl, ordinal_type iface,
172 ordinal_type hgradCardinality, ordinal_type numFaces,
173 ordinal_type numFaceDofs, ordinal_type numEdgeDofs,
174 ordinal_type faceDim, ordinal_type dim):
175 basisCoeffs_(basisCoeffs),
176 negPartialProjTan_(negPartialProjTan), negPartialProjCurlNormal_(negPartialProjCurlNormal),
177 hgradBasisGradAtBasisEPoints_(hgradBasisGradAtBasisEPoints), wHgradBasisGradAtBasisEPoints_(wHgradBasisGradAtBasisEPoints),
178 basisCurlAtBasisCurlEPoints_(basisCurlAtBasisCurlEPoints), basisCurlNormalAtBasisCurlEPoints_(basisCurlNormalAtBasisCurlEPoints),
179 basisAtBasisEPoints_(basisAtBasisEPoints),
180 normalTargetCurlAtTargetEPoints_(normalTargetCurlAtTargetEPoints), basisTanAtBasisEPoints_(basisTanAtBasisEPoints),
181 hgradBasisGradAtTargetEPoints_(hgradBasisGradAtTargetEPoints), wHgradBasisGradAtTargetEPoints_(wHgradBasisGradAtTargetEPoints),
182 wNormalBasisCurlAtBasisCurlEPoints_(wNormalBasisCurlAtBasisCurlEPoints), basisCurlAtTargetCurlEPoints_(basisCurlAtTargetCurlEPoints),
183 wNormalBasisCurlBasisAtTargetCurlEPoints_(wNormalBasisCurlBasisAtTargetCurlEPoints), targetAtTargetEPoints_(targetAtTargetEPoints),
184 targetTanAtTargetEPoints_(targetTanAtTargetEPoints), targetCurlAtTargetCurlEPoints_(targetCurlAtTargetCurlEPoints),
185 basisEWeights_(basisEWeights), targetEWeights_(targetEWeights),
186 basisCurlEWeights_(basisCurlEWeights), targetCurlEWeights_(targetCurlEWeights), tagToOrdinal_(tagToOrdinal),
187 hGradTagToOrdinal_(hGradTagToOrdinal), computedDofs_(computedDofs),
188 refFaceNormal_(refFaceNormal), offsetBasis_(offsetBasis),
189 offsetBasisCurl_(offsetBasisCurl), offsetTarget_(offsetTarget),
190 offsetTargetCurl_(offsetTargetCurl), iface_(iface),
191 hgradCardinality_(hgradCardinality), numFaces_(numFaces),
192 numFaceDofs_(numFaceDofs), numEdgeDofs_(numEdgeDofs),
193 faceDim_(faceDim), dim_(dim){}
194
195 void
196 KOKKOS_INLINE_FUNCTION
197 operator()(const ordinal_type ic) const {
198
199 typename ViewType1::value_type n[3] = {refFaceNormal_(0), refFaceNormal_(1), refFaceNormal_(2)};
200 typename ViewType1::value_type tmp=0;
201
202 ordinal_type numBasisEPoints = basisEWeights_.extent(0);
203 ordinal_type numTargetEPoints = targetEWeights_.extent(0);
204 for(ordinal_type j=0; j <hgradCardinality_; ++j) {
205 ordinal_type face_dof = hGradTagToOrdinal_(faceDim_, iface_, j);
206 for(ordinal_type iq=0; iq <numBasisEPoints; ++iq) {
207 for(ordinal_type d=0; d <dim_; ++d) {
208 ordinal_type dp1 = (d+1) % dim_;
209 ordinal_type dp2 = (d+2) % dim_;
210 // basis \times n
211 wHgradBasisGradAtBasisEPoints_(ic,j,iq,d) = (hgradBasisGradAtBasisEPoints_(face_dof,iq,dp1)*n[dp2] - hgradBasisGradAtBasisEPoints_(face_dof,iq,dp2)*n[dp1]) * basisEWeights_(iq);
212 }
213 }
214
215 for(ordinal_type iq=0; iq <numTargetEPoints; ++iq) {
216 for(ordinal_type d=0; d <dim_; ++d) {
217 ordinal_type dp1 = (d+1) % dim_;
218 ordinal_type dp2 = (d+2) % dim_;
219 wHgradBasisGradAtTargetEPoints_(ic,j,iq,d) = (hgradBasisGradAtTargetEPoints_(face_dof,iq,dp1)*n[dp2] - hgradBasisGradAtTargetEPoints_(face_dof,iq,dp2)*n[dp1]) * targetEWeights_(iq);
220 }
221 }
222 }
223
224 ordinal_type numBasisCurlEPoints = basisCurlEWeights_.extent(0);
225 for(ordinal_type j=0; j <numFaceDofs_; ++j) {
226 ordinal_type jdof = tagToOrdinal_(faceDim_, iface_, j);
227 for(ordinal_type iq=0; iq <numBasisEPoints; ++iq) {
228 for(ordinal_type d=0; d <dim_; ++d) {
229 ordinal_type dp1 = (d+1) % dim_;
230 ordinal_type dp2 = (d+2) % dim_;
231 // basis \times n
232 basisTanAtBasisEPoints_(0,j,iq,d) = basisAtBasisEPoints_(jdof,offsetBasis_+iq,dp1)*n[dp2] - basisAtBasisEPoints_(jdof,offsetBasis_+iq,dp2)*n[dp1];
233 }
234 }
235 //Note: we are not considering the jacobian of the orientation map for normals since it is simply a scalar term for the integrals and it does not affect the projection
236 for(ordinal_type iq=0; iq <numBasisCurlEPoints; ++iq) {
237 tmp=0;
238 for(ordinal_type d=0; d <dim_; ++d)
239 tmp += n[d]*basisCurlAtBasisCurlEPoints_(jdof,offsetBasisCurl_+iq,d);
240 basisCurlNormalAtBasisCurlEPoints_(0,j,iq) = tmp;
241 wNormalBasisCurlAtBasisCurlEPoints_(ic,j,iq) = tmp * basisCurlEWeights_(iq);
242 }
243
244 ordinal_type numTargetCurlEPoints = targetCurlEWeights_.extent(0);
245 for(ordinal_type iq=0; iq <numTargetCurlEPoints; ++iq) {
246 tmp=0;
247 // target \cdot n
248 for(ordinal_type d=0; d <dim_; ++d)
249 tmp += n[d]*basisCurlAtTargetCurlEPoints_(jdof,offsetTargetCurl_+iq,d);
250 wNormalBasisCurlBasisAtTargetCurlEPoints_(ic,j,iq) = tmp*targetCurlEWeights_(iq);
251 }
252 }
253
254 for(ordinal_type j=0; j <numEdgeDofs_; ++j) {
255 ordinal_type jdof = computedDofs_(j);
256 for(ordinal_type iq=0; iq <numBasisEPoints; ++iq)
257 for(ordinal_type d=0; d <dim_; ++d) {
258 ordinal_type dp1 = (d+1) % dim_;
259 ordinal_type dp2 = (d+2) % dim_;
260 negPartialProjCurlNormal_(ic,iq) -= n[d]*basisCoeffs_(ic,jdof)*basisCurlAtBasisCurlEPoints_(jdof,offsetBasisCurl_+iq,d);
261 // basis \times n
262 negPartialProjTan_(ic,iq,d) -= (basisAtBasisEPoints_(jdof,offsetBasis_+iq,dp1)*n[dp2] - basisAtBasisEPoints_(jdof,offsetBasis_+iq,dp2)*n[dp1])*basisCoeffs_(ic,jdof);
263 }
264 }
265
266 ordinal_type numTargetCurlEPoints = targetCurlEWeights_.extent(0);
267 for(ordinal_type iq=0; iq <numTargetEPoints; ++iq)
268 for(ordinal_type d=0; d <dim_; ++d) {
269 ordinal_type dp1 = (d+1) % dim_;
270 ordinal_type dp2 = (d+2) % dim_;
271 // target \times n
272 targetTanAtTargetEPoints_(ic,iq,d) = (targetAtTargetEPoints_(ic,offsetTarget_+iq,dp1)*n[dp2] - targetAtTargetEPoints_(ic,offsetTarget_+iq,dp2)*n[dp1]);
273 }
274
275 for(ordinal_type iq=0; iq <numTargetCurlEPoints; ++iq)
276 for(ordinal_type d=0; d <dim_; ++d) {
277 // target \cdot n
278 normalTargetCurlAtTargetEPoints_(ic,iq) += n[d]*targetCurlAtTargetCurlEPoints_(ic,offsetTargetCurl_+iq,d);
279 }
280 }
281};
282
283
284template<typename ViewType1, typename ViewType2, typename ViewType3,
285typename ViewType4>
287 const ViewType1 basisCoeffs_;
288 const ViewType1 negPartialProj_;
289 const ViewType1 negPartialProjCurl_;
290 const ViewType1 cellBasisAtBasisEPoints_;
291 const ViewType1 cellBasisCurlAtBasisCurlEPoints_;
292 const ViewType1 basisAtBasisEPoints_;
293 const ViewType1 hgradBasisGradAtBasisEPoints_;
294 const ViewType1 basisCurlAtBasisCurlEPoints_;
295 const ViewType1 hgradBasisGradAtTargetEPoints_;
296 const ViewType1 basisCurlAtTargetCurlEPoints_;
297 const ViewType2 basisEWeights_;
298 const ViewType2 basisCurlEWeights_;
299 const ViewType1 wHgradBasisGradAtBasisEPoints_;
300 const ViewType1 wBasisCurlAtBasisCurlEPoints_;
301 const ViewType2 targetEWeights_;
302 const ViewType2 targetCurlEWeights_;
303 const ViewType1 wHgradBasisGradAtTargetEPoints_;
304 const ViewType1 wBasisCurlAtTargetCurlEPoints_;
305 const ViewType3 computedDofs_;
306 const ViewType4 tagToOrdinal_;
307 const ViewType4 hGradTagToOrdinal_;
308 ordinal_type numCellDofs_;
309 ordinal_type hgradCardinality_;
310 ordinal_type offsetBasis_;
311 ordinal_type offsetBasisCurl_;
312 ordinal_type offsetTargetCurl_;
313 ordinal_type numEdgeFaceDofs_;
314 ordinal_type dim_;
315 ordinal_type derDim_;
316
317 ComputeBasisCoeffsOnCell_HCurl(const ViewType1 basisCoeffs, ViewType1 negPartialProj, ViewType1 negPartialProjCurl,
318 const ViewType1 cellBasisAtBasisEPoints, const ViewType1 cellBasisCurlAtBasisCurlEPoints,
319 const ViewType1 basisAtBasisEPoints, const ViewType1 hgradBasisGradAtBasisEPoints, const ViewType1 basisCurlAtBasisCurlEPoints,
320 const ViewType1 hgradBasisGradAtTargetEPoints, const ViewType1 basisCurlAtTargetCurlEPoints,
321 const ViewType2 basisEWeights, const ViewType2 basisCurlEWeights,
322 const ViewType1 wHgradBasisGradAtBasisEPoints, const ViewType1 wBasisCurlAtBasisCurlEPoints,
323 const ViewType2 targetEWeights, const ViewType2 targetCurlEWeights,
324 const ViewType1 wHgradBasisGradAtTargetEPoints,
325 const ViewType1 wBasisCurlAtTargetCurlEPoints, const ViewType3 computedDofs,
326 const ViewType4 tagToOrdinal, const ViewType4 hGradTagToOrdinal,
327 ordinal_type numCellDofs, ordinal_type hgradCardinality,
328 ordinal_type offsetBasis, ordinal_type offsetBasisCurl, ordinal_type offsetTargetCurl,
329 ordinal_type numEdgeFaceDofs, ordinal_type dim, ordinal_type derDim) :
330 basisCoeffs_(basisCoeffs), negPartialProj_(negPartialProj), negPartialProjCurl_(negPartialProjCurl),
331 cellBasisAtBasisEPoints_(cellBasisAtBasisEPoints), cellBasisCurlAtBasisCurlEPoints_(cellBasisCurlAtBasisCurlEPoints),
332 basisAtBasisEPoints_(basisAtBasisEPoints), hgradBasisGradAtBasisEPoints_(hgradBasisGradAtBasisEPoints),
333 basisCurlAtBasisCurlEPoints_(basisCurlAtBasisCurlEPoints),
334 hgradBasisGradAtTargetEPoints_(hgradBasisGradAtTargetEPoints),
335 basisCurlAtTargetCurlEPoints_(basisCurlAtTargetCurlEPoints),
336 basisEWeights_(basisEWeights), basisCurlEWeights_(basisCurlEWeights),
337 wHgradBasisGradAtBasisEPoints_(wHgradBasisGradAtBasisEPoints),
338 wBasisCurlAtBasisCurlEPoints_(wBasisCurlAtBasisCurlEPoints),
339 targetEWeights_(targetEWeights), targetCurlEWeights_(targetCurlEWeights),
340 wHgradBasisGradAtTargetEPoints_(wHgradBasisGradAtTargetEPoints),
341 wBasisCurlAtTargetCurlEPoints_(wBasisCurlAtTargetCurlEPoints),
342 computedDofs_(computedDofs), tagToOrdinal_(tagToOrdinal), hGradTagToOrdinal_(hGradTagToOrdinal),
343 numCellDofs_(numCellDofs), hgradCardinality_(hgradCardinality),
344 offsetBasis_(offsetBasis), offsetBasisCurl_(offsetBasisCurl), offsetTargetCurl_(offsetTargetCurl),
345 numEdgeFaceDofs_(numEdgeFaceDofs), dim_(dim), derDim_(derDim) {}
346
347 void
348 KOKKOS_INLINE_FUNCTION
349 operator()(const ordinal_type ic) const {
350
351 ordinal_type numBasisPoints = basisEWeights_.extent(0);
352 ordinal_type numBasisCurlPoints = basisCurlEWeights_.extent(0);
353 ordinal_type numTargetPoints = targetEWeights_.extent(0);
354 ordinal_type numTargetCurlPoints = targetCurlEWeights_.extent(0);
355 for(ordinal_type j=0; j <hgradCardinality_; ++j) {
356 ordinal_type idof = hGradTagToOrdinal_(dim_, 0, j);
357 for(ordinal_type d=0; d <dim_; ++d) {
358 for(ordinal_type iq=0; iq <numBasisPoints; ++iq)
359 wHgradBasisGradAtBasisEPoints_(ic,j,iq,d) = hgradBasisGradAtBasisEPoints_(idof,iq,d)*basisEWeights_(iq);
360 for(ordinal_type iq=0; iq <numTargetPoints; ++iq)
361 wHgradBasisGradAtTargetEPoints_(ic,j,iq,d) = hgradBasisGradAtTargetEPoints_(idof,iq,d)*targetEWeights_(iq);
362 }
363 }
364 for(ordinal_type j=0; j <numCellDofs_; ++j) {
365 ordinal_type idof = tagToOrdinal_(dim_, 0, j);
366 for(ordinal_type d=0; d <dim_; ++d)
367 for(ordinal_type iq=0; iq <numBasisPoints; ++iq)
368 cellBasisAtBasisEPoints_(0,j,iq,d)=basisAtBasisEPoints_(idof,offsetBasis_+iq,d);
369
370 for(ordinal_type d=0; d <derDim_; ++d) {
371 for(ordinal_type iq=0; iq <numBasisCurlPoints; ++iq) {
372 cellBasisCurlAtBasisCurlEPoints_(0,j,iq,d)=basisCurlAtBasisCurlEPoints_.access(idof,offsetBasisCurl_+iq,d);
373 wBasisCurlAtBasisCurlEPoints_(ic,j,iq,d)=cellBasisCurlAtBasisCurlEPoints_(0,j,iq,d)*basisCurlEWeights_(iq);
374 }
375 for(ordinal_type iq=0; iq <numTargetCurlPoints; ++iq)
376 wBasisCurlAtTargetCurlEPoints_(ic,j,iq,d) = basisCurlAtTargetCurlEPoints_.access(idof,offsetTargetCurl_+iq,d)*targetCurlEWeights_(iq);
377 }
378 }
379 for(ordinal_type j=0; j < numEdgeFaceDofs_; ++j) {
380 ordinal_type jdof = computedDofs_(j);
381 for(ordinal_type d=0; d <derDim_; ++d)
382 for(ordinal_type iq=0; iq <numBasisCurlPoints; ++iq)
383 negPartialProjCurl_(ic,iq,d) -= basisCoeffs_(ic,jdof)*basisCurlAtBasisCurlEPoints_.access(jdof,offsetBasisCurl_+iq,d);
384 for(ordinal_type d=0; d <dim_; ++d)
385 for(ordinal_type iq=0; iq <numBasisPoints; ++iq)
386 negPartialProj_(ic,iq,d) -= basisCoeffs_(ic,jdof)*basisAtBasisEPoints_(jdof,offsetBasis_+iq,d);
387 }
388 }
389};
390
391} // FunctorsProjectionTools namespace
392
393
394template<typename DeviceType>
395template<typename basisCoeffsValueType, class ...basisCoeffsProperties,
396typename funValsValueType, class ...funValsProperties,
397typename BasisType,
398typename ortValueType,class ...ortProperties>
399void
400ProjectionTools<DeviceType>::getHCurlBasisCoeffs(Kokkos::DynRankView<basisCoeffsValueType,basisCoeffsProperties...> basisCoeffs,
401 const Kokkos::DynRankView<funValsValueType,funValsProperties...> targetAtTargetEPoints,
402 const Kokkos::DynRankView<funValsValueType,funValsProperties...> targetCurlAtTargetCurlEPoints,
403 const Kokkos::DynRankView<ortValueType, ortProperties...> orts,
404 const BasisType* cellBasis,
406
407 typedef typename BasisType::scalarType scalarType;
408 typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
409 typedef Kokkos::pair<ordinal_type,ordinal_type> range_type;
410 const auto cellTopo = cellBasis->getBaseCellTopology();
411 ordinal_type dim = cellTopo.getDimension();
412 ordinal_type basisCardinality = cellBasis->getCardinality();
413 ordinal_type numCells = targetAtTargetEPoints.extent(0);
414 const ordinal_type edgeDim = 1;
415 const ordinal_type faceDim = 2;
416 const ordinal_type derDim = dim == 3 ? dim : 1;
417
418 const Kokkos::RangePolicy<ExecSpaceType> policy(0, numCells);
419
420 const std::string& name = cellBasis->getName();
421
422 ordinal_type numEdges = (cellBasis->getDofCount(1, 0) > 0) ? cellTopo.getEdgeCount() : 0;
423 ordinal_type numFaces = (cellBasis->getDofCount(2, 0) > 0) ? cellTopo.getFaceCount() : 0;
424
425 ScalarViewType refEdgeTangent("refEdgeTangent", dim);
426 ScalarViewType refFaceNormal("refFaceNormal", dim);
427
428 ordinal_type numEdgeDofs(0);
429 for(ordinal_type ie=0; ie<numEdges; ++ie)
430 numEdgeDofs += cellBasis->getDofCount(edgeDim,ie);
431
432 ordinal_type numTotalFaceDofs(0);
433 for(ordinal_type iface=0; iface<numFaces; ++iface)
434 numTotalFaceDofs += cellBasis->getDofCount(faceDim,iface);
435
436 auto tagToOrdinal = Kokkos::create_mirror_view_and_copy(MemSpaceType(), cellBasis->getAllDofOrdinal());
437
438 Kokkos::View<ordinal_type*, DeviceType> computedDofs("computedDofs",numEdgeDofs+numTotalFaceDofs);
439
440 auto targetEPointsRange = projStruct->getTargetPointsRange();
441 auto targetCurlEPointsRange = projStruct->getTargetDerivPointsRange();
442
443 auto basisEPointsRange = projStruct->getBasisPointsRange();
444 auto basisCurlEPointsRange = projStruct->getBasisDerivPointsRange();
445
446 ordinal_type numTotalBasisEPoints = projStruct->getNumBasisEvalPoints(), numTotalBasisCurlEPoints = projStruct->getNumBasisDerivEvalPoints();
447 auto basisEPoints = projStruct->getAllEvalPoints(EvalPointsType::BASIS);
448 auto basisCurlEPoints = projStruct->getAllDerivEvalPoints(EvalPointsType::BASIS);
449
450 ordinal_type numTotalTargetEPoints = projStruct->getNumTargetEvalPoints(), numTotalTargetCurlEPoints = projStruct->getNumTargetDerivEvalPoints();
451 auto targetEPoints = projStruct->getAllEvalPoints(EvalPointsType::TARGET);
452 auto targetCurlEPoints = projStruct->getAllDerivEvalPoints(EvalPointsType::TARGET);
453
454 ScalarViewType basisAtBasisEPoints("basisAtBasisEPoints",basisCardinality, numTotalBasisEPoints, dim);
455 ScalarViewType basisAtTargetEPoints("basisAtTargetEPoints",basisCardinality, numTotalTargetEPoints, dim);
456 cellBasis->getValues(basisAtTargetEPoints, targetEPoints);
457 cellBasis->getValues(basisAtBasisEPoints, basisEPoints);
458
459 ScalarViewType basisCurlAtBasisCurlEPoints;
460 ScalarViewType basisCurlAtTargetCurlEPoints;
461 if(numTotalBasisCurlEPoints>0) {
462 ScalarViewType nonOrientedBasisCurlAtTargetCurlEPoints, nonOrientedBasisCurlAtBasisCurlEPoints;
463 if (dim == 3) {
464 basisCurlAtBasisCurlEPoints = ScalarViewType ("basisCurlAtBasisCurlEPoints",basisCardinality, numTotalBasisCurlEPoints, dim);
465 basisCurlAtTargetCurlEPoints = ScalarViewType("basisCurlAtTargetCurlEPoints",basisCardinality, numTotalTargetCurlEPoints, dim);
466 } else {
467 basisCurlAtBasisCurlEPoints = ScalarViewType ("basisCurlAtBasisCurlEPoints",basisCardinality, numTotalBasisCurlEPoints);
468 basisCurlAtTargetCurlEPoints = ScalarViewType("basisCurlAtTargetCurlEPoints",basisCardinality, numTotalTargetCurlEPoints);
469 }
470
471 cellBasis->getValues(basisCurlAtBasisCurlEPoints, basisCurlEPoints,OPERATOR_CURL);
472 cellBasis->getValues(basisCurlAtTargetCurlEPoints, targetCurlEPoints,OPERATOR_CURL);
473 }
474
475 ScalarViewType refBasisCoeffs("refBasisCoeffs", basisCoeffs.extent(0), basisCoeffs.extent(1));
476
477 ordinal_type computedDofsCount = 0;
478 for(ordinal_type ie=0; ie<numEdges; ++ie) {
479
480 ordinal_type edgeCardinality = cellBasis->getDofCount(edgeDim,ie);
481 ordinal_type numBasisEPoints = range_size(basisEPointsRange(edgeDim, ie));
482 ordinal_type numTargetEPoints = range_size(targetEPointsRange(edgeDim, ie));
483
484 CellTools<DeviceType>::getReferenceEdgeTangent(refEdgeTangent, ie, cellTopo);
485
486 ScalarViewType basisTanAtBasisEPoints("basisTanAtBasisEPoints",1,edgeCardinality, numBasisEPoints);
487 ScalarViewType weightedTanBasisAtBasisEPoints("weightedTanBasisAtBasisEPoints",1,edgeCardinality, numBasisEPoints);
488 ScalarViewType weightedTanBasisAtTargetEPoints("weightedTanBasisAtTargetEPoints",numCells,edgeCardinality, numTargetEPoints);
489 ScalarViewType targetTanAtTargetEPoints("normalTargetAtTargetEPoints",numCells, numTargetEPoints);
490
491 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(edgeDim,ie));
492 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(edgeDim,ie));
493
494 ordinal_type offsetBasis = basisEPointsRange(edgeDim, ie).first;
495 ordinal_type offsetTarget = targetEPointsRange(edgeDim, ie).first;
496
497 using functorTypeWBasisEdge = FunctorsProjectionTools::ComputeWBasisEdge_HCurl<ScalarViewType, decltype(basisEWeights), decltype(tagToOrdinal)>;
498 Kokkos::parallel_for(Kokkos::MDRangePolicy<ExecSpaceType, Kokkos::Rank<2> >({0,0}, {edgeCardinality,numBasisEPoints}),
499 functorTypeWBasisEdge(basisTanAtBasisEPoints,weightedTanBasisAtBasisEPoints,basisAtBasisEPoints,basisEWeights,refEdgeTangent,tagToOrdinal,edgeDim,ie,offsetBasis));
500
501 using functorTypeEdge = FunctorsProjectionTools::ComputeBasisCoeffsOnEdges_HCurl<ScalarViewType, decltype(targetEWeights), decltype(tagToOrdinal), decltype(targetAtTargetEPoints)>;
502 Kokkos::parallel_for(policy, functorTypeEdge(targetEWeights,
503 basisAtTargetEPoints, weightedTanBasisAtTargetEPoints, tagToOrdinal,
504 targetAtTargetEPoints, targetTanAtTargetEPoints,
505 refEdgeTangent, edgeCardinality,
506 offsetTarget, edgeDim,
507 dim, ie));
508
509 ScalarViewType edgeMassMat_("edgeMassMat_", 1, edgeCardinality+1, edgeCardinality+1),
510 edgeRhsMat_("rhsMat_", numCells, edgeCardinality+1);
511
512 ScalarViewType eWeights_("eWeights_", 1, 1, basisEWeights.extent(0)), targetEWeights_("targetEWeights", numCells, 1, targetEWeights.extent(0));
513 RealSpaceTools<DeviceType>::clone(eWeights_, basisEWeights);
514 RealSpaceTools<DeviceType>::clone(targetEWeights_, targetEWeights);
515
516 range_type range_H(0, edgeCardinality);
517 range_type range_B(edgeCardinality, edgeCardinality+1);
518 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(edgeMassMat_,Kokkos::ALL(),range_H,range_H), basisTanAtBasisEPoints, weightedTanBasisAtBasisEPoints);
519 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(edgeMassMat_,Kokkos::ALL(),range_H,range_B), basisTanAtBasisEPoints, eWeights_);
520 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(edgeRhsMat_,Kokkos::ALL(),range_H), targetTanAtTargetEPoints, weightedTanBasisAtTargetEPoints);
521 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(edgeRhsMat_,Kokkos::ALL(),range_B), targetTanAtTargetEPoints, targetEWeights_);
522
523 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
524 ScalarViewType t_("t",numCells, edgeCardinality+1);
525 WorkArrayViewType w_("w",numCells, edgeCardinality+1);
526
527 auto edgeDofs = Kokkos::subview(tagToOrdinal, edgeDim, ie, range_type(0,edgeCardinality));
528 ElemSystem edgeSystem("edgeSystem", true);
529 edgeSystem.solve(refBasisCoeffs, edgeMassMat_, edgeRhsMat_, t_, w_, edgeDofs, edgeCardinality, 1);
530
531 auto computedEdgeDofs = Kokkos::subview(computedDofs, range_type(computedDofsCount,computedDofsCount+edgeCardinality));
532 deep_copy(computedEdgeDofs, edgeDofs);
533 computedDofsCount += edgeCardinality;
534
535 }
536
537 typename RefSubcellParametrization<DeviceType>::ConstViewType subcellParamFace;
539 if(numFaces>0) {
540 subcellParamFace = RefSubcellParametrization<DeviceType>::get(faceDim, cellBasis->getBaseCellTopology().getKey());
541 if(cellTopo.getKey() == shards::getCellTopologyData<shards::Hexahedron<8> >()->key)
542 hgradBasis = new Basis_HGRAD_HEX_Cn_FEM<DeviceType,scalarType,scalarType>(cellBasis->getDegree(),POINTTYPE_WARPBLEND);
543 else if(cellTopo.getKey() == shards::getCellTopologyData<shards::Tetrahedron<4> >()->key)
544 hgradBasis = new Basis_HGRAD_TET_Cn_FEM<DeviceType,scalarType,scalarType>(cellBasis->getDegree(),POINTTYPE_WARPBLEND);
545 else if(cellTopo.getKey() == shards::getCellTopologyData<shards::Wedge<6> >()->key) {
546 hgradBasis = new typename DerivedNodalBasisFamily<DeviceType,scalarType,scalarType>::HGRAD_WEDGE(cellBasis->getDegree(),POINTTYPE_WARPBLEND);
547 }
548 else {
549 std::stringstream ss;
550 ss << ">>> ERROR (Intrepid2::ProjectionTools::getHCurlBasisCoeffs): "
551 << "Method not implemented for basis " << name;
552 INTREPID2_TEST_FOR_EXCEPTION( true, std::runtime_error, ss.str().c_str() );
553 }
554 }
555 for(ordinal_type iface=0; iface<numFaces; ++iface) {
556
557 ordinal_type numTargetEPoints = range_size(targetEPointsRange(faceDim, iface));
558 ordinal_type numTargetCurlEPoints = range_size(targetCurlEPointsRange(faceDim, iface));
559 ordinal_type numBasisEPoints = range_size(basisEPointsRange(faceDim, iface));
560 ordinal_type numBasisCurlEPoints = range_size(basisCurlEPointsRange(faceDim, iface));
561
562 ordinal_type numFaceDofs = cellBasis->getDofCount(faceDim,iface);
563
564 ScalarViewType hgradBasisGradAtBasisEPoints("hgradBasisGradAtBasisEPoints", hgradBasis->getCardinality(), numBasisEPoints, dim);
565 ScalarViewType hgradBasisGradAtTargetEPoints("hgradBasisGradAtTargetEPoints", hgradBasis->getCardinality(), numTargetEPoints, dim);
566
567 ordinal_type hgradCardinality = hgradBasis->getDofCount(faceDim,iface);
568
569 hgradBasis->getValues(hgradBasisGradAtBasisEPoints, Kokkos::subview(basisEPoints, basisEPointsRange(faceDim, iface), Kokkos::ALL()), OPERATOR_GRAD);
570 hgradBasis->getValues(hgradBasisGradAtTargetEPoints, Kokkos::subview(targetEPoints, targetEPointsRange(faceDim, iface), Kokkos::ALL()), OPERATOR_GRAD);
571
572 //no need to orient these basis as they act locally as test functions.
573
574 auto hGradTagToOrdinal = Kokkos::create_mirror_view_and_copy(MemSpaceType(), hgradBasis->getAllDofOrdinal());
575
576 ScalarViewType basisTanAtBasisEPoints("basisTanAtBasisEPoints",1,numFaceDofs, numBasisEPoints,dim);
577 ScalarViewType basisCurlNormalAtBasisCurlEPoints("normaBasisCurlAtBasisEPoints",1,numFaceDofs, numBasisCurlEPoints);
578 ScalarViewType wNormalBasisCurlAtBasisCurlEPoints("weightedNormalBasisCurlAtBasisEPoints",numCells,numFaceDofs, numBasisCurlEPoints);
579
580 ScalarViewType targetTanAtTargetEPoints("targetTanAtTargetEPoints",numCells, numTargetEPoints, dim);
581 ScalarViewType normalTargetCurlAtTargetEPoints("normalTargetCurlAtTargetEPoints",numCells, numTargetCurlEPoints);
582 ScalarViewType wNormalBasisCurlBasisAtTargetCurlEPoints("weightedNormalBasisCurlAtTargetCurlEPoints",numCells,numFaceDofs, numTargetCurlEPoints);
583
584 ScalarViewType wHgradBasisGradAtBasisEPoints("wHgradBasisGradAtBasisEPoints",numCells, hgradCardinality, numBasisEPoints, dim);
585 ScalarViewType wHgradBasisGradAtTargetEPoints("wHgradBasisGradAtTargetEPoints",numCells, hgradCardinality, numTargetEPoints, dim);
586
587 ScalarViewType negPartialProjCurlNormal("mNormalComputedProjection", numCells,numBasisEPoints);
588 ScalarViewType negPartialProjTan("negPartialProjTan", numCells,numBasisEPoints,dim);
589
590 ordinal_type offsetBasis = basisEPointsRange(faceDim, iface).first;
591 ordinal_type offsetBasisCurl = basisCurlEPointsRange(faceDim, iface).first;
592 ordinal_type offsetTarget = targetEPointsRange(faceDim, iface).first;
593 ordinal_type offsetTargetCurl = targetCurlEPointsRange(faceDim, iface).first;
594
595 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(faceDim,iface));
596 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(faceDim,iface));
597 auto targetCurlEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetDerivEvalWeights(faceDim,iface));
598 auto basisCurlEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisDerivEvalWeights(faceDim,iface));
599
600 CellTools<DeviceType>::getReferenceFaceNormal(refFaceNormal, iface, cellTopo);
601
602 using functorTypeFaces = FunctorsProjectionTools::ComputeBasisCoeffsOnFaces_HCurl<ScalarViewType, decltype(targetAtTargetEPoints), decltype(basisEWeights),
603 decltype(tagToOrdinal), decltype(computedDofs)>;
604 Kokkos::parallel_for(policy, functorTypeFaces(refBasisCoeffs,
605 negPartialProjTan, negPartialProjCurlNormal,
606 hgradBasisGradAtBasisEPoints, wHgradBasisGradAtBasisEPoints,
607 basisCurlAtBasisCurlEPoints, basisCurlNormalAtBasisCurlEPoints,
608 basisAtBasisEPoints,
609 normalTargetCurlAtTargetEPoints, basisTanAtBasisEPoints,
610 hgradBasisGradAtTargetEPoints, wHgradBasisGradAtTargetEPoints,
611 wNormalBasisCurlAtBasisCurlEPoints, basisCurlAtTargetCurlEPoints,
612 wNormalBasisCurlBasisAtTargetCurlEPoints, targetAtTargetEPoints,
613 targetTanAtTargetEPoints, targetCurlAtTargetCurlEPoints,
614 basisEWeights, targetEWeights,
615 basisCurlEWeights, targetCurlEWeights, tagToOrdinal,
616 hGradTagToOrdinal,
617 computedDofs, refFaceNormal, offsetBasis,
618 offsetBasisCurl, offsetTarget,
619 offsetTargetCurl, iface,
620 hgradCardinality, numFaces,
621 numFaceDofs, numEdgeDofs,
622 faceDim, dim));
623
624
625 ScalarViewType faceMassMat_("faceMassMat_", 1, numFaceDofs+hgradCardinality, numFaceDofs+hgradCardinality),
626 faceRhsMat_("rhsMat_", numCells, numFaceDofs+hgradCardinality);
627 range_type range_H(0, numFaceDofs);
628 range_type range_B(numFaceDofs, numFaceDofs+hgradCardinality);
629 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(faceMassMat_,Kokkos::ALL(),range_H,range_H), basisCurlNormalAtBasisCurlEPoints, Kokkos::subview(wNormalBasisCurlAtBasisCurlEPoints, std::make_pair(0,1), Kokkos::ALL(), Kokkos::ALL()) );
630 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(faceMassMat_,Kokkos::ALL(),range_H,range_B), basisTanAtBasisEPoints, Kokkos::subview(wHgradBasisGradAtBasisEPoints, std::make_pair(0,1), Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL()) );
631
632 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(faceRhsMat_,Kokkos::ALL(),range_H), normalTargetCurlAtTargetEPoints, wNormalBasisCurlBasisAtTargetCurlEPoints);
633 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(faceRhsMat_,Kokkos::ALL(),range_H), negPartialProjCurlNormal, wNormalBasisCurlAtBasisCurlEPoints,true);
634
635 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(faceRhsMat_,Kokkos::ALL(),range_B), targetTanAtTargetEPoints, wHgradBasisGradAtTargetEPoints);
636 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(faceRhsMat_,Kokkos::ALL(),range_B), negPartialProjTan, wHgradBasisGradAtBasisEPoints,true);
637
638
639 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
640 ScalarViewType t_("t",numCells, numFaceDofs+hgradCardinality);
641 WorkArrayViewType w_("w",numCells, numFaceDofs+hgradCardinality);
642
643 auto faceDofs = Kokkos::subview(tagToOrdinal, faceDim, iface, range_type(0,numFaceDofs));
644 ElemSystem faceSystem( "faceSystem", true);
645 faceSystem.solve(refBasisCoeffs, faceMassMat_, faceRhsMat_, t_, w_, faceDofs, numFaceDofs, hgradCardinality);
646
647 auto computedFaceDofs = Kokkos::subview(computedDofs, range_type(computedDofsCount,computedDofsCount+numFaceDofs));
648 deep_copy(computedFaceDofs, faceDofs);
649 computedDofsCount += numFaceDofs;
650
651 }
652 delete hgradBasis;
653
654 ordinal_type numCellDofs = cellBasis->getDofCount(dim,0);
655 if(numCellDofs>0) {
656 if(cellTopo.getKey() == shards::getCellTopologyData<shards::Hexahedron<8> >()->key)
657 hgradBasis = new Basis_HGRAD_HEX_Cn_FEM<DeviceType,scalarType,scalarType>(cellBasis->getDegree());
658 else if(cellTopo.getKey() == shards::getCellTopologyData<shards::Tetrahedron<4> >()->key)
659 hgradBasis = new Basis_HGRAD_TET_Cn_FEM<DeviceType,scalarType,scalarType>(cellBasis->getDegree(),POINTTYPE_WARPBLEND);
660 else if(cellTopo.getKey() == shards::getCellTopologyData<shards::Wedge<6> >()->key)
661 hgradBasis = new typename DerivedNodalBasisFamily<DeviceType,scalarType,scalarType>::HGRAD_WEDGE(cellBasis->getDegree(),POINTTYPE_WARPBLEND);
662 else if(cellTopo.getKey() == shards::getCellTopologyData<shards::Triangle<3> >()->key)
663 hgradBasis = new Basis_HGRAD_TRI_Cn_FEM<DeviceType,scalarType,scalarType>(cellBasis->getDegree(),POINTTYPE_WARPBLEND);
664 else if(cellTopo.getKey() == shards::getCellTopologyData<shards::Quadrilateral<4> >()->key)
665 hgradBasis = new Basis_HGRAD_QUAD_Cn_FEM<DeviceType,scalarType,scalarType>(cellBasis->getDegree(),POINTTYPE_WARPBLEND);
666 else {
667 std::stringstream ss;
668 ss << ">>> ERROR (Intrepid2::ProjectionTools::getHCurlBasisCoeffs): "
669 << "Method not implemented for basis " << name;
670 INTREPID2_TEST_FOR_EXCEPTION( true, std::runtime_error, ss.str().c_str() );
671 }
672
673 range_type cellPointsRange = targetEPointsRange(dim, 0);
674 range_type cellCurlPointsRange = targetCurlEPointsRange(dim, 0);
675
676 ordinal_type numTargetCurlEPoints = range_size(targetCurlEPointsRange(dim,0));
677 ordinal_type numBasisCurlEPoints = range_size(basisCurlEPointsRange(dim,0));
678 ordinal_type numBasisEPoints = range_size(basisEPointsRange(dim,0));
679 ordinal_type numTargetEPoints = range_size(targetEPointsRange(dim,0));
680
681 ScalarViewType hgradBasisGradAtBasisEPoints("hgradBasisGradAtBasisEPoints",hgradBasis->getCardinality(), numBasisEPoints, dim);
682 ScalarViewType hgradBasisGradAtTargetEPoints("hgradBasisGradAtTargetEPoints",hgradBasis->getCardinality(), numTargetEPoints, dim);
683
684 ordinal_type hgradCardinality = hgradBasis->getDofCount(dim,0);
685 ScalarViewType wHgradBasisGradAtTargetEPoints("wHgradBasisGradAtTargetEPoints",numCells, hgradCardinality, numTargetEPoints, dim);
686 ScalarViewType wHgradBasisGradAtBasisEPoints("wHgradBasisGradAtBasisEPoints",numCells, hgradCardinality, numBasisEPoints, dim);
687
688 hgradBasis->getValues(hgradBasisGradAtBasisEPoints,Kokkos::subview(basisEPoints, basisEPointsRange(dim, 0), Kokkos::ALL()), OPERATOR_GRAD);
689 hgradBasis->getValues(hgradBasisGradAtTargetEPoints,Kokkos::subview(targetEPoints, targetEPointsRange(dim, 0), Kokkos::ALL()),OPERATOR_GRAD);
690
691 ScalarViewType cellBasisAtBasisEPoints("basisCellAtEPoints",1,numCellDofs, numBasisEPoints, dim);
692 ScalarViewType cellBasisCurlAtCurlEPoints("cellBasisCurlAtCurlEPoints",1,numCellDofs, numBasisCurlEPoints, derDim);
693 ScalarViewType negPartialProjCurl("negPartialProjCurl", numCells, numBasisEPoints, derDim);
694 ScalarViewType negPartialProj("negPartialProj", numCells, numBasisEPoints, dim);
695 ScalarViewType wBasisCurlAtCurlEPoints("weightedBasisCurlAtBasisEPoints",numCells,numCellDofs, numBasisCurlEPoints,derDim);
696 ScalarViewType wBasisCurlBasisAtTargetCurlEPoints("weightedBasisCurlAtTargetCurlEPoints",numCells,numCellDofs, numTargetCurlEPoints,derDim);
697
698 auto targetEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetEvalWeights(dim,0));
699 auto basisEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisEvalWeights(dim,0));
700 auto targetCurlEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getTargetDerivEvalWeights(dim,0));
701 auto basisCurlEWeights = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getBasisDerivEvalWeights(dim,0));
702 ordinal_type offsetBasis = basisEPointsRange(dim, 0).first;
703 ordinal_type offsetBasisCurl = basisCurlEPointsRange(dim, 0).first;
704 ordinal_type offsetTargetCurl = targetCurlEPointsRange(dim, 0).first;
705
706
707 auto hGradTagToOrdinal = Kokkos::create_mirror_view_and_copy(MemSpaceType(), hgradBasis->getAllDofOrdinal());
708
709 using functorTypeCell = FunctorsProjectionTools::ComputeBasisCoeffsOnCell_HCurl<ScalarViewType, decltype(basisEWeights),
710 decltype(computedDofs), decltype(tagToOrdinal)>;
711 Kokkos::parallel_for(policy, functorTypeCell(refBasisCoeffs, negPartialProj, negPartialProjCurl,
712 cellBasisAtBasisEPoints, cellBasisCurlAtCurlEPoints,
713 basisAtBasisEPoints, hgradBasisGradAtBasisEPoints, basisCurlAtBasisCurlEPoints,
714 hgradBasisGradAtTargetEPoints, basisCurlAtTargetCurlEPoints,
715 basisEWeights, basisCurlEWeights, wHgradBasisGradAtBasisEPoints,
716 wBasisCurlAtCurlEPoints, targetEWeights, targetCurlEWeights,
717 wHgradBasisGradAtTargetEPoints,
718 wBasisCurlBasisAtTargetCurlEPoints, computedDofs,
719 tagToOrdinal, hGradTagToOrdinal,
720 numCellDofs, hgradCardinality,
721 offsetBasis, offsetBasisCurl, offsetTargetCurl,
722 numEdgeDofs+numTotalFaceDofs, dim, derDim));
723
724 ScalarViewType cellMassMat_("cellMassMat_", 1, numCellDofs+hgradCardinality, numCellDofs+hgradCardinality),
725 cellRhsMat_("rhsMat_", numCells, numCellDofs+hgradCardinality);
726
727 range_type range_H(0, numCellDofs);
728 range_type range_B(numCellDofs, numCellDofs+hgradCardinality);
729 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(cellMassMat_,Kokkos::ALL(),range_H,range_H), cellBasisCurlAtCurlEPoints, Kokkos::subview(wBasisCurlAtCurlEPoints, std::make_pair(0,1), Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL()) );
730 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(cellMassMat_,Kokkos::ALL(),range_H,range_B), cellBasisAtBasisEPoints, Kokkos::subview(wHgradBasisGradAtBasisEPoints, std::make_pair(0,1), Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL()) );
731 if(dim==3)
732 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(cellRhsMat_,Kokkos::ALL(),range_H), Kokkos::subview(targetCurlAtTargetCurlEPoints,Kokkos::ALL(),cellCurlPointsRange,Kokkos::ALL()), wBasisCurlBasisAtTargetCurlEPoints);
733 else
734 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(cellRhsMat_,Kokkos::ALL(),range_H), Kokkos::subview(targetCurlAtTargetCurlEPoints,Kokkos::ALL(),cellCurlPointsRange), Kokkos::subview(wBasisCurlBasisAtTargetCurlEPoints,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),0));
735 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(cellRhsMat_,Kokkos::ALL(),range_H), negPartialProjCurl, wBasisCurlAtCurlEPoints, true);
736
737 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(cellRhsMat_,Kokkos::ALL(),range_B), Kokkos::subview(targetAtTargetEPoints,Kokkos::ALL(),cellPointsRange,Kokkos::ALL()), wHgradBasisGradAtTargetEPoints);
738 FunctionSpaceTools<DeviceType >::integrate(Kokkos::subview(cellRhsMat_,Kokkos::ALL(),range_B), negPartialProj, wHgradBasisGradAtBasisEPoints, true);
739
740 typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
741 ScalarViewType t_("t",numCells, numCellDofs+hgradCardinality);
742 WorkArrayViewType w_("w",numCells, numCellDofs+hgradCardinality);
743
744 auto cellDofs = Kokkos::subview(tagToOrdinal, dim, 0, Kokkos::ALL());
745 ElemSystem cellSystem( "cellSystem", true);
746 cellSystem.solve(refBasisCoeffs, cellMassMat_, cellRhsMat_, t_, w_, cellDofs, numCellDofs, hgradCardinality);
747
748 delete hgradBasis;
749 }
750
751 OrientationTools<DeviceType>::modifyBasisByOrientationInverse(basisCoeffs, refBasisCoeffs, orts, cellBasis, true);
752}
753
754} // Intrepid2 namespace
755
756#endif
757
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.
Stateless class that acts as a factory for a family of nodal bases (hypercube topologies only at this...
Implementation of the default H(grad)-compatible FEM basis of degree 2 on Hexahedron cell.
Implementation of the default H(grad)-compatible FEM basis of degree n on Quadrilateral cell Implemen...
Implementation of the default H(grad)-compatible Lagrange basis of arbitrary degree on Tetrahedron ce...
Implementation of the default H(grad)-compatible Lagrange basis of arbitrary degree on Triangle cell.
An abstract base class that defines interface for concrete basis implementations for Finite Element (...
virtual KOKKOS_INLINE_FUNCTION void getValues(OutputViewType, const PointViewType, const EOperator, const typename Kokkos::TeamPolicy< ExecutionSpace >::member_type &teamMember, const typename ExecutionSpace::scratch_memory_space &scratchStorage, const ordinal_type subcellDim=-1, const ordinal_type subcellOrdinal=-1) const
Team-level evaluation of basis functions on a reference cell.
ordinal_type getDofCount(const ordinal_type subcDim, const ordinal_type subcOrd) const
DoF count for specified subcell.
ordinal_type getCardinality() const
Returns cardinality of the basis.
const OrdinalTypeArray3DHost getAllDofOrdinal() const
DoF tag to ordinal data structure.
static void getReferenceFaceNormal(RefFaceNormalViewType refFaceNormal, const ordinal_type faceOrd, const shards::CellTopology parentCell)
Computes constant normal vectors to faces of 3D reference cell.
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.
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 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...
An helper class to compute the evaluation points and weights needed for performing projections.
const range_tag getTargetDerivPointsRange() const
Returns the range tag of the target function derivative evaluation points on subcells.
ordinal_type getNumBasisDerivEvalPoints()
Returns number of evaluation points for basis derivatives.
ordinal_type getNumTargetDerivEvalPoints()
Returns number of points where to evaluate the derivatives of the target function.
host_view_type getTargetDerivEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the function derivatives evaluation weights on a subcell.
const range_tag getBasisPointsRange() const
Returns the range tag of the basis evaluation points subcells.
host_view_type getBasisDerivEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the basis derivatives evaluation weights on a subcell.
ordinal_type getNumBasisEvalPoints()
Returns number of basis evaluation points.
host_view_type getTargetEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the function evaluation weights on a subcell.
host_view_type getBasisEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the basis evaluation weights on a subcell.
const range_tag getBasisDerivPointsRange() const
Returns the range tag of the derivative evaluation points on subcell.
const range_tag getTargetPointsRange() const
Returns the range of the target function evaluation points on subcells.
view_type getAllDerivEvalPoints(EvalPointsType type=TARGET) const
Returns all the evaluation points for basis/target derivatives in the cell.
ordinal_type getNumTargetEvalPoints()
Returns number of points where to evaluate the target function.
view_type getAllEvalPoints(EvalPointsType type=TARGET) const
Returns the basis/target evaluation points in the cell.
static void getHCurlBasisCoeffs(Kokkos::DynRankView< basisCoeffsValueType, basisCoeffsProperties... > basisCoeffs, const Kokkos::DynRankView< funValsValueType, funValsProperties... > targetAtEvalPoints, const Kokkos::DynRankView< funValsValueType, funValsProperties... > targetCurlAtCurlEvalPoints, const Kokkos::DynRankView< ortValueType, ortProperties... > cellOrientations, const BasisType *cellBasis, ProjectionStruct< DeviceType, typename BasisType::scalarType > *projStruct)
Computes the basis coefficients of the HCurl projection of the target function.
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...