Intrepid2
Intrepid2_OrientationToolsDefModifyBasis.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 
43 
48 #ifndef __INTREPID2_ORIENTATIONTOOLS_DEF_MODIFY_BASIS_HPP__
49 #define __INTREPID2_ORIENTATIONTOOLS_DEF_MODIFY_BASIS_HPP__
50 
51 // disable clang warnings
52 #if defined (__clang__) && !defined (__INTEL_COMPILER)
53 #pragma clang system_header
54 #endif
55 
57 
58 namespace Intrepid2 {
59 
60  template<typename DT>
61  template<typename elemOrtValueType, class ...elemOrtProperties,
62  typename elemNodeValueType, class ...elemNodeProperties>
63  void
65  getOrientation( Kokkos::DynRankView<elemOrtValueType,elemOrtProperties...> elemOrts,
66  const Kokkos::DynRankView<elemNodeValueType,elemNodeProperties...> elemNodes,
67  const shards::CellTopology cellTopo,
68  bool isSide) {
69  // small meta data modification and it uses shards; let's do this on host
70  auto elemOrtsHost = Kokkos::create_mirror_view(elemOrts);
71  auto elemNodesHost = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), elemNodes);
72 
73  const ordinal_type numCells = elemNodes.extent(0);
74  for (auto cell=0;cell<numCells;++cell) {
75  const auto nodes = Kokkos::subview(elemNodesHost, cell, Kokkos::ALL());
76  elemOrtsHost(cell) = Orientation::getOrientation(cellTopo, nodes, isSide);
77  }
78 
79  Kokkos::deep_copy(elemOrts, elemOrtsHost);
80  }
81 
82  template<typename ortViewType,
83  typename OutputViewType,
84  typename inputViewType,
85  typename o2tViewType,
86  typename t2oViewType,
87  typename dataViewType>
89  ortViewType orts;
90  OutputViewType output;
91  inputViewType input;
92  o2tViewType ordinalToTag;
93  t2oViewType tagToOrdinal;
94 
95  const dataViewType matData;
96  const ordinal_type cellDim, numVerts, numEdges, numFaces, numPoints, dimBasis;
97  const bool leftMultiply;
98  // for simple left-multiplied basis value modification, numPoints is the dimension after the field dimension
99  // for matrix value modification (C,F1,F2), numPoints is F2 when left multiplied, and F1 when right multiplied
100  const bool transpose; // when true, multiply by the transpose of the matrix
101 
102  F_modifyBasisByOrientation(ortViewType orts_,
103  OutputViewType output_,
104  inputViewType input_,
105  o2tViewType ordinalToTag_,
106  t2oViewType tagToOrdinal_,
107  const dataViewType matData_,
108  const ordinal_type cellDim_,
109  const ordinal_type numVerts_,
110  const ordinal_type numEdges_,
111  const ordinal_type numFaces_,
112  const ordinal_type numPoints_,
113  const ordinal_type dimBasis_,
114  const bool leftMultiply_ = true,
115  const bool transpose_ = false)
116  : orts(orts_),
117  output(output_),
118  input(input_),
119  ordinalToTag(ordinalToTag_),
120  tagToOrdinal(tagToOrdinal_),
121  matData(matData_),
122  cellDim(cellDim_),
123  numVerts(numVerts_),
124  numEdges(numEdges_),
125  numFaces(numFaces_),
126  numPoints(numPoints_),
127  dimBasis(dimBasis_),
128  leftMultiply(leftMultiply_),
129  transpose(transpose_)
130  {}
131 
132  KOKKOS_INLINE_FUNCTION
133  void operator()(const ordinal_type cell) const {
134  typedef typename inputViewType::non_const_value_type input_value_type;
135 
136  auto out = Kokkos::subview(output, cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
137  auto in = (input.rank() == output.rank()) ?
138  Kokkos::subview(input, cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL())
139  : Kokkos::subview(input, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
140 
141  // edge transformation
142  ordinal_type existEdgeDofs = 0;
143  if (numEdges > 0) {
144  ordinal_type ortEdges[12];
145  orts(cell).getEdgeOrientation(ortEdges, numEdges);
146 
147  // apply coeff matrix
148  for (ordinal_type edgeId=0;edgeId<numEdges;++edgeId) {
149  const ordinal_type ordEdge = (1 < tagToOrdinal.extent(0) ? (static_cast<size_type>(edgeId) < tagToOrdinal.extent(1) ? tagToOrdinal(1, edgeId, 0) : -1) : -1);
150 
151  if (ordEdge != -1) {
152  existEdgeDofs = 1;
153  const ordinal_type ndofEdge = ordinalToTag(ordEdge, 3);
154  const auto mat = Kokkos::subview(matData,
155  edgeId, ortEdges[edgeId],
156  Kokkos::ALL(), Kokkos::ALL());
157 
158  for (ordinal_type j=0;j<numPoints;++j)
159  for (ordinal_type i=0;i<ndofEdge;++i) {
160  const ordinal_type ii = tagToOrdinal(1, edgeId, i);
161 
162  for (ordinal_type k=0;k<dimBasis;++k) {
163  input_value_type temp = 0.0;
164  for (ordinal_type l=0;l<ndofEdge;++l) {
165  const ordinal_type ll = tagToOrdinal(1, edgeId, l);
166  auto & input_ = leftMultiply ? in(ll, j, k) : in(j, ll, k);
167  auto & mat_il = transpose ? mat(l,i) : mat(i,l);
168  temp += mat_il*input_;
169  }
170  auto & output_ = leftMultiply ? out(ii, j, k) : out(j, ii, k);
171  output_ = temp;
172  }
173  }
174  }
175  }
176  }
177 
178  // face transformation
179  if (numFaces > 0) {
180  ordinal_type ortFaces[12];
181  orts(cell).getFaceOrientation(ortFaces, numFaces);
182 
183  // apply coeff matrix
184  for (ordinal_type faceId=0;faceId<numFaces;++faceId) {
185  const ordinal_type ordFace = (2 < tagToOrdinal.extent(0) ? (static_cast<size_type>(faceId) < tagToOrdinal.extent(1) ? tagToOrdinal(2, faceId, 0) : -1) : -1);
186 
187  if (ordFace != -1) {
188  const ordinal_type ndofFace = ordinalToTag(ordFace, 3);
189  const auto mat = Kokkos::subview(matData,
190  numEdges*existEdgeDofs+faceId, ortFaces[faceId],
191  Kokkos::ALL(), Kokkos::ALL());
192 
193  for (ordinal_type j=0;j<numPoints;++j)
194  for (ordinal_type i=0;i<ndofFace;++i) {
195  const ordinal_type ii = tagToOrdinal(2, faceId, i);
196 
197  for (ordinal_type k=0;k<dimBasis;++k) {
198  input_value_type temp = 0.0;
199  for (ordinal_type l=0;l<ndofFace;++l) {
200  const ordinal_type ll = tagToOrdinal(2, faceId, l);
201  auto & input_ = leftMultiply ? in(ll, j, k) : in(j, ll, k);
202  auto & mat_il = transpose ? mat(l,i) : mat(i,l);
203  temp += mat_il*input_;
204  }
205 
206  auto & output_ = leftMultiply ? out(ii, j, k) : out(j, ii, k);
207  output_ = temp;
208  }
209  }
210  }
211  }
212  }
213 
214  //side orientations
215  ordinal_type faceOrt(0), edgeOrt(0);
216  if(cellDim == 2) orts(cell).getFaceOrientation(&faceOrt, 1);
217  if (faceOrt != 0) {
218  const ordinal_type ordFace = (2 < tagToOrdinal.extent(0) ? (static_cast<size_type>(0) < tagToOrdinal.extent(1) ? tagToOrdinal(2, 0, 0) : -1) : -1);
219 
220  if (ordFace != -1) {
221  const ordinal_type ndofFace = ordinalToTag(ordFace, 3);
222  const auto mat = Kokkos::subview(matData,
223  numEdges*existEdgeDofs, faceOrt,
224  Kokkos::ALL(), Kokkos::ALL());
225 
226  for (ordinal_type j=0;j<numPoints;++j)
227  for (ordinal_type i=0;i<ndofFace;++i) {
228  const ordinal_type ii = tagToOrdinal(2, 0, i);
229 
230  for (ordinal_type k=0;k<dimBasis;++k) {
231  input_value_type temp = 0.0;
232  for (ordinal_type l=0;l<ndofFace;++l) {
233  const ordinal_type ll = tagToOrdinal(2, 0, l);
234  auto & input_ = leftMultiply ? in(ll, j, k) : in(j, ll, k);
235  auto & mat_il = transpose ? mat(l,i) : mat(i,l);
236  temp += mat_il*input_;
237  }
238  auto & output_ = leftMultiply ? out(ii, j, k) : out(j, ii, k);
239  output_ = temp;
240  }
241  }
242  }
243  }
244 
245  if(cellDim == 1) orts(cell).getEdgeOrientation(&edgeOrt, 1);
246  if (edgeOrt != 0) {
247  const ordinal_type ordEdge = (1 < tagToOrdinal.extent(0) ? (static_cast<size_type>(0) < tagToOrdinal.extent(1) ? tagToOrdinal(1, 0, 0) : -1) : -1);
248 
249  if (ordEdge != -1) {
250  const ordinal_type ndofEdge = ordinalToTag(ordEdge, 3);
251  const auto mat = Kokkos::subview(matData,
252  0, edgeOrt,
253  Kokkos::ALL(), Kokkos::ALL());
254 
255  for (ordinal_type j=0;j<numPoints;++j)
256  for (ordinal_type i=0;i<ndofEdge;++i) {
257  const ordinal_type ii = tagToOrdinal(1, 0, i);
258 
259  for (ordinal_type k=0;k<dimBasis;++k) {
260  input_value_type temp = 0.0;
261  for (ordinal_type l=0;l<ndofEdge;++l) {
262  const ordinal_type ll = tagToOrdinal(1, 0, l);
263  auto & input_ = leftMultiply ? in(ll, j, k) : in(j, ll, k);
264  auto & mat_il = transpose ? mat(l,i) : mat(i,l);
265  temp += mat_il*input_;
266  }
267  auto & output_ = leftMultiply ? out(ii, j, k) : out(j, ii, k);
268  output_ = temp;
269  }
270  }
271  }
272  }
273  }
274  };
275 
276  template<typename DT>
277  template<typename outputValueType, class ...outputProperties,
278  typename inputValueType, class ...inputProperties,
279  typename OrientationViewType,
280  typename BasisType>
281  void
283  modifyBasisByOrientation( Kokkos::DynRankView<outputValueType,outputProperties...> output,
284  const Kokkos::DynRankView<inputValueType, inputProperties...> input,
285  const OrientationViewType orts,
286  const BasisType* basis,
287  const bool transpose) {
288 #ifdef HAVE_INTREPID2_DEBUG
289  {
290  if (input.rank() == output.rank())
291  {
292  for (size_type i=0;i<input.rank();++i)
293  INTREPID2_TEST_FOR_EXCEPTION( input.extent(i) != output.extent(i), std::invalid_argument,
294  ">>> ERROR (OrientationTools::modifyBasisByOrientation): Input and output dimensions do not match.");
295  }
296  else if (input.rank() == output.rank() - 1)
297  {
298  for (size_type i=0;i<input.rank();++i)
299  INTREPID2_TEST_FOR_EXCEPTION( input.extent(i) != output.extent(i+1), std::invalid_argument,
300  ">>> ERROR (OrientationTools::modifyBasisByOrientation): Input dimensions must match output dimensions exactly, or else match all but the first dimension (in the case that input does not have a 'cell' dimension).");
301  }
302  else
303  {
304  INTREPID2_TEST_FOR_EXCEPTION(true, std::invalid_argument,
305  ">>> ERROR (OrientationTools::modifyBasisByOrientation): input and output ranks must either match, or input rank must be one less than that of output.")
306  }
307 
308  INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(output.extent(1)) != basis->getCardinality(), std::invalid_argument,
309  ">>> ERROR (OrientationTools::modifyBasisByOrientation): Field dimension of input/output does not match to basis cardinality.");
310  }
311 #endif
312 
313  const shards::CellTopology cellTopo = basis->getBaseCellTopology();
314  const ordinal_type cellDim = cellTopo.getDimension();
315 
316  //Initialize output with values from input
317  if(input.rank() == output.rank())
318  Kokkos::deep_copy(output, input);
319  else
320  RealSpaceTools<DT>::clone(output, input);
321 
322  if ((cellDim < 3) || basis->requireOrientation()) {
323  auto ordinalToTag = Kokkos::create_mirror_view_and_copy(typename DT::memory_space(), basis->getAllDofTags());
324  auto tagToOrdinal = Kokkos::create_mirror_view_and_copy(typename DT::memory_space(), basis->getAllDofOrdinal());
325 
326  const ordinal_type
327  numCells = output.extent(0),
328  //numBasis = output.extent(1),
329  numPoints = output.extent(2),
330  dimBasis = output.extent(3); //returns 1 when output.rank() < 4;
331 
332  const CoeffMatrixDataViewType matData = createCoeffMatrix(basis);
333 
334  ordinal_type numVerts(0), numEdges(0), numFaces(0);
335 
336  if (basis->requireOrientation()) {
337  numVerts = cellTopo.getVertexCount()*ordinal_type(basis->getDofCount(0, 0) > 0);
338  numEdges = cellTopo.getEdgeCount()*ordinal_type(basis->getDofCount(1, 0) > 0);
339  numFaces = cellTopo.getFaceCount()*ordinal_type(basis->getDofCount(2, 0) > 0);
340  }
341 
342  bool leftMultiply = true;
343 
344  const Kokkos::RangePolicy<typename DT::execution_space> policy(0, numCells);
346  <decltype(orts),
347  decltype(output),decltype(input),
348  decltype(ordinalToTag),decltype(tagToOrdinal),
349  decltype(matData)> FunctorType;
350  Kokkos::parallel_for
351  (policy,
352  FunctorType(orts,
353  output, input,
354  ordinalToTag, tagToOrdinal,
355  matData,
356  cellDim, numVerts, numEdges, numFaces,
357  numPoints, dimBasis, leftMultiply, transpose));
358  }
359  }
360 
361  template<typename DT>
362  template<typename outputValueType, class ...outputProperties,
363  typename inputValueType, class ...inputProperties,
364  typename OrientationViewType,
365  typename BasisType>
366  void
368  modifyBasisByOrientationTranspose( Kokkos::DynRankView<outputValueType,outputProperties...> output,
369  const Kokkos::DynRankView<inputValueType, inputProperties...> input,
370  const OrientationViewType orts,
371  const BasisType* basis ) {
372  bool transpose = true;
373  modifyBasisByOrientation(output, input, orts, basis, transpose);
374  }
375 
376  template<typename DT>
377  template<typename outputValueType, class ...outputProperties,
378  typename inputValueType, class ...inputProperties,
379  typename OrientationViewType,
380  typename BasisType>
381  void
383  modifyBasisByOrientationInverse( Kokkos::DynRankView<outputValueType,outputProperties...> output,
384  const Kokkos::DynRankView<inputValueType, inputProperties...> input,
385  const OrientationViewType orts,
386  const BasisType* basis,
387  const bool transpose ) {
388  #ifdef HAVE_INTREPID2_DEBUG
389  {
390  if (input.rank() == output.rank())
391  {
392  for (size_type i=0;i<input.rank();++i)
393  INTREPID2_TEST_FOR_EXCEPTION( input.extent(i) != output.extent(i), std::invalid_argument,
394  ">>> ERROR (OrientationTools::modifyBasisByOrientation): Input and output dimensions do not match.");
395  }
396  else if (input.rank() == output.rank() - 1)
397  {
398  for (size_type i=0;i<input.rank();++i)
399  INTREPID2_TEST_FOR_EXCEPTION( input.extent(i) != output.extent(i+1), std::invalid_argument,
400  ">>> ERROR (OrientationTools::modifyBasisByOrientation): Input dimensions must match output dimensions exactly, or else match all but the first dimension (in the case that input does not have a 'cell' dimension).");
401  }
402  else
403  {
404  INTREPID2_TEST_FOR_EXCEPTION(true, std::invalid_argument,
405  ">>> ERROR (OrientationTools::modifyBasisByOrientation): input and output ranks must either match, or input rank must be one less than that of output.")
406  }
407 
408  INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(output.extent(1)) != basis->getCardinality(), std::invalid_argument,
409  ">>> ERROR (OrientationTools::modifyBasisByOrientation): Field dimension of input/output does not match to basis cardinality.");
410  }
411  #endif
412 
413  const shards::CellTopology cellTopo = basis->getBaseCellTopology();
414  const ordinal_type cellDim = cellTopo.getDimension();
415 
416  //Initialize output with values from input
417  if(input.rank() == output.rank())
418  Kokkos::deep_copy(output, input);
419  else
420  RealSpaceTools<DT>::clone(output, input);
421 
422  if ((cellDim < 3) || basis->requireOrientation()) {
423  auto ordinalToTag = Kokkos::create_mirror_view_and_copy(typename DT::memory_space(), basis->getAllDofTags());
424  auto tagToOrdinal = Kokkos::create_mirror_view_and_copy(typename DT::memory_space(), basis->getAllDofOrdinal());
425 
426  const ordinal_type
427  numCells = output.extent(0),
428  //numBasis = output.extent(1),
429  numPoints = output.extent(2),
430  dimBasis = output.extent(3); //returns 1 when output.rank() < 4;
431 
432  const CoeffMatrixDataViewType matData = createInvCoeffMatrix(basis);
433 
434  ordinal_type numVerts(0), numEdges(0), numFaces(0);
435 
436  if (basis->requireOrientation()) {
437  numVerts = cellTopo.getVertexCount()*ordinal_type(basis->getDofCount(0, 0) > 0);
438  numEdges = cellTopo.getEdgeCount()*ordinal_type(basis->getDofCount(1, 0) > 0);
439  numFaces = cellTopo.getFaceCount()*ordinal_type(basis->getDofCount(2, 0) > 0);
440  }
441 
442  bool leftMultiply = true;
443 
444  const Kokkos::RangePolicy<typename DT::execution_space> policy(0, numCells);
446  <decltype(orts),
447  decltype(output),decltype(input),
448  decltype(ordinalToTag),decltype(tagToOrdinal),
449  decltype(matData)> FunctorType;
450  Kokkos::parallel_for
451  (policy,
452  FunctorType(orts,
453  output, input,
454  ordinalToTag, tagToOrdinal,
455  matData,
456  cellDim, numVerts, numEdges, numFaces,
457  numPoints, dimBasis, leftMultiply, transpose));
458  }
459  }
460 
461  template<typename DT>
462  template<typename outputValueType, class ...outputProperties,
463  typename inputValueType, class ...inputProperties,
464  typename OrientationViewType,
465  typename BasisTypeLeft,
466  typename BasisTypeRight>
467  void
469  modifyMatrixByOrientation(Kokkos::DynRankView<outputValueType,outputProperties...> output,
470  const Kokkos::DynRankView<inputValueType, inputProperties...> input,
471  const OrientationViewType orts,
472  const BasisTypeLeft* basisLeft,
473  const BasisTypeRight* basisRight)
474  {
475  const ordinal_type numCells = output.extent(0);
476  const ordinal_type numFieldsLeft = basisLeft->getCardinality();
477  const ordinal_type numFieldsRight = basisRight->getCardinality();
478 #ifdef HAVE_INTREPID2_DEBUG
479  {
480  if (input.rank() == output.rank())
481  {
482  for (size_type i=0;i<input.rank();++i)
483  INTREPID2_TEST_FOR_EXCEPTION( input.extent(i) != output.extent(i), std::invalid_argument,
484  ">>> ERROR (OrientationTools::modifyMatrixByOrientation): Input and output dimensions do not match.");
485  }
486  else if (input.rank() == output.rank() - 1)
487  {
488  for (size_type i=0;i<input.rank();++i)
489  INTREPID2_TEST_FOR_EXCEPTION( input.extent(i) != output.extent(i+1), std::invalid_argument,
490  ">>> ERROR (OrientationTools::modifyMatrixByOrientation): Input dimensions must match output dimensions exactly, or else match all but the first dimension (in the case that input does not have a 'cell' dimension).");
491  }
492  else
493  {
494  INTREPID2_TEST_FOR_EXCEPTION(true, std::invalid_argument,
495  ">>> ERROR (OrientationTools::modifyMatrixByOrientation): input and output ranks must either match, or input rank must be one less than that of output.")
496  }
497 
498  INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(output.extent(1)) != numFieldsLeft, std::invalid_argument,
499  ">>> ERROR (OrientationTools::modifyMatrixByOrientation): First field dimension of input/output does not match left basis cardinality.");
500  INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(output.extent(2)) != numFieldsRight, std::invalid_argument,
501  ">>> ERROR (OrientationTools::modifyMatrixByOrientation): Second field dimension of input/output does not match right basis cardinality.");
502  INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(output.extent(3)) != 1, std::invalid_argument,
503  ">>> ERROR (OrientationTools::modifyMatrixByOrientation): Third dimension of output must be 1.");
504 
505  }
506 #endif
507  const shards::CellTopology cellTopo = basisLeft->getBaseCellTopology();
508  const ordinal_type cellDim = cellTopo.getDimension();
509 
510  // apply orientations on left
511  decltype(output) outputLeft("temp view - output from left application", numCells, numFieldsLeft, numFieldsRight);
512 
513  //Initialize outputLeft with values from input
514  if(input.rank() == output.rank())
515  Kokkos::deep_copy(outputLeft, input);
516  else
517  RealSpaceTools<DT>::clone(outputLeft, input);
518 
519  if ((cellDim < 3) || basisLeft->requireOrientation()) {
520  bool leftMultiply = true;
521  auto ordinalToTag = Kokkos::create_mirror_view_and_copy(typename DT::memory_space(), basisLeft->getAllDofTags());
522  auto tagToOrdinal = Kokkos::create_mirror_view_and_copy(typename DT::memory_space(), basisLeft->getAllDofOrdinal());
523 
524  const ordinal_type
525  numOtherFields = output.extent(2),
526  dimBasis = output.extent(3); //returns 1 when output.rank() < 4;
527 
528  const CoeffMatrixDataViewType matData = createCoeffMatrix(basisLeft);
529 
530  ordinal_type numVerts(0), numEdges(0), numFaces(0);
531 
532  if (basisLeft->requireOrientation()) {
533  numVerts = cellTopo.getVertexCount()*ordinal_type(basisLeft->getDofCount(0, 0) > 0);
534  numEdges = cellTopo.getEdgeCount()*ordinal_type(basisLeft->getDofCount(1, 0) > 0);
535  numFaces = cellTopo.getFaceCount()*ordinal_type(basisLeft->getDofCount(2, 0) > 0);
536  }
537 
538  const Kokkos::RangePolicy<typename DT::execution_space> policy(0, numCells);
540  <decltype(orts),
541  decltype(outputLeft),decltype(input),
542  decltype(ordinalToTag),decltype(tagToOrdinal),
543  decltype(matData)> FunctorType;
544  Kokkos::parallel_for
545  (policy,
546  FunctorType(orts,
547  outputLeft, input,
548  ordinalToTag, tagToOrdinal,
549  matData,
550  cellDim, numVerts, numEdges, numFaces,
551  numOtherFields, dimBasis, leftMultiply));
552  }
553 
554  // apply orientations on right
555  //Initialize output with values from outputLeft
556  Kokkos::deep_copy(output, outputLeft);
557  if ((cellDim < 3) || basisRight->requireOrientation()) {
558  bool leftMultiply = false;
559  auto ordinalToTag = Kokkos::create_mirror_view_and_copy(typename DT::memory_space(), basisRight->getAllDofTags());
560  auto tagToOrdinal = Kokkos::create_mirror_view_and_copy(typename DT::memory_space(), basisRight->getAllDofOrdinal());
561 
562  const ordinal_type
563  numOtherFields = output.extent(1),
564  dimBasis = output.extent(3); //returns 1 when output.rank() < 4;
565 
566  const CoeffMatrixDataViewType matData = createCoeffMatrix(basisRight);
567 
568  ordinal_type numVerts(0), numEdges(0), numFaces(0);
569 
570  if (basisRight->requireOrientation()) {
571  numVerts = cellTopo.getVertexCount()*ordinal_type(basisRight->getDofCount(0, 0) > 0);
572  numEdges = cellTopo.getEdgeCount()*ordinal_type(basisRight->getDofCount(1, 0) > 0);
573  numFaces = cellTopo.getFaceCount()*ordinal_type(basisRight->getDofCount(2, 0) > 0);
574  }
575 
576  const Kokkos::RangePolicy<typename DT::execution_space> policy(0, numCells);
578  <decltype(orts),
579  decltype(output),decltype(outputLeft),
580  decltype(ordinalToTag),decltype(tagToOrdinal),
581  decltype(matData)> FunctorType;
582  Kokkos::parallel_for
583  (policy,
584  FunctorType(orts,
585  output, outputLeft,
586  ordinalToTag, tagToOrdinal,
587  matData,
588  cellDim, numVerts, numEdges, numFaces,
589  numOtherFields, dimBasis, leftMultiply));
590  }
591  }
592 
593 } // namespace Intrepid2
594 
595 #endif
Kokkos::View< double ****, DeviceType > CoeffMatrixDataViewType
subcell ordinal, orientation, matrix m x n
static void clone(Kokkos::DynRankView< outputValueType, outputProperties...> output, const Kokkos::DynRankView< inputValueType, inputProperties...> input)
Clone input array.
static void modifyMatrixByOrientation(Kokkos::DynRankView< outputValueType, outputProperties...> output, const Kokkos::DynRankView< inputValueType, inputProperties...> input, const OrientationViewType orts, const BasisTypeLeft *basisLeft, const BasisTypeRight *basisRight)
Modify an assembled (C,F1,F2) matrix according to orientation of the cells.
Header file for the Intrepid2::Orientation class.
static void modifyBasisByOrientation(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.
static void getOrientation(Kokkos::DynRankView< elemOrtValueType, elemOrtProperties...> elemOrts, const Kokkos::DynRankView< elemNodeValueType, elemNodeProperties...> elemNodes, const shards::CellTopology cellTopo, bool isSide=false)
Compute orientations of cells in a workset.
static void modifyBasisByOrientationTranspose(Kokkos::DynRankView< outputValueType, outputProperties...> output, const Kokkos::DynRankView< inputValueType, inputProperties...> input, const OrientationViewType orts, const BasisType *basis)
Modify basis due to orientation, applying the transpose of the operator applied in modifyBasisByOrien...
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...