Intrepid2
Intrepid2_FunctionSpaceToolsDef.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Intrepid2 Package
5 // Copyright (2007) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Kyungjoo Kim (kyukim@sandia.gov), or
38 // Mauro Perego (mperego@sandia.gov)
39 //
40 // ************************************************************************
41 // @HEADER
42 
49 #ifndef __INTREPID2_FUNCTIONSPACETOOLS_DEF_HPP__
50 #define __INTREPID2_FUNCTIONSPACETOOLS_DEF_HPP__
51 
52 namespace Intrepid2 {
53 
54  // ------------------------------------------------------------------------------------
55  template<typename SpT>
56  template<typename outputValueType, class ...outputProperties,
57  typename inputValueType, class ...inputProperties>
58  void
60  HGRADtransformVALUE( Kokkos::DynRankView<outputValueType,outputProperties...> output,
61  const Kokkos::DynRankView<inputValueType, inputProperties...> input ) {
62  ArrayTools<SpT>::cloneFields(output, input);
63  }
64 
65  // ------------------------------------------------------------------------------------
66 
67  namespace FunctorFunctionSpaceTools {
71  template <typename outputViewType,
72  typename jacInverseViewType,
73  typename inputViewType,
74  ordinal_type spaceDim>
76  outputViewType _output;
77  const jacInverseViewType _jacInverse;
78  const inputViewType _input;
79 
80  // output CPDD, left CPDD or PDD, right FPD
81  KOKKOS_INLINE_FUNCTION
82  F_HGRADtransformGRAD(outputViewType output_,
83  jacInverseViewType jacInverse_,
84  inputViewType input_)
85  : _output(output_),
86  _jacInverse(jacInverse_),
87  _input(input_) {}
88 
89  KOKKOS_INLINE_FUNCTION
90  void operator()(const ordinal_type cl,
91  const ordinal_type bf,
92  const ordinal_type pt) const {
93  /* */ auto y = Kokkos::subview(_output, cl, bf, pt, Kokkos::ALL());
94  const auto A = Kokkos::subview(_jacInverse, cl, pt, Kokkos::ALL(), Kokkos::ALL());
95  const auto x = Kokkos::subview(_input, bf, pt, Kokkos::ALL());
96 
97  if (spaceDim == 2) {
98  Kernels::Serial::matvec_trans_product_d2( y, A, x );
99  } else {
100  Kernels::Serial::matvec_trans_product_d3( y, A, x );
101  }
102  }
103  };
104  }
105 
106  template<typename SpT>
107  template<typename outputValValueType, class ...outputValProperties,
108  typename jacobianInverseValueType, class ...jacobianInverseProperties,
109  typename inputValValueType, class ...inputValProperties>
110  void
112  HGRADtransformGRAD( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
113  const Kokkos::DynRankView<jacobianInverseValueType,jacobianInverseProperties...> jacobianInverse,
114  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
115 #ifdef HAVE_INTREPID2_DEBUG
116  {
117  INTREPID2_TEST_FOR_EXCEPTION( inputVals.rank() != 3 ||
118  jacobianInverse.rank() != 4 ||
119  outputVals.rank() != 4, std::invalid_argument,
120  ">>> ERROR (FunctionSpaceTools::HGRADtransformGRAD): Ranks are not compatible.");
121  INTREPID2_TEST_FOR_EXCEPTION( outputVals.extent(0) != jacobianInverse.extent(0), std::invalid_argument,
122  ">>> ERROR (FunctionSpaceTools::HGRADtransformGRAD): Cell dimension does not match.");
123  INTREPID2_TEST_FOR_EXCEPTION( outputVals.extent(1) != inputVals.extent(0), std::invalid_argument,
124  ">>> ERROR (FunctionSpaceTools::HGRADtransformGRAD): Field dimension does not match.");
125  INTREPID2_TEST_FOR_EXCEPTION( outputVals.extent(2) != inputVals.extent(1) ||
126  jacobianInverse.extent(1) != inputVals.extent(1), std::invalid_argument,
127  ">>> ERROR (FunctionSpaceTools::HGRADtransformGRAD): Point dimension does not match.");
128  const auto spaceDim = outputVals.extent(3);
129  INTREPID2_TEST_FOR_EXCEPTION( jacobianInverse.extent(2) != spaceDim ||
130  jacobianInverse.extent(3) != spaceDim ||
131  inputVals.extent(2) != spaceDim , std::invalid_argument,
132  ">>> ERROR (FunctionSpaceTools::HGRADtransformGRAD): matvec dimensions are not compatible.");
133  }
134 #endif
135  ArrayTools<SpT>::matvecProductDataField(outputVals, jacobianInverse, inputVals, 'T');
136 
137  // this modification is for 2d and 3d (not 1d)
138  // this is an attempt to measure the overhead of subview of dynrankview.
139 
140  // typedef Kokkos::DynRankView<outputValValueType,outputValProperties...> outputViewType;
141  // typedef const Kokkos::DynRankView<jacobianInverseValueType,jacobianInverseProperties...> jacInverseViewType;
142  // typedef const Kokkos::DynRankView<inputValValueType,inputValProperties...> inputViewType;
143 
144  // const ordinal_type
145  // C = outputVals.extent(0),
146  // F = outputVals.extent(1),
147  // P = outputVals.extent(2);
148 
149  // using range_policy_type = Kokkos::Experimental::MDRangePolicy
150  // < SpT, Kokkos::Experimental::Rank<3>, Kokkos::IndexType<ordinal_type> >;
151  // range_policy_type policy( { 0, 0, 0 },
152  // { C, F, P } );
153 
154  // const ordinal_type spaceDim = inputVals.extent(2);
155  // switch (spaceDim) {
156  // case 2: {
157  // typedef FunctorFunctionSpaceTools::F_HGRADtransformGRAD<outputViewType, jacInverseViewType, inputViewType, 2> FunctorType;
158  // Kokkos::parallel_for( policy, FunctorType(outputVals, jacobianInverse, inputVals) );
159  // break;
160  // }
161  // case 3: {
162  // typedef FunctorFunctionSpaceTools::F_HGRADtransformGRAD<outputViewType, jacInverseViewType, inputViewType, 3> FunctorType;
163  // Kokkos::parallel_for( policy, FunctorType(outputVals, jacobianInverse, inputVals) );
164  // break;
165  // }
166  // default: {
167  // INTREPID2_TEST_FOR_EXCEPTION( true, std::invalid_argument,
168  // ">>> ERROR (FunctionSpaceTools::HGRADtransformGRAD): spaceDim is not 2 or 3.");
169  // break;
170  // }
171  // }
172  }
173 
174  // ------------------------------------------------------------------------------------
175 
176  template<typename SpT>
177  template<typename outputValValueType, class ...outputValProperties,
178  typename jacobianInverseValueType, class ...jacobianInverseProperties,
179  typename inputValValueType, class ...inputValProperties>
180  void
182  HCURLtransformVALUE( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
183  const Kokkos::DynRankView<jacobianInverseValueType,jacobianInverseProperties...> jacobianInverse,
184  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
185  ArrayTools<SpT>::matvecProductDataField(outputVals, jacobianInverse, inputVals, 'T');
186  }
187 
188  // ------------------------------------------------------------------------------------
189 
190  template<typename SpT>
191  template<typename outputValValueType, class ...outputValProperties,
192  typename jacobianValueType, class ...jacobianProperties,
193  typename jacobianDetValueType, class ...jacobianDetProperties,
194  typename inputValValueType, class ...inputValProperties>
195  void
197  HCURLtransformCURL( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
198  const Kokkos::DynRankView<jacobianValueType, jacobianProperties...> jacobian,
199  const Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
200  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
201  ArrayTools<SpT>::matvecProductDataField(outputVals, jacobian, inputVals, 'N');
202  ArrayTools<SpT>::scalarMultiplyDataField(outputVals, jacobianDet, outputVals, true);
203  }
204 
205  // ------------------------------------------------------------------------------------
206 
207  template<typename SpT>
208  template<typename outputValValueType, class ...outputValProperties,
209  typename jacobianValueType, class ...jacobianProperties,
210  typename jacobianDetValueType, class ...jacobianDetProperties,
211  typename inputValValueType, class ...inputValProperties>
212  void
214  HDIVtransformVALUE( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
215  const Kokkos::DynRankView<jacobianValueType, jacobianProperties...> jacobian,
216  const Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
217  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
218  ArrayTools<SpT>::matvecProductDataField(outputVals, jacobian, inputVals, 'N');
219  ArrayTools<SpT>::scalarMultiplyDataField(outputVals, jacobianDet, outputVals, true);
220  }
221 
222  // ------------------------------------------------------------------------------------
223 
224  template<typename SpT>
225  template<typename outputValValueType, class ...outputValProperties,
226  typename jacobianDetValueType, class ...jacobianDetProperties,
227  typename inputValValueType, class ...inputValProperties>
228  void
230  HDIVtransformDIV( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
231  const Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
232  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
233  ArrayTools<SpT>::scalarMultiplyDataField(outputVals, jacobianDet, inputVals, true);
234  }
235 
236  // ------------------------------------------------------------------------------------
237 
238  template<typename SpT>
239  template<typename outputValValueType, class ...outputValProperties,
240  typename jacobianDetValueType, class ...jacobianDetProperties,
241  typename inputValValueType, class ...inputValProperties>
242  void
244  HVOLtransformVALUE( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
245  const Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
246  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
247  ArrayTools<SpT>::scalarMultiplyDataField(outputVals, jacobianDet, inputVals, true);
248  }
249 
250  // ------------------------------------------------------------------------------------
251 
252  template<typename SpT>
253  template<typename outputValueValueType, class ...outputValueProperties,
254  typename leftValueValueType, class ...leftValueProperties,
255  typename rightValueValueType, class ...rightValueProperties>
256  void
258  integrate( Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
259  const Kokkos::DynRankView<leftValueValueType, leftValueProperties...> leftValues,
260  const Kokkos::DynRankView<rightValueValueType, rightValueProperties...> rightValues,
261  const bool sumInto ) {
262 
263 #ifdef HAVE_INTREPID2_DEBUG
264  {
265  INTREPID2_TEST_FOR_EXCEPTION( leftValues.rank() < 2 ||
266  leftValues.rank() > 4, std::invalid_argument,
267  ">>> ERROR (FunctionSpaceTools::integrate): Left data must have rank 2, 3 or 4.");
268  INTREPID2_TEST_FOR_EXCEPTION( outputValues.rank() < 1 ||
269  outputValues.rank() > 3, std::invalid_argument,
270  ">>> ERROR (FunctionSpaceTools::integrate): Output values must have rank 1, 2 or 3.");
271  }
272 #endif
273 
274  const ordinal_type outRank = outputValues.rank();
275  const ordinal_type leftRank = leftValues.rank();
276  const ordinal_type mode = outRank*10 + leftRank;
277 
278  switch (mode) {
279  // DataData
280  case 12:
282  leftValues,
283  rightValues,
284  sumInto );
285  break;
286  case 13:
288  leftValues,
289  rightValues,
290  sumInto );
291  break;
292  case 14:
294  leftValues,
295  rightValues,
296  sumInto );
297  break;
298 
299  // DataField
300  case 22:
302  leftValues,
303  rightValues,
304  sumInto );
305  break;
306  case 23:
308  leftValues,
309  rightValues,
310  sumInto );
311  break;
312  case 24:
314  leftValues,
315  rightValues,
316  sumInto );
317  break;
318 
319  // FieldField
320  case 33:
322  leftValues,
323  rightValues,
324  sumInto );
325  break;
326  case 34:
328  leftValues,
329  rightValues,
330  sumInto );
331  break;
332  case 35:
334  leftValues,
335  rightValues,
336  sumInto );
337  break;
338  default: {
339  INTREPID2_TEST_FOR_EXCEPTION( outRank < 1 || outRank > 3, std::runtime_error,
340  ">>> ERROR (FunctionSpaceTools::integrate): outRank must be 1,2, or 3.");
341  INTREPID2_TEST_FOR_EXCEPTION( leftRank < 2 || leftRank > 4, std::runtime_error,
342  ">>> ERROR (FunctionSpaceTools::integrate): leftRank must be 1,2, 3 or 4.");
343  }
344  }
345  }
346 
347  // ------------------------------------------------------------------------------------
348  namespace FunctorFunctionSpaceTools {
352  template<typename outputValViewType,
353  typename inputDetViewType,
354  typename inputWeightViewType>
356  outputValViewType _outputVals;
357  const inputDetViewType _inputDet;
358  const inputWeightViewType _inputWeight;
359 
360  KOKKOS_INLINE_FUNCTION
361  F_computeCellMeasure( outputValViewType outputVals_,
362  inputDetViewType inputDet_,
363  inputWeightViewType inputWeight_)
364  : _outputVals(outputVals_),
365  _inputDet(inputDet_),
366  _inputWeight(inputWeight_) {}
367 
368  typedef ordinal_type value_type;
369 
370 // KOKKOS_INLINE_FUNCTION
371 // void init( value_type &dst ) const {
372 // dst = false;
373 // }
374 
375 // KOKKOS_INLINE_FUNCTION
376 // void join( volatile value_type &dst,
377 // const volatile value_type &src ) const {
378 // dst |= src;
379 // }
380 
381  KOKKOS_INLINE_FUNCTION
382  void operator()(const size_type cl, value_type &dst) const {
383  // negative jacobian check
384  const bool hasNegativeDet = (_inputDet(cl, 0) < 0.0);
385  dst |= hasNegativeDet;
386 
387  // make it absolute
388  const auto sign = (hasNegativeDet ? -1.0 : 1.0);
389  const ordinal_type pt_end = _outputVals.extent(1);
390  for (ordinal_type pt=0;pt<pt_end;++pt)
391  _outputVals(cl, pt) = sign*_inputDet(cl, pt)*_inputWeight(pt);
392  }
393  };
394  }
395 
396  template<typename SpT>
397  template<typename outputValValueType, class ...outputValProperties,
398  typename inputDetValueType, class ...inputDetProperties,
399  typename inputWeightValueType, class ...inputWeightProperties>
400  bool
402  computeCellMeasure( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
403  const Kokkos::DynRankView<inputDetValueType, inputDetProperties...> inputDet,
404  const Kokkos::DynRankView<inputWeightValueType,inputWeightProperties...> inputWeights ) {
405 #ifdef HAVE_INTREPID2_DEBUG
406  {
407  INTREPID2_TEST_FOR_EXCEPTION( inputDet.rank() != 2 ||
408  inputWeights.rank() != 1 ||
409  outputVals.rank() != 2, std::invalid_argument,
410  ">>> ERROR (FunctionSpaceTools::computeCellMeasure): Ranks are not compatible.");
411  INTREPID2_TEST_FOR_EXCEPTION( outputVals.extent(0) != inputDet.extent(0), std::invalid_argument,
412  ">>> ERROR (FunctionSpaceTools::computeCellMeasure): Cell dimension does not match.");
413  INTREPID2_TEST_FOR_EXCEPTION( inputDet.extent(1) != outputVals.extent(1) ||
414  inputWeights.extent(0) != outputVals.extent(1), std::invalid_argument,
415  ">>> ERROR (FunctionSpaceTools::computeCellMeasure): Point dimension does not match.");
416  }
417 #endif
418  typedef Kokkos::DynRankView<outputValValueType, outputValProperties...> outputValViewType;
419  typedef Kokkos::DynRankView<inputDetValueType, inputDetProperties...> inputDetViewType;
420  typedef Kokkos::DynRankView<inputWeightValueType,inputWeightProperties...> inputWeightViewType;
422  <outputValViewType,inputDetViewType,inputWeightViewType> FunctorType;
423 
424  const ordinal_type C = inputDet.extent(0);
425  Kokkos::RangePolicy<SpT,Kokkos::Schedule<Kokkos::Static> > policy(0, C);
426 
427  typename FunctorType::value_type hasNegativeDet = false;
428  Kokkos::parallel_reduce( policy, FunctorType(outputVals, inputDet, inputWeights), hasNegativeDet );
429 
430  return hasNegativeDet;
431  }
432 
433  // ------------------------------------------------------------------------------------
434 
435  template<typename SpT>
436  template<typename outputValValueType, class ...outputValProperties,
437  typename inputJacValueType, class ...inputJacProperties,
438  typename inputWeightValueType, class ...inputWeightProperties,
439  typename scratchValueType, class ...scratchProperties>
440  void
442  computeFaceMeasure( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
443  const Kokkos::DynRankView<inputJacValueType, inputJacProperties...> inputJac,
444  const Kokkos::DynRankView<inputWeightValueType,inputWeightProperties...> inputWeights,
445  const ordinal_type whichFace,
446  const shards::CellTopology parentCell,
447  const Kokkos::DynRankView<scratchValueType, scratchProperties...> scratch ) {
448 #ifdef HAVE_INTREPID2_DEBUG
449  INTREPID2_TEST_FOR_EXCEPTION( inputJac.rank() != 4, std::invalid_argument,
450  ">>> ERROR (FunctionSpaceTools::computeFaceMeasure): Input Jacobian container must have rank 4.");
451  INTREPID2_TEST_FOR_EXCEPTION( scratch.rank() != 1, std::invalid_argument,
452  ">>> ERROR (FunctionSpaceTools::computeFaceMeasure): Scratch view imust have rank 1.");
453  INTREPID2_TEST_FOR_EXCEPTION( scratch.span() < inputJac.span(), std::invalid_argument,
454  ">>> ERROR (FunctionSpaceTools::computeFaceMeasure): Scratch storage must be greater than or equal to inputJac's one.");
455 #endif
456 
457  // face normals (reshape scratch)
458  // Kokkos::DynRankView<scratchValueType,scratchProperties...> faceNormals(scratch.data(),
459  // inputJac.extent(0),
460  // inputJac.extent(1),
461  // inputJac.extent(2));
462  auto vcprop = Kokkos::common_view_alloc_prop(scratch);
463  //typedef Kokkos::DynRankView<scratchValueType, typename decltype(scratch)::memory_space> viewType;
464  typedef Kokkos::DynRankView<scratchValueType, SpT> viewType;
465  viewType faceNormals(Kokkos::view_wrap(scratch.data(), vcprop),
466  inputJac.extent(0),
467  inputJac.extent(1),
468  inputJac.extent(2));
469 
470  // compute normals
471  CellTools<SpT>::getPhysicalFaceNormals(faceNormals, inputJac, whichFace, parentCell);
472 
473  // compute lenghts of normals
474  RealSpaceTools<SpT>::vectorNorm(outputVals, faceNormals, NORM_TWO);
475 
476  // multiply with weights
477  ArrayTools<SpT>::scalarMultiplyDataData(outputVals, outputVals, inputWeights);
478  }
479 
480  // ------------------------------------------------------------------------------------
481 
482  template<typename SpT>
483  template<typename outputValValueType, class ...outputValProperties,
484  typename inputJacValueType, class ...inputJacProperties,
485  typename inputWeightValueType, class ...inputWeightProperties,
486  typename scratchValueType, class ...scratchProperties>
487  void
489  computeEdgeMeasure( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
490  const Kokkos::DynRankView<inputJacValueType, inputJacProperties...> inputJac,
491  const Kokkos::DynRankView<inputWeightValueType,inputWeightProperties...> inputWeights,
492  const ordinal_type whichEdge,
493  const shards::CellTopology parentCell,
494  const Kokkos::DynRankView<scratchValueType, scratchProperties...> scratch ) {
495 #ifdef HAVE_INTREPID2_DEBUG
496  INTREPID2_TEST_FOR_EXCEPTION( (inputJac.rank() != 4), std::invalid_argument,
497  ">>> ERROR (FunctionSpaceTools::computeEdgeMeasure): Input Jacobian container must have rank 4.");
498  INTREPID2_TEST_FOR_EXCEPTION( scratch.rank() != 1, std::invalid_argument,
499  ">>> ERROR (FunctionSpaceTools::computeEdgeMeasure): Scratch view must have a rank 1.");
500  INTREPID2_TEST_FOR_EXCEPTION( scratch.span() < inputJac.span(), std::invalid_argument,
501  ">>> ERROR (FunctionSpaceTools::computeEdgeMeasure): Scratch storage must be greater than or equal to inputJac'one.");
502 #endif
503 
504  // edge tangents (reshape scratch)
505  // Kokkos::DynRankView<scratchValueType,scratchProperties...> edgeTangents(scratch.data(),
506  // inputJac.extent(0),
507  // inputJac.extent(1),
508  // inputJac.extent(2));
509  auto vcprop = Kokkos::common_view_alloc_prop(scratch);
510  //typedef Kokkos::DynRankView<scratchValueType, typename decltype(scratch)::memory_space> viewType;
511  typedef Kokkos::DynRankView<scratchValueType, SpT> viewType;
512  viewType edgeTangents(Kokkos::view_wrap(scratch.data(), vcprop),
513  inputJac.extent(0),
514  inputJac.extent(1),
515  inputJac.extent(2));
516 
517  // compute normals
518  CellTools<SpT>::getPhysicalEdgeTangents(edgeTangents, inputJac, whichEdge, parentCell);
519 
520  // compute lenghts of tangents
521  RealSpaceTools<SpT>::vectorNorm(outputVals, edgeTangents, NORM_TWO);
522 
523  // multiply with weights
524  ArrayTools<SpT>::scalarMultiplyDataData(outputVals, outputVals, inputWeights);
525  }
526 
527  // ------------------------------------------------------------------------------------
528 
529  template<typename SpT>
530  template<typename outputValValueType, class ...outputValProperties,
531  typename inputMeasureValueType, class ...inputMeasureProperties,
532  typename inputValValueType, class ...inputValProperteis>
533  void
535  multiplyMeasure( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
536  const Kokkos::DynRankView<inputMeasureValueType,inputMeasureProperties...> inputMeasure,
537  const Kokkos::DynRankView<inputValValueType, inputValProperteis...> inputVals ) {
538  scalarMultiplyDataField( outputVals,
539  inputMeasure,
540  inputVals );
541  }
542 
543  // ------------------------------------------------------------------------------------
544 
545  template<typename SpT>
546  template<typename outputFieldValueType, class ...outputFieldProperties,
547  typename inputDataValueType, class ...inputDataProperties,
548  typename inputFieldValueType, class ...inputFieldProperties>
549  void
551  scalarMultiplyDataField( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
552  const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
553  const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields,
554  const bool reciprocal ) {
556  inputData,
557  inputFields,
558  reciprocal );
559  }
560 
561  // ------------------------------------------------------------------------------------
562 
563  template<typename SpT>
564  template<typename outputDataValuetype, class ...outputDataProperties,
565  typename inputDataLeftValueType, class ...inputDataLeftProperties,
566  typename inputDataRightValueType, class ...inputDataRightProperties>
567  void
569  scalarMultiplyDataData( Kokkos::DynRankView<outputDataValuetype, outputDataProperties...> outputData,
570  const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
571  const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight,
572  const bool reciprocal ) {
574  inputDataLeft,
575  inputDataRight,
576  reciprocal );
577  }
578 
579  // ------------------------------------------------------------------------------------
580 
581  template<typename SpT>
582  template<typename outputFieldValueType, class ...outputFieldProperties,
583  typename inputDataValueType, class ...inputDataProperties,
584  typename inputFieldValueType, class ...inputFieldProperties>
585  void
587  dotMultiplyDataField( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
588  const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
589  const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields ) {
591  inputData,
592  inputFields );
593  }
594 
595  // ------------------------------------------------------------------------------------
596 
597  template<typename SpT>
598  template<typename outputDataValueType, class ...outputDataProperties,
599  typename inputDataLeftValueType, class ...inputDataLeftProperties,
600  typename inputDataRightValueType, class ...inputDataRightProperties>
601  void
603  dotMultiplyDataData( Kokkos::DynRankView<outputDataValueType, outputDataProperties...> outputData,
604  const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
605  const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight ) {
607  inputDataLeft,
608  inputDataRight );
609  }
610 
611  // ------------------------------------------------------------------------------------
612 
613  template<typename SpT>
614  template<typename outputFieldValueType, class ...outputFieldProperties,
615  typename inputDataValueType, class ...inputDataProperties,
616  typename inputFieldValueType, class ...inputFieldProperties>
617  void
619  vectorMultiplyDataField( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
620  const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
621  const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields ) {
622  const auto outRank = outputFields.rank();
623  switch (outRank) {
624  case 3:
625  case 4:
627  inputData,
628  inputFields );
629  break;
630  case 5:
632  inputData,
633  inputFields );
634  break;
635  default: {
636  INTREPID2_TEST_FOR_EXCEPTION( outRank < 3 && outRank > 5, std::runtime_error,
637  ">>> ERROR (FunctionSpaceTools::vectorMultiplyDataField): Output container must have rank 3, 4 or 5.");
638  }
639  }
640  }
641 
642  // ------------------------------------------------------------------------------------
643 
644  template<typename SpT>
645  template<typename outputDataValueType, class ...outputDataProperties,
646  typename inputDataLeftValueType, class ...inputDataLeftProperties,
647  typename inputDataRightValueType, class ...inputDataRightProperties>
648  void
650  vectorMultiplyDataData( Kokkos::DynRankView<outputDataValueType, outputDataProperties...> outputData,
651  const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
652  const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight ) {
653  const auto outRank = outputData.rank();
654  switch (outRank) {
655  case 2:
656  case 3:
658  inputDataLeft,
659  inputDataRight );
660  break;
661  case 4:
663  inputDataLeft,
664  inputDataRight );
665  break;
666  default: {
667  INTREPID2_TEST_FOR_EXCEPTION( outRank < 2 && outRank > 4, std::runtime_error,
668  ">>> ERROR (FunctionSpaceTools::vectorMultiplyDataData): Output container must have rank 2, 3 or 4.");
669  }
670  }
671  }
672 
673  // ------------------------------------------------------------------------------------
674 
675  template<typename SpT>
676  template<typename outputFieldValueType, class ...outputFieldProperties,
677  typename inputDataValueType, class ...inputDataProperties,
678  typename inputFieldValueType, class ...inputFieldProperties>
679  void
681  tensorMultiplyDataField( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
682  const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
683  const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields,
684  const char transpose ) {
685 
686  const auto outRank = outputFields.rank();
687  switch (outRank) {
688  case 4:
690  inputData,
691  inputFields,
692  transpose );
693  break;
694  case 5:
696  inputData,
697  inputFields,
698  transpose );
699  break;
700  default: {
701  INTREPID2_TEST_FOR_EXCEPTION( outRank < 4 && outRank > 5, std::runtime_error,
702  ">>> ERROR (FunctionSpaceTools::tensorMultiplyDataField): Output container must have rank 4 or 5.");
703  }
704  }
705  }
706 
707  // ------------------------------------------------------------------------------------
708 
709  template<typename SpT>
710  template<typename outputDataValueType, class ...outputDataProperties,
711  typename inputDataLeftValueType, class ...inputDataLeftProperties,
712  typename inputDataRightValueType, class ...inputDataRightProperties>
713  void
715  tensorMultiplyDataData( Kokkos::DynRankView<outputDataValueType, outputDataProperties...> outputData,
716  const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
717  const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight,
718  const char transpose ) {
719  const auto outRank = outputData.rank();
720  switch (outRank) {
721  case 3:
723  inputDataLeft,
724  inputDataRight,
725  transpose );
726  break;
727  case 4:
729  inputDataLeft,
730  inputDataRight,
731  transpose );
732  break;
733  default: {
734  INTREPID2_TEST_FOR_EXCEPTION( outRank < 4 && outRank > 5, std::runtime_error,
735  ">>> ERROR (FunctionSpaceTools::tensorMultiplyDataField): Output container must have rank 4 or 5.");
736  }
737  }
738  }
739 
740  // ------------------------------------------------------------------------------------
741 
742  namespace FunctorFunctionSpaceTools {
743 
747  template<typename inoutOperatorViewType,
748  typename fieldSignViewType>
750  inoutOperatorViewType _inoutOperator;
751  const fieldSignViewType _fieldSigns;
752 
753  KOKKOS_INLINE_FUNCTION
754  F_applyLeftFieldSigns( inoutOperatorViewType inoutOperator_,
755  fieldSignViewType fieldSigns_ )
756  : _inoutOperator(inoutOperator_), _fieldSigns(fieldSigns_) {}
757 
758  KOKKOS_INLINE_FUNCTION
759  void operator()(const ordinal_type cl) const {
760  const ordinal_type nlbf = _inoutOperator.extent(1);
761  const ordinal_type nrbf = _inoutOperator.extent(2);
762 
763  for (ordinal_type lbf=0;lbf<nlbf;++lbf)
764  for (ordinal_type rbf=0;rbf<nrbf;++rbf)
765  _inoutOperator(cl, lbf, rbf) *= _fieldSigns(cl, lbf);
766  }
767  };
768  }
769 
770  template<typename SpT>
771  template<typename inoutOperatorValueType, class ...inoutOperatorProperties,
772  typename fieldSignValueType, class ...fieldSignProperties>
773  void
775  applyLeftFieldSigns( Kokkos::DynRankView<inoutOperatorValueType,inoutOperatorProperties...> inoutOperator,
776  const Kokkos::DynRankView<fieldSignValueType, fieldSignProperties...> fieldSigns ) {
777 
778 #ifdef HAVE_INTREPID2_DEBUG
779  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.rank() != 3, std::invalid_argument,
780  ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): Input operator container must have rank 3.");
781  INTREPID2_TEST_FOR_EXCEPTION( fieldSigns.rank() != 2, std::invalid_argument,
782  ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): Input field signs container must have rank 2.");
783  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.extent(0) != fieldSigns.extent(0), std::invalid_argument,
784  ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): Zeroth dimensions (number of cells) of the operator and field signs containers must agree!");
785  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.extent(1) != fieldSigns.extent(1), std::invalid_argument,
786  ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): First dimensions (number of left fields) of the operator and field signs containers must agree!");
787 #endif
788 
789  typedef Kokkos::DynRankView<inoutOperatorValueType,inoutOperatorProperties...> inoutOperatorViewType;
790  typedef Kokkos::DynRankView<fieldSignValueType, fieldSignProperties...> fieldSignViewType;
792  <inoutOperatorViewType,fieldSignViewType> FunctorType;
793  typedef typename ExecSpace<typename inoutOperatorViewType::execution_space,SpT>::ExecSpaceType ExecSpaceType;
794 
795  const ordinal_type C = inoutOperator.extent(0);
796  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, C);
797  Kokkos::parallel_for( policy, FunctorType(inoutOperator, fieldSigns) );
798  }
799 
800  // ------------------------------------------------------------------------------------
801 
802  namespace FunctorFunctionSpaceTools {
806  template<typename inoutOperatorViewType,
807  typename fieldSignViewType>
809  inoutOperatorViewType _inoutOperator;
810  const fieldSignViewType _fieldSigns;
811 
812  KOKKOS_INLINE_FUNCTION
813  F_applyRightFieldSigns( inoutOperatorViewType inoutOperator_,
814  fieldSignViewType fieldSigns_ )
815  : _inoutOperator(inoutOperator_), _fieldSigns(fieldSigns_) {}
816 
817  KOKKOS_INLINE_FUNCTION
818  void operator()(const ordinal_type cl) const {
819  const ordinal_type nlbf = _inoutOperator.extent(1);
820  const ordinal_type nrbf = _inoutOperator.extent(2);
821 
822  for (ordinal_type lbf=0;lbf<nlbf;++lbf)
823  for (ordinal_type rbf=0;rbf<nrbf;++rbf)
824  _inoutOperator(cl, lbf, rbf) *= _fieldSigns(cl, rbf);
825  }
826  };
827  }
828 
829  template<typename SpT>
830  template<typename inoutOperatorValueType, class ...inoutOperatorProperties,
831  typename fieldSignValueType, class ...fieldSignProperties>
832  void
834  applyRightFieldSigns( Kokkos::DynRankView<inoutOperatorValueType,inoutOperatorProperties...> inoutOperator,
835  const Kokkos::DynRankView<fieldSignValueType, fieldSignProperties...> fieldSigns ) {
836 
837 #ifdef HAVE_INTREPID2_DEBUG
838  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.rank() != 3, std::invalid_argument,
839  ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Input operator container must have rank 3.");
840  INTREPID2_TEST_FOR_EXCEPTION( fieldSigns.rank() != 2, std::invalid_argument,
841  ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Input field signs container must have rank 2.");
842  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.extent(0) != fieldSigns.extent(0), std::invalid_argument,
843  ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Zeroth dimensions (number of cells) of the operator and field signs containers must agree!");
844  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.extent(2) != fieldSigns.extent(1), std::invalid_argument,
845  ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Second dimension of the operator container and first dimension of the field signs container (number of right fields) must agree!");
846 #endif
847 
848  typedef Kokkos::DynRankView<inoutOperatorValueType,inoutOperatorProperties...> inoutOperatorViewType;
849  typedef Kokkos::DynRankView<fieldSignValueType, fieldSignProperties...> fieldSignViewType;
851  <inoutOperatorViewType,fieldSignViewType> FunctorType;
852  typedef typename ExecSpace<typename inoutOperatorViewType::execution_space,SpT>::ExecSpaceType ExecSpaceType;
853 
854  const ordinal_type C = inoutOperator.extent(0);
855  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, C);
856  Kokkos::parallel_for( policy, FunctorType(inoutOperator, fieldSigns) );
857  }
858 
859  // ------------------------------------------------------------------------------------
860 
861  namespace FunctorFunctionSpaceTools {
865  template<typename inoutFunctionViewType,
866  typename fieldSignViewType>
868  inoutFunctionViewType _inoutFunction;
869  const fieldSignViewType _fieldSigns;
870 
871  KOKKOS_INLINE_FUNCTION
872  F_applyFieldSigns( inoutFunctionViewType inoutFunction_,
873  fieldSignViewType fieldSigns_)
874  : _inoutFunction(inoutFunction_), _fieldSigns(fieldSigns_) {}
875 
876  KOKKOS_INLINE_FUNCTION
877  void operator()(const ordinal_type cl) const {
878  const ordinal_type nbfs = _inoutFunction.extent(1);
879  const ordinal_type npts = _inoutFunction.extent(2);
880  const ordinal_type iend = _inoutFunction.extent(3);
881  const ordinal_type jend = _inoutFunction.extent(4);
882 
883  for (ordinal_type bf=0;bf<nbfs;++bf)
884  for (ordinal_type pt=0;pt<npts;++pt)
885  for (ordinal_type i=0;i<iend;++i)
886  for (ordinal_type j=0;j<jend;++j)
887  _inoutFunction(cl, bf, pt, i, j) *= _fieldSigns(cl, bf);
888  }
889  };
890  }
891 
892  template<typename SpT>
893  template<typename inoutFunctionValueType, class ...inoutFunctionProperties,
894  typename fieldSignValueType, class ...fieldSignProperties>
895  void
897  applyFieldSigns( Kokkos::DynRankView<inoutFunctionValueType,inoutFunctionProperties...> inoutFunction,
898  const Kokkos::DynRankView<fieldSignValueType, fieldSignProperties...> fieldSigns ) {
899 
900 #ifdef HAVE_INTREPID2_DEBUG
901  INTREPID2_TEST_FOR_EXCEPTION( inoutFunction.rank() < 2 || inoutFunction.rank() > 5, std::invalid_argument,
902  ">>> ERROR (FunctionSpaceTools::applyFieldSigns): Input function container must have rank 2, 3, 4, or 5.");
903  INTREPID2_TEST_FOR_EXCEPTION( fieldSigns.rank() != 2, std::invalid_argument,
904  ">>> ERROR (FunctionSpaceTools::applyFieldSigns): Input field signs container must have rank 2.");
905  INTREPID2_TEST_FOR_EXCEPTION( inoutFunction.extent(0) != fieldSigns.extent(0), std::invalid_argument,
906  ">>> ERROR (FunctionSpaceTools::applyFieldSigns): Zeroth dimensions (number of integration domains) of the function and field signs containers must agree!");
907  INTREPID2_TEST_FOR_EXCEPTION( inoutFunction.extent(1) != fieldSigns.extent(1), std::invalid_argument,
908  ">>> ERROR (FunctionSpaceTools::applyFieldSigns): First dimensions (number of fields) of the function and field signs containers must agree!");
909 
910 #endif
911 
912  typedef Kokkos::DynRankView<inoutFunctionValueType,inoutFunctionProperties...> inoutFunctionViewType;
913  typedef Kokkos::DynRankView<fieldSignValueType, fieldSignProperties...> fieldSignViewType;
915  <inoutFunctionViewType,fieldSignViewType> FunctorType;
916  typedef typename ExecSpace<typename inoutFunctionViewType::execution_space,SpT>::ExecSpaceType ExecSpaceType;
917 
918  const ordinal_type C = inoutFunction.extent(0);
919  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, C);
920  Kokkos::parallel_for( policy, FunctorType(inoutFunction, fieldSigns) );
921  }
922 
923  // ------------------------------------------------------------------------------------
924 
925  namespace FunctorFunctionSpaceTools {
926 
930  template<typename outputPointViewType,
931  typename inputCoeffViewType,
932  typename inputFieldViewType>
933  struct F_evaluate {
934  outputPointViewType _outputPointVals;
935  const inputCoeffViewType _inputCoeffs;
936  const inputFieldViewType _inputFields;
937 
938  KOKKOS_INLINE_FUNCTION
939  F_evaluate( outputPointViewType outputPointVals_,
940  inputCoeffViewType inputCoeffs_,
941  inputFieldViewType inputFields_ )
942  : _outputPointVals(outputPointVals_), _inputCoeffs(inputCoeffs_), _inputFields(inputFields_) {}
943 
944  KOKKOS_INLINE_FUNCTION
945  void operator()(const ordinal_type cl) const {
946  const ordinal_type nbfs = _inputFields.extent(1);
947  const ordinal_type npts = _inputFields.extent(2);
948 
949  const ordinal_type iend = _inputFields.extent(3);
950  const ordinal_type jend = _inputFields.extent(4);
951 
952  for (ordinal_type bf=0;bf<nbfs;++bf)
953  for (ordinal_type pt=0;pt<npts;++pt)
954  for (ordinal_type i=0;i<iend;++i)
955  for (ordinal_type j=0;j<jend;++j)
956  _outputPointVals(cl, pt, i, j) += _inputCoeffs(cl, bf) * _inputFields(cl, bf, pt, i, j);
957  }
958  };
959  }
960 
961  template<typename SpT>
962  template<typename outputPointValueType, class ...outputPointProperties,
963  typename inputCoeffValueType, class ...inputCoeffProperties,
964  typename inputFieldValueType, class ...inputFieldProperties>
965  void
967  evaluate( Kokkos::DynRankView<outputPointValueType,outputPointProperties...> outputPointVals,
968  const Kokkos::DynRankView<inputCoeffValueType, inputCoeffProperties...> inputCoeffs,
969  const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields ) {
970 
971 #ifdef HAVE_INTREPID2_DEBUG
972  INTREPID2_TEST_FOR_EXCEPTION( inputFields.rank() < 3 || inputFields.rank() > 5, std::invalid_argument,
973  ">>> ERROR (FunctionSpaceTools::evaluate): Input fields container must have rank 3, 4, or 5.");
974  INTREPID2_TEST_FOR_EXCEPTION( inputCoeffs.rank() != 2, std::invalid_argument,
975  ">>> ERROR (FunctionSpaceTools::evaluate): Input coefficient container must have rank 2.");
976  INTREPID2_TEST_FOR_EXCEPTION( outputPointVals.rank() != (inputFields.rank()-1), std::invalid_argument,
977  ">>> ERROR (FunctionSpaceTools::evaluate): Output values container must have rank one less than the rank of the input fields container.");
978  INTREPID2_TEST_FOR_EXCEPTION( inputCoeffs.extent(0) != inputFields.extent(0), std::invalid_argument,
979  ">>> ERROR (FunctionSpaceTools::evaluate): Zeroth dimensions (number of cells) of the coefficient and fields input containers must agree!");
980  INTREPID2_TEST_FOR_EXCEPTION( inputCoeffs.extent(1) != inputFields.extent(1), std::invalid_argument,
981  ">>> ERROR (FunctionSpaceTools::evaluate): First dimensions (number of fields) of the coefficient and fields input containers must agree!");
982  INTREPID2_TEST_FOR_EXCEPTION( outputPointVals.extent(0) != inputFields.extent(0), std::invalid_argument,
983  ">>> ERROR (FunctionSpaceTools::evaluate): Zeroth dimensions (number of cells) of the input fields container and the output values container must agree!");
984  for (size_type i=1;i<outputPointVals.rank();++i)
985  INTREPID2_TEST_FOR_EXCEPTION( outputPointVals.extent(i) != inputFields.extent(i+1), std::invalid_argument,
986  ">>> ERROR (FunctionSpaceTools::evaluate): outputPointVals dimension(i) does not match to inputFields dimension(i+1).");
987 #endif
988 
989  typedef Kokkos::DynRankView<outputPointValueType,outputPointProperties...> outputPointValViewType;
990  typedef Kokkos::DynRankView<inputCoeffValueType, inputCoeffProperties...> inputCoeffViewType;
991  typedef Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFieldViewType;
993  <outputPointValViewType,inputCoeffViewType,inputFieldViewType> FunctorType;
994  typedef typename ExecSpace<typename inputCoeffViewType::execution_space,SpT>::ExecSpaceType ExecSpaceType;
995 
996  const ordinal_type C = inputFields.extent(0);
997  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, C);
998  Kokkos::parallel_for( policy, FunctorType(outputPointVals, inputCoeffs, inputFields) );
999  }
1000 
1001  // ------------------------------------------------------------------------------------
1002 
1003 } // end namespace Intrepid2
1004 
1005 #endif
Functor for calculation HGRADtransformGRAD, see Intrepid2::FunctionSpaceTools for more...
Functor for applyLeftFieldSigns, see Intrepid2::FunctionSpaceTools for more.
static void contractDataFieldVector(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields, const bool sumInto=false)
Contracts the "point" and "space" dimensions P and D of a rank-4 container and a rank-3 container wit...
static void scalarMultiplyDataData(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight, const bool reciprocal=false)
There are two use cases: (1) multiplies a rank-2, 3, or 4 container inputDataRight with dimensions (C...
static void contractDataDataVector(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight, const bool sumInto=false)
Contracts the "point" and "space" dimensions P and D of rank-3 containers with dimensions (C...
static void HGRADtransformGRAD(Kokkos::DynRankView< outputValValueType, outputValProperties... > outputVals, const Kokkos::DynRankView< jacobianInverseValueType, jacobianInverseProperties... > jacobianInverse, const Kokkos::DynRankView< inputValValueType, inputValProperties... > inputVals)
Transformation of a gradient field in the H-grad space, defined at points on a reference cell...
static void matmatProductDataData(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight, const char transpose='N')
There are two use cases: (1) matrix-matrix product of a rank-4 container inputDataRight with dimensio...
static void HDIVtransformVALUE(Kokkos::DynRankView< outputValValueType, outputValProperties... > outputVals, const Kokkos::DynRankView< jacobianValueType, jacobianProperties... > jacobian, const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties... > jacobianDet, const Kokkos::DynRankView< inputValValueType, inputValProperties... > inputVals)
Transformation of a (vector) value field in the H-div space, defined at points on a reference cell...
static void dotMultiplyDataData(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight)
There are two use cases: (1) dot product of a rank-2, 3 or 4 container inputDataRight with dimensions...
static void HCURLtransformCURL(Kokkos::DynRankView< outputValValueType, outputValProperties... > outputVals, const Kokkos::DynRankView< jacobianValueType, jacobianProperties... > jacobian, const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties... > jacobianDet, const Kokkos::DynRankView< inputValValueType, inputValProperties... > inputVals)
Transformation of a curl field in the H-curl space, defined at points on a reference cell...
static void scalarMultiplyDataField(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields, const bool reciprocal=false)
There are two use cases: (1) multiplies a rank-3, 4, or 5 container inputFields with dimensions (C...
static void outerProductDataData(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValuetype, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight)
There are two use cases: (1) outer product of a rank-3 container inputDataRight with dimensions (C...
static void tensorMultiplyDataData(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight, const char transpose='N')
Matrix-vector or matrix-matrix product of data and data; please read the description below...
static void vectorNorm(Kokkos::DynRankView< normArrayValueType, normArrayProperties... > normArray, const Kokkos::DynRankView< inVecValueType, inVecProperties... > inVecs, const ENorm normType)
Computes norms (1, 2, infinity) of vectors stored in a array of total rank 2 (array of vectors)...
static void dotMultiplyDataField(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields)
Dot product of data and fields; please read the description below.
static void scalarMultiplyDataField(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataPropertes... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields, const bool reciprocal=false)
Scalar multiplication of data and fields; please read the description below.
Functor for applyRightFieldSigns, see Intrepid2::FunctionSpaceTools for more.
static void contractFieldFieldTensor(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< leftFieldValueType, leftFieldProperties... > leftFields, const Kokkos::DynRankView< rightFieldValueType, rightFieldProperties... > rightFields, const bool sumInto=false)
Contracts the "point" and "space" dimensions P, D1, and D2 of two rank-5 containers with dimensions (...
static void scalarMultiplyDataData(Kokkos::DynRankView< outputDataValuetype, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight, const bool reciprocal=false)
Scalar multiplication of data and data; please read the description below.
static void getPhysicalFaceNormals(Kokkos::DynRankView< faceNormalValueType, faceNormalProperties... > faceNormals, const Kokkos::DynRankView< worksetJacobianValueType, worksetJacobianProperties... > worksetJacobians, const ordinal_type worksetFaceOrd, const shards::CellTopology parentCell)
Computes non-normalized normal vectors to physical faces in a face workset ; (see Subcell worksets fo...
static void applyRightFieldSigns(Kokkos::DynRankView< inoutOperatorValueType, inoutOperatorProperties... > inoutOperator, const Kokkos::DynRankView< fieldSignValueType, fieldSignProperties... > fieldSigns)
Applies right (column) signs, stored in the user-provided container fieldSigns and indexed by (C...
static void contractDataDataTensor(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight, const bool sumInto=false)
Contracts the "point" and "space" dimensions P, D1 and D2 of rank-4 containers with dimensions (C...
static void HVOLtransformVALUE(Kokkos::DynRankView< outputValValueType, outputValProperties... > outputVals, const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties... > jacobianDet, const Kokkos::DynRankView< inputValValueType, inputValProperties... > inputVals)
Transformation of a (scalar) value field in the H-vol space, defined at points on a reference cell...
static void matmatProductDataField(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields, const char transpose='N')
There are two use cases: (1) matrix-matrix product of a rank-5 container inputFields with dimensions ...
static void contractFieldFieldVector(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< leftFieldValueType, leftFieldProperties... > leftFields, const Kokkos::DynRankView< rightFieldValueType, rightFieldProperties... > rightFields, const bool sumInto=false)
Contracts the "point" and "space" dimensions P and D1 of two rank-4 containers with dimensions (C...
static void vectorMultiplyDataField(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields)
Cross or outer product of data and fields; please read the description below.
static bool computeCellMeasure(Kokkos::DynRankView< outputValValueType, outputValProperties... > outputVals, const Kokkos::DynRankView< inputDetValueType, inputDetPropertes... > inputDet, const Kokkos::DynRankView< inputWeightValueType, inputWeightPropertes... > inputWeights)
Returns the weighted integration measures outputVals with dimensions (C,P) used for the computation o...
static void getPhysicalEdgeTangents(Kokkos::DynRankView< edgeTangentValueType, edgeTangentProperties... > edgeTangents, const Kokkos::DynRankView< worksetJacobianValueType, worksetJacobianProperties... > worksetJacobians, const ordinal_type worksetEdgeOrd, const shards::CellTopology parentCell)
Computes non-normalized tangent vectors to physical edges in an edge workset ; (see Subcell worksets ...
static void vectorMultiplyDataData(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight)
Cross or outer product of data and data; please read the description below.
static void HCURLtransformVALUE(Kokkos::DynRankView< outputValValueType, outputValProperties... > outputVals, const Kokkos::DynRankView< jacobianInverseValueType, jacobianInverseProperties... > jacobianInverse, const Kokkos::DynRankView< inputValValueType, inputValProperties... > inputVals)
Transformation of a (vector) value field in the H-curl space, defined at points on a reference cell...
static void tensorMultiplyDataField(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields, const char transpose='N')
Matrix-vector or matrix-matrix product of data and fields; please read the description below...
static void crossProductDataField(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields)
There are two use cases: (1) cross product of a rank-4 container inputFields with dimensions (C...
static void matvecProductDataData(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight, const char transpose='N')
There are two use cases: (1) matrix-vector product of a rank-3 container inputDataRight with dimensio...
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 crossProductDataData(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight)
There are two use cases: (1) cross product of a rank-3 container inputDataRight with dimensions (C...
static void applyFieldSigns(Kokkos::DynRankView< inoutFunctionValueType, inoutFunctionProperties... > inoutFunction, const Kokkos::DynRankView< fieldSignValueType, fieldSignProperties... > fieldSigns)
Applies field signs, stored in the user-provided container fieldSigns and indexed by (C...
static void computeEdgeMeasure(Kokkos::DynRankView< outputValValueType, outputValProperties... > outputVals, const Kokkos::DynRankView< inputJacValueType, inputJacProperties... > inputJac, const Kokkos::DynRankView< inputWeightValueType, inputWeightProperties... > inputWeights, const int whichEdge, const shards::CellTopology parentCell, const Kokkos::DynRankView< scratchValueType, scratchProperties... > scratch)
Returns the weighted integration measures outVals with dimensions (C,P) used for the computation of e...
static void HGRADtransformVALUE(Kokkos::DynRankView< outputValueType, outputProperties... > output, const Kokkos::DynRankView< inputValueType, inputProperties... > input)
Transformation of a (scalar) value field in the H-grad space, defined at points on a reference cell...
static void contractDataDataScalar(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight, const bool sumInto=false)
Contracts the "point" dimensions P of rank-2 containers with dimensions (C,P), and returns the result...
static void applyLeftFieldSigns(Kokkos::DynRankView< inoutOperatorValueType, inoutOperatorProperties... > inoutOperator, const Kokkos::DynRankView< fieldSignValueType, fieldSignProperties... > fieldSigns)
Applies left (row) signs, stored in the user-provided container fieldSigns and indexed by (C...
static void contractFieldFieldScalar(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< leftFieldValueType, leftFieldProperties... > leftFields, const Kokkos::DynRankView< rightFieldValueType, rightFieldProperties... > rightFields, const bool sumInto=false)
Contracts the "point" dimension P of two rank-3 containers with dimensions (C,L,P) and (C...
Functor to evaluate functions, see Intrepid2::FunctionSpaceTools for more.
static void computeFaceMeasure(Kokkos::DynRankView< outputValValueType, outputValProperties... > outputVals, const Kokkos::DynRankView< inputJacValueType, inputJacProperties... > inputJac, const Kokkos::DynRankView< inputWeightValueType, inputWeightPropertes... > inputWeights, const int whichFace, const shards::CellTopology parentCell, const Kokkos::DynRankView< scratchValueType, scratchProperties... > scratch)
Returns the weighted integration measures outputVals with dimensions (C,P) used for the computation o...
static void evaluate(Kokkos::DynRankView< outputPointValueType, outputPointProperties... > outputPointVals, const Kokkos::DynRankView< inputCoeffValueType, inputCoeffProperties... > inputCoeffs, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields)
Computes point values outPointVals of a discrete function specified by the basis inFields and coeffic...
static void cloneFields(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields)
Replicates a rank-2, 3, or 4 container with dimensions (F,P), (F,P,D1) or (F,P,D1,D2), representing the values of a scalar, vector or a tensor field, into an output value container of size (C,F,P), (C,F,P,D1) or (C,F,P,D1,D2).
static void HDIVtransformDIV(Kokkos::DynRankView< outputValValueType, outputValProperties... > outputVals, const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties... > jacobianDet, const Kokkos::DynRankView< inputValValueType, inputValProperties... > inputVals)
Transformation of a divergence field in the H-div space, defined at points on a reference cell...
Functor for calculation of cell measure, see Intrepid2::FunctionSpaceTools for more.
static void outerProductDataField(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields)
There are two use cases: (1) outer product of a rank-4 container inputFields with dimensions (C...
static void dotMultiplyDataData(Kokkos::DynRankView< outputDataValueType, outputDataProperties... > outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties... > inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties... > inputDataRight)
Dot product of data and data; please read the description below.
static void multiplyMeasure(Kokkos::DynRankView< outputValValueType, outputValProperties... > outputVals, const Kokkos::DynRankView< inputMeasureValueType, inputMeasureProperties... > inputMeasure, const Kokkos::DynRankView< inputValValueType, inputValProperteis... > inputVals)
Multiplies fields inputVals by weighted measures inputMeasure and returns the field array outputVals;...
static void contractDataFieldScalar(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields, const bool sumInto=false)
Contracts the "point" dimensions P of a rank-3 containers and a rank-2 container with dimensions (C...
static void contractDataFieldTensor(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields, const bool sumInto=false)
Contracts the "point" and "space" dimensions P, D1 and D2 of a rank-5 container and a rank-4 containe...
Functor for applyFieldSigns, see Intrepid2::FunctionSpaceTools for more.
static void dotMultiplyDataField(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields)
There are two use cases: (1) dot product of a rank-3, 4 or 5 container inputFields with dimensions (C...
static void matvecProductDataField(Kokkos::DynRankView< outputFieldValueType, outputFieldProperties... > outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties... > inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties... > inputFields, const char transpose='N')
There are two use cases: (1) matrix-vector product of a rank-4 container inputFields with dimensions ...