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 
54 
55 #include "Teuchos_TimeMonitor.hpp"
56 
57 namespace Intrepid2 {
58  // ------------------------------------------------------------------------------------
59  template<typename DeviceType>
60  template<typename outputValueType, class ...outputProperties,
61  typename inputValueType, class ...inputProperties>
62  void
64  HGRADtransformVALUE( Kokkos::DynRankView<outputValueType,outputProperties...> output,
65  const Kokkos::DynRankView<inputValueType, inputProperties...> input ) {
66  if(output.rank() == input.rank()) {
67 #ifdef HAVE_INTREPID2_DEBUG
68  {
69  for (size_type i=0;i< input.rank();++i) {
70  INTREPID2_TEST_FOR_EXCEPTION( (input.extent(i) != output.extent(i)), std::invalid_argument,
71  ">>> ERROR (FunctionSpaceTools::HGRADtransformVALUE): Dimensions of input and output fields containers do not match.");
72  }
73  }
74 #endif
75  RealSpaceTools<DeviceType>::clone(output, input);
76  }
77  else
79  }
80 
81  // ------------------------------------------------------------------------------------
82 
83  namespace FunctorFunctionSpaceTools {
87  template <typename OutputViewType,
88  typename jacInverseViewType,
89  typename inputViewType,
90  ordinal_type spaceDim>
92  OutputViewType _output;
93  const jacInverseViewType _jacInverse;
94  const inputViewType _input;
95 
96  // output CPDD, left CPDD or PDD, right FPD
97  KOKKOS_INLINE_FUNCTION
98  F_HGRADtransformGRAD(OutputViewType output_,
99  jacInverseViewType jacInverse_,
100  inputViewType input_)
101  : _output(output_),
102  _jacInverse(jacInverse_),
103  _input(input_) {}
104 
105  KOKKOS_INLINE_FUNCTION
106  void operator()(const ordinal_type cl,
107  const ordinal_type bf,
108  const ordinal_type pt) const {
109  /* */ auto y = Kokkos::subview(_output, cl, bf, pt, Kokkos::ALL());
110  const auto A = Kokkos::subview(_jacInverse, cl, pt, Kokkos::ALL(), Kokkos::ALL());
111  const auto x = Kokkos::subview(_input, bf, pt, Kokkos::ALL());
112 
113  if (spaceDim == 2) {
114  Kernels::Serial::matvec_trans_product_d2( y, A, x );
115  } else {
116  Kernels::Serial::matvec_trans_product_d3( y, A, x );
117  }
118  }
119  };
120  }
121 
122  template<typename DeviceType>
123  template<typename outputValValueType, class ...outputValProperties,
124  typename jacobianInverseValueType, class ...jacobianInverseProperties,
125  typename inputValValueType, class ...inputValProperties>
126  void
128  HGRADtransformGRAD( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
129  const Kokkos::DynRankView<jacobianInverseValueType,jacobianInverseProperties...> jacobianInverse,
130  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
131  return HCURLtransformVALUE(outputVals, jacobianInverse, inputVals);
132 
133  // this modification is for 2d and 3d (not 1d)
134  // this is an attempt to measure the overhead of subview of dynrankview.
135 
136  // typedef Kokkos::DynRankView<outputValValueType,outputValProperties...> OutputViewType;
137  // typedef const Kokkos::DynRankView<jacobianInverseValueType,jacobianInverseProperties...> jacInverseViewType;
138  // typedef const Kokkos::DynRankView<inputValValueType,inputValProperties...> inputViewType;
139 
140  // const ordinal_type
141  // C = outputVals.extent(0),
142  // F = outputVals.extent(1),
143  // P = outputVals.extent(2);
144 
145  // using range_policy_type = Kokkos::Experimental::MDRangePolicy
146  // < DeviceType, Kokkos::Experimental::Rank<3>, Kokkos::IndexType<ordinal_type> >;
147  // range_policy_type policy( { 0, 0, 0 },
148  // { C, F, P } );
149 
150  // const ordinal_type spaceDim = inputVals.extent(2);
151  // switch (spaceDim) {
152  // case 2: {
153  // typedef FunctorFunctionSpaceTools::F_HGRADtransformGRAD<OutputViewType, jacInverseViewType, inputViewType, 2> FunctorType;
154  // Kokkos::parallel_for( policy, FunctorType(outputVals, jacobianInverse, inputVals) );
155  // break;
156  // }
157  // case 3: {
158  // typedef FunctorFunctionSpaceTools::F_HGRADtransformGRAD<OutputViewType, jacInverseViewType, inputViewType, 3> FunctorType;
159  // Kokkos::parallel_for( policy, FunctorType(outputVals, jacobianInverse, inputVals) );
160  // break;
161  // }
162  // default: {
163  // INTREPID2_TEST_FOR_EXCEPTION( true, std::invalid_argument,
164  // ">>> ERROR (FunctionSpaceTools::HGRADtransformGRAD): spaceDim is not 2 or 3.");
165  // break;
166  // }
167  // }
168  }
169 
170  // ------------------------------------------------------------------------------------
171 
172  template<typename DeviceType>
173  template<typename outputValValueType, class ...outputValProperties,
174  typename jacobianInverseValueType, class ...jacobianInverseProperties,
175  typename inputValValueType, class ...inputValProperties>
176  void
178  HCURLtransformVALUE( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
179  const Kokkos::DynRankView<jacobianInverseValueType,jacobianInverseProperties...> jacobianInverse,
180  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
181  ArrayTools<DeviceType>::matvecProductDataField(outputVals, jacobianInverse, inputVals, 'T');
182  }
183 
184  // ------------------------------------------------------------------------------------
185 
186  template<typename DeviceType>
187  template<typename outputValValueType, class ...outputValProperties,
188  typename jacobianValueType, class ...jacobianProperties,
189  typename jacobianDetValueType, class ...jacobianDetProperties,
190  typename inputValValueType, class ...inputValProperties>
191  void
193  HCURLtransformCURL( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
194  const Kokkos::DynRankView<jacobianValueType, jacobianProperties...> jacobian,
195  const Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
196  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
197  if(jacobian.data()==NULL || jacobian.extent(2)==2) //2D case
198  return HVOLtransformVALUE(outputVals, jacobianDet, inputVals);
199  else
200  return HDIVtransformVALUE(outputVals, jacobian, jacobianDet, inputVals);
201  }
202 
203  // ------------------------------------------------------------------------------------
204 
205  template<typename DeviceType>
206  template<typename outputValValueType, class ...outputValProperties,
207  typename jacobianDetValueType, class ...jacobianDetProperties,
208  typename inputValValueType, class ...inputValProperties>
209  void
211  HCURLtransformCURL( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
212  const Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
213  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
214 #ifdef HAVE_INTREPID2_DEBUG
215  {
216  INTREPID2_TEST_FOR_EXCEPTION( outputVals.rank() == 4, std::invalid_argument,
217  ">>> ERROR (FunctionSpaceTools::HCURLtransformCURL): Output rank must have rank 3.\n If these are 3D fields, then use the appropriate overload of this function.");
218  }
219 #endif
220  return HVOLtransformVALUE(outputVals, jacobianDet, inputVals);
221  }
222 
223  // ------------------------------------------------------------------------------------
224 
225  template<typename DeviceType>
226  template<typename outputValValueType, class ...outputValProperties,
227  typename jacobianValueType, class ...jacobianProperties,
228  typename jacobianDetValueType, class ...jacobianDetProperties,
229  typename inputValValueType, class ...inputValProperties>
230  void
232  HGRADtransformCURL( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
233  const Kokkos::DynRankView<jacobianValueType, jacobianProperties...> jacobian,
234  const Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
235  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
236 #ifdef HAVE_INTREPID2_DEBUG
237  {
238  INTREPID2_TEST_FOR_EXCEPTION( outputVals.extent(3)!=2, std::invalid_argument,
239  ">>> ERROR (FunctionSpaceTools::HGRADtransformCURL):\n output field is 3D by the function is meant for 2D fields");
240  }
241 #endif
242  return HDIVtransformVALUE(outputVals, jacobian, jacobianDet, inputVals);
243  }
244 
245  // ------------------------------------------------------------------------------------
246 
247  template<typename DeviceType>
248  template<typename outputValValueType, class ...outputValProperties,
249  typename jacobianValueType, class ...jacobianProperties,
250  typename jacobianDetValueType, class ...jacobianDetProperties,
251  typename inputValValueType, class ...inputValProperties>
252  void
254  HDIVtransformVALUE( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
255  const Kokkos::DynRankView<jacobianValueType, jacobianProperties...> jacobian,
256  const Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
257  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
258  ArrayTools<DeviceType>::matvecProductDataField(outputVals, jacobian, inputVals, 'N');
259  ArrayTools<DeviceType>::scalarMultiplyDataField(outputVals, jacobianDet, outputVals, true);
260  }
261 
262  // ------------------------------------------------------------------------------------
263 
264  template<typename DeviceType>
265  template<typename outputValValueType, class ...outputValProperties,
266  typename jacobianDetValueType, class ...jacobianDetProperties,
267  typename inputValValueType, class ...inputValProperties>
268  void
270  HDIVtransformDIV( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
271  const Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
272  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
273  return HVOLtransformVALUE(outputVals, jacobianDet, inputVals);
274  }
275 
276  // ------------------------------------------------------------------------------------
277 
278  template<typename DeviceType>
279  template<typename outputValValueType, class ...outputValProperties,
280  typename jacobianDetValueType, class ...jacobianDetProperties,
281  typename inputValValueType, class ...inputValProperties>
282  void
284  HVOLtransformVALUE( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
285  const Kokkos::DynRankView<jacobianDetValueType,jacobianDetProperties...> jacobianDet,
286  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
287  ArrayTools<DeviceType>::scalarMultiplyDataField(outputVals, jacobianDet, inputVals, true);
288  }
289 
290  // ------------------------------------------------------------------------------------
291 
292  template<typename DeviceType>
293  template<typename outputValueValueType, class ...outputValueProperties,
294  typename leftValueValueType, class ...leftValueProperties,
295  typename rightValueValueType, class ...rightValueProperties>
296  void
298  integrate( Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
299  const Kokkos::DynRankView<leftValueValueType, leftValueProperties...> leftValues,
300  const Kokkos::DynRankView<rightValueValueType, rightValueProperties...> rightValues,
301  const bool sumInto ) {
302 
303 #ifdef HAVE_INTREPID2_DEBUG
304  {
305  INTREPID2_TEST_FOR_EXCEPTION( leftValues.rank() < 2 ||
306  leftValues.rank() > 4, std::invalid_argument,
307  ">>> ERROR (FunctionSpaceTools::integrate): Left data must have rank 2, 3 or 4.");
308  INTREPID2_TEST_FOR_EXCEPTION( outputValues.rank() < 1 ||
309  outputValues.rank() > 3, std::invalid_argument,
310  ">>> ERROR (FunctionSpaceTools::integrate): Output values must have rank 1, 2 or 3.");
311  }
312 #endif
313 
314  const ordinal_type outRank = outputValues.rank();
315  const ordinal_type leftRank = leftValues.rank();
316  const ordinal_type mode = outRank*10 + leftRank;
317 
318  switch (mode) {
319  // DataData
320  case 12:
322  leftValues,
323  rightValues,
324  sumInto );
325  break;
326  case 13:
328  leftValues,
329  rightValues,
330  sumInto );
331  break;
332  case 14:
334  leftValues,
335  rightValues,
336  sumInto );
337  break;
338 
339  // DataField
340  case 22:
342  leftValues,
343  rightValues,
344  sumInto );
345  break;
346  case 23:
348  leftValues,
349  rightValues,
350  sumInto );
351  break;
352  case 24:
354  leftValues,
355  rightValues,
356  sumInto );
357  break;
358 
359  // FieldField
360  case 33:
362  leftValues,
363  rightValues,
364  sumInto );
365  break;
366  case 34:
368  leftValues,
369  rightValues,
370  sumInto );
371  break;
372  case 35:
374  leftValues,
375  rightValues,
376  sumInto );
377  break;
378  default: {
379  INTREPID2_TEST_FOR_EXCEPTION( outRank < 1 || outRank > 3, std::runtime_error,
380  ">>> ERROR (FunctionSpaceTools::integrate): outRank must be 1,2, or 3.");
381  INTREPID2_TEST_FOR_EXCEPTION( leftRank < 2 || leftRank > 4, std::runtime_error,
382  ">>> ERROR (FunctionSpaceTools::integrate): leftRank must be 1,2, 3 or 4.");
383  }
384  }
385  }
386 
387  // ------------------------------------------------------------------------------------
388  namespace FunctorFunctionSpaceTools {
392  template<typename outputValViewType,
393  typename inputDetViewType,
394  typename inputWeightViewType>
396  outputValViewType _outputVals;
397  const inputDetViewType _inputDet;
398  const inputWeightViewType _inputWeight;
399 
400  KOKKOS_INLINE_FUNCTION
401  F_computeCellMeasure( outputValViewType outputVals_,
402  inputDetViewType inputDet_,
403  inputWeightViewType inputWeight_)
404  : _outputVals(outputVals_),
405  _inputDet(inputDet_),
406  _inputWeight(inputWeight_) {}
407 
408  typedef ordinal_type value_type;
409 
410 // KOKKOS_INLINE_FUNCTION
411 // void init( value_type &dst ) const {
412 // dst = false;
413 // }
414 
415 // KOKKOS_INLINE_FUNCTION
416 // void join( volatile value_type &dst,
417 // const volatile value_type &src ) const {
418 // dst |= src;
419 // }
420 
421  KOKKOS_INLINE_FUNCTION
422  void operator()(const size_type cl, value_type &dst) const {
423  // negative jacobian check
424  const bool hasNegativeDet = (_inputDet(cl, 0) < 0.0);
425  dst |= hasNegativeDet;
426 
427  // make it absolute
428  const auto sign = (hasNegativeDet ? -1.0 : 1.0);
429  const ordinal_type pt_end = _outputVals.extent(1);
430  for (ordinal_type pt=0;pt<pt_end;++pt)
431  _outputVals(cl, pt) = sign*_inputDet(cl, pt)*_inputWeight(pt);
432  }
433  };
434  }
435 
436  template<typename DeviceType>
437  template<typename outputValValueType, class ...outputValProperties,
438  typename inputDetValueType, class ...inputDetProperties,
439  typename inputWeightValueType, class ...inputWeightProperties>
440  bool
442  computeCellMeasure( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
443  const Kokkos::DynRankView<inputDetValueType, inputDetProperties...> inputDet,
444  const Kokkos::DynRankView<inputWeightValueType,inputWeightProperties...> inputWeights ) {
445 #ifdef HAVE_INTREPID2_DEBUG
446  {
447  INTREPID2_TEST_FOR_EXCEPTION( inputDet.rank() != 2 ||
448  inputWeights.rank() != 1 ||
449  outputVals.rank() != 2, std::invalid_argument,
450  ">>> ERROR (FunctionSpaceTools::computeCellMeasure): Ranks are not compatible.");
451  INTREPID2_TEST_FOR_EXCEPTION( outputVals.extent(0) != inputDet.extent(0), std::invalid_argument,
452  ">>> ERROR (FunctionSpaceTools::computeCellMeasure): Cell dimension does not match.");
453  INTREPID2_TEST_FOR_EXCEPTION( inputDet.extent(1) != outputVals.extent(1) ||
454  inputWeights.extent(0) != outputVals.extent(1), std::invalid_argument,
455  ">>> ERROR (FunctionSpaceTools::computeCellMeasure): Point dimension does not match.");
456  }
457 #endif
458  constexpr bool are_accessible =
459  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
460  typename decltype(outputVals)::memory_space>::accessible &&
461  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
462  typename decltype(inputDet)::memory_space>::accessible &&
463  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
464  typename decltype(inputWeights)::memory_space>::accessible;
465  static_assert(are_accessible, "FunctionSpaceTools<DeviceType>::computeCellMeasure(..): input/output views' memory spaces are not compatible with DeviceType");
466 
468  <decltype(outputVals),decltype(inputDet),decltype(inputWeights)>;
469 
470  const ordinal_type C = inputDet.extent(0);
471  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, C);
472 
473  typename FunctorType::value_type hasNegativeDet = false;
474  Kokkos::parallel_reduce( policy, FunctorType(outputVals, inputDet, inputWeights), hasNegativeDet );
475 
476  return hasNegativeDet;
477  }
478 
479  // ------------------------------------------------------------------------------------
480 
481  template<typename DeviceType>
482  template<typename outputValValueType, class ...outputValProperties,
483  typename inputJacValueType, class ...inputJacProperties,
484  typename inputWeightValueType, class ...inputWeightProperties,
485  typename scratchValueType, class ...scratchProperties>
486  void
488  computeFaceMeasure( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
489  const Kokkos::DynRankView<inputJacValueType, inputJacProperties...> inputJac,
490  const Kokkos::DynRankView<inputWeightValueType,inputWeightProperties...> inputWeights,
491  const ordinal_type whichFace,
492  const shards::CellTopology parentCell,
493  const Kokkos::DynRankView<scratchValueType, scratchProperties...> scratch ) {
494 #ifdef HAVE_INTREPID2_DEBUG
495  INTREPID2_TEST_FOR_EXCEPTION( inputJac.rank() != 4, std::invalid_argument,
496  ">>> ERROR (FunctionSpaceTools::computeFaceMeasure): Input Jacobian container must have rank 4.");
497  INTREPID2_TEST_FOR_EXCEPTION( scratch.rank() != 1, std::invalid_argument,
498  ">>> ERROR (FunctionSpaceTools::computeFaceMeasure): Scratch view imust have rank 1.");
499  INTREPID2_TEST_FOR_EXCEPTION( scratch.span() < inputJac.span(), std::invalid_argument,
500  ">>> ERROR (FunctionSpaceTools::computeFaceMeasure): Scratch storage must be greater than or equal to inputJac's one.");
501 #endif
502 
503  // face normals (reshape scratch)
504  // Kokkos::DynRankView<scratchValueType,scratchProperties...> faceNormals(scratch.data(),
505  // inputJac.extent(0),
506  // inputJac.extent(1),
507  // inputJac.extent(2));
508  auto vcprop = Kokkos::common_view_alloc_prop(scratch);
509  //typedef Kokkos::DynRankView<scratchValueType, typename decltype(scratch)::memory_space> viewType;
510  typedef Kokkos::DynRankView<scratchValueType, DeviceType> viewType;
511  viewType faceNormals(Kokkos::view_wrap(scratch.data(), vcprop),
512  inputJac.extent(0),
513  inputJac.extent(1),
514  inputJac.extent(2));
515 
516  // compute normals
517  CellTools<DeviceType>::getPhysicalFaceNormals(faceNormals, inputJac, whichFace, parentCell);
518 
519  // compute lenghts of normals
520  RealSpaceTools<DeviceType>::vectorNorm(outputVals, faceNormals, NORM_TWO);
521 
522  // multiply with weights
523  ArrayTools<DeviceType>::scalarMultiplyDataData(outputVals, outputVals, inputWeights);
524  }
525 
526  // ------------------------------------------------------------------------------------
527 
528  template<typename DeviceType>
529  template<typename outputValValueType, class ...outputValProperties,
530  typename inputJacValueType, class ...inputJacProperties,
531  typename inputWeightValueType, class ...inputWeightProperties,
532  typename scratchValueType, class ...scratchProperties>
533  void
535  computeEdgeMeasure( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
536  const Kokkos::DynRankView<inputJacValueType, inputJacProperties...> inputJac,
537  const Kokkos::DynRankView<inputWeightValueType,inputWeightProperties...> inputWeights,
538  const ordinal_type whichEdge,
539  const shards::CellTopology parentCell,
540  const Kokkos::DynRankView<scratchValueType, scratchProperties...> scratch ) {
541 #ifdef HAVE_INTREPID2_DEBUG
542  INTREPID2_TEST_FOR_EXCEPTION( (inputJac.rank() != 4), std::invalid_argument,
543  ">>> ERROR (FunctionSpaceTools::computeEdgeMeasure): Input Jacobian container must have rank 4.");
544  INTREPID2_TEST_FOR_EXCEPTION( scratch.rank() != 1, std::invalid_argument,
545  ">>> ERROR (FunctionSpaceTools::computeEdgeMeasure): Scratch view must have a rank 1.");
546  INTREPID2_TEST_FOR_EXCEPTION( scratch.span() < inputJac.span(), std::invalid_argument,
547  ">>> ERROR (FunctionSpaceTools::computeEdgeMeasure): Scratch storage must be greater than or equal to inputJac'one.");
548 #endif
549 
550  // edge tangents (reshape scratch)
551  // Kokkos::DynRankView<scratchValueType,scratchProperties...> edgeTangents(scratch.data(),
552  // inputJac.extent(0),
553  // inputJac.extent(1),
554  // inputJac.extent(2));
555  auto vcprop = Kokkos::common_view_alloc_prop(scratch);
556  //typedef Kokkos::DynRankView<scratchValueType, typename decltype(scratch)::memory_space> viewType;
557  typedef Kokkos::DynRankView<scratchValueType, DeviceType> viewType;
558  viewType edgeTangents(Kokkos::view_wrap(scratch.data(), vcprop),
559  inputJac.extent(0),
560  inputJac.extent(1),
561  inputJac.extent(2));
562 
563  // compute normals
564  CellTools<DeviceType>::getPhysicalEdgeTangents(edgeTangents, inputJac, whichEdge, parentCell);
565 
566  // compute lenghts of tangents
567  RealSpaceTools<DeviceType>::vectorNorm(outputVals, edgeTangents, NORM_TWO);
568 
569  // multiply with weights
570  ArrayTools<DeviceType>::scalarMultiplyDataData(outputVals, outputVals, inputWeights);
571  }
572 
573  // ------------------------------------------------------------------------------------
574 
575  template<typename DeviceType>
576  template<typename outputValValueType, class ...outputValProperties,
577  typename inputMeasureValueType, class ...inputMeasureProperties,
578  typename inputValValueType, class ...inputValProperties>
579  void
581  multiplyMeasure( Kokkos::DynRankView<outputValValueType, outputValProperties...> outputVals,
582  const Kokkos::DynRankView<inputMeasureValueType,inputMeasureProperties...> inputMeasure,
583  const Kokkos::DynRankView<inputValValueType, inputValProperties...> inputVals ) {
584  scalarMultiplyDataField( outputVals,
585  inputMeasure,
586  inputVals );
587  }
588 
589  // ------------------------------------------------------------------------------------
590 
591  template<typename DeviceType>
592  template<typename outputFieldValueType, class ...outputFieldProperties,
593  typename inputDataValueType, class ...inputDataProperties,
594  typename inputFieldValueType, class ...inputFieldProperties>
595  void
597  scalarMultiplyDataField( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
598  const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
599  const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields,
600  const bool reciprocal ) {
602  inputData,
603  inputFields,
604  reciprocal );
605  }
606 
607  // ------------------------------------------------------------------------------------
608 
609  template<typename DeviceType>
610  template<typename outputDataValuetype, class ...outputDataProperties,
611  typename inputDataLeftValueType, class ...inputDataLeftProperties,
612  typename inputDataRightValueType, class ...inputDataRightProperties>
613  void
615  scalarMultiplyDataData( Kokkos::DynRankView<outputDataValuetype, outputDataProperties...> outputData,
616  const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
617  const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight,
618  const bool reciprocal ) {
620  inputDataLeft,
621  inputDataRight,
622  reciprocal );
623  }
624 
625  // ------------------------------------------------------------------------------------
626 
627  template<typename DeviceType>
628  template<typename outputFieldValueType, class ...outputFieldProperties,
629  typename inputDataValueType, class ...inputDataProperties,
630  typename inputFieldValueType, class ...inputFieldProperties>
631  void
633  dotMultiplyDataField( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
634  const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
635  const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields ) {
637  inputData,
638  inputFields );
639  }
640 
641  // ------------------------------------------------------------------------------------
642 
643  template<typename DeviceType>
644  template<typename outputDataValueType, class ...outputDataProperties,
645  typename inputDataLeftValueType, class ...inputDataLeftProperties,
646  typename inputDataRightValueType, class ...inputDataRightProperties>
647  void
649  dotMultiplyDataData( Kokkos::DynRankView<outputDataValueType, outputDataProperties...> outputData,
650  const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
651  const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight ) {
653  inputDataLeft,
654  inputDataRight );
655  }
656 
657  // ------------------------------------------------------------------------------------
658 
659  template<typename DeviceType>
660  template<typename outputFieldValueType, class ...outputFieldProperties,
661  typename inputDataValueType, class ...inputDataProperties,
662  typename inputFieldValueType, class ...inputFieldProperties>
663  void
665  vectorMultiplyDataField( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
666  const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
667  const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields ) {
668  const auto outRank = outputFields.rank();
669  switch (outRank) {
670  case 3:
671  case 4:
673  inputData,
674  inputFields );
675  break;
676  case 5:
678  inputData,
679  inputFields );
680  break;
681  default: {
682  INTREPID2_TEST_FOR_EXCEPTION( outRank < 3 && outRank > 5, std::runtime_error,
683  ">>> ERROR (FunctionSpaceTools::vectorMultiplyDataField): Output container must have rank 3, 4 or 5.");
684  }
685  }
686  }
687 
688  // ------------------------------------------------------------------------------------
689 
690  template<typename DeviceType>
691  template<typename outputDataValueType, class ...outputDataProperties,
692  typename inputDataLeftValueType, class ...inputDataLeftProperties,
693  typename inputDataRightValueType, class ...inputDataRightProperties>
694  void
696  vectorMultiplyDataData( Kokkos::DynRankView<outputDataValueType, outputDataProperties...> outputData,
697  const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
698  const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight ) {
699  const auto outRank = outputData.rank();
700  switch (outRank) {
701  case 2:
702  case 3:
704  inputDataLeft,
705  inputDataRight );
706  break;
707  case 4:
709  inputDataLeft,
710  inputDataRight );
711  break;
712  default: {
713  INTREPID2_TEST_FOR_EXCEPTION( outRank < 2 && outRank > 4, std::runtime_error,
714  ">>> ERROR (FunctionSpaceTools::vectorMultiplyDataData): Output container must have rank 2, 3 or 4.");
715  }
716  }
717  }
718 
719  // ------------------------------------------------------------------------------------
720 
721  template<typename DeviceType>
722  template<typename outputFieldValueType, class ...outputFieldProperties,
723  typename inputDataValueType, class ...inputDataProperties,
724  typename inputFieldValueType, class ...inputFieldProperties>
725  void
727  tensorMultiplyDataField( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
728  const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
729  const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields,
730  const char transpose ) {
731 
732  const auto outRank = outputFields.rank();
733  switch (outRank) {
734  case 4:
736  inputData,
737  inputFields,
738  transpose );
739  break;
740  case 5:
742  inputData,
743  inputFields,
744  transpose );
745  break;
746  default: {
747  INTREPID2_TEST_FOR_EXCEPTION( outRank < 4 && outRank > 5, std::runtime_error,
748  ">>> ERROR (FunctionSpaceTools::tensorMultiplyDataField): Output container must have rank 4 or 5.");
749  }
750  }
751  }
752 
753  // ------------------------------------------------------------------------------------
754 
755  template<typename DeviceType>
756  template<typename outputDataValueType, class ...outputDataProperties,
757  typename inputDataLeftValueType, class ...inputDataLeftProperties,
758  typename inputDataRightValueType, class ...inputDataRightProperties>
759  void
761  tensorMultiplyDataData( Kokkos::DynRankView<outputDataValueType, outputDataProperties...> outputData,
762  const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
763  const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight,
764  const char transpose ) {
765  const auto outRank = outputData.rank();
766  switch (outRank) {
767  case 3:
769  inputDataLeft,
770  inputDataRight,
771  transpose );
772  break;
773  case 4:
775  inputDataLeft,
776  inputDataRight,
777  transpose );
778  break;
779  default: {
780  INTREPID2_TEST_FOR_EXCEPTION( outRank < 4 && outRank > 5, std::runtime_error,
781  ">>> ERROR (FunctionSpaceTools::tensorMultiplyDataField): Output container must have rank 4 or 5.");
782  }
783  }
784  }
785 
786  // ------------------------------------------------------------------------------------
787 
788  namespace FunctorFunctionSpaceTools {
789 
793  template<typename inoutOperatorViewType,
794  typename fieldSignViewType>
796  inoutOperatorViewType _inoutOperator;
797  const fieldSignViewType _fieldSigns;
798 
799  KOKKOS_INLINE_FUNCTION
800  F_applyLeftFieldSigns( inoutOperatorViewType inoutOperator_,
801  fieldSignViewType fieldSigns_ )
802  : _inoutOperator(inoutOperator_), _fieldSigns(fieldSigns_) {}
803 
804  KOKKOS_INLINE_FUNCTION
805  void operator()(const ordinal_type cl) const {
806  const ordinal_type nlbf = _inoutOperator.extent(1);
807  const ordinal_type nrbf = _inoutOperator.extent(2);
808 
809  for (ordinal_type lbf=0;lbf<nlbf;++lbf)
810  for (ordinal_type rbf=0;rbf<nrbf;++rbf)
811  _inoutOperator(cl, lbf, rbf) *= _fieldSigns(cl, lbf);
812  }
813  };
814  }
815 
816  template<typename DeviceType>
817  template<typename inoutOperatorValueType, class ...inoutOperatorProperties,
818  typename fieldSignValueType, class ...fieldSignProperties>
819  void
821  applyLeftFieldSigns( Kokkos::DynRankView<inoutOperatorValueType,inoutOperatorProperties...> inoutOperator,
822  const Kokkos::DynRankView<fieldSignValueType, fieldSignProperties...> fieldSigns ) {
823 
824 #ifdef HAVE_INTREPID2_DEBUG
825  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.rank() != 3, std::invalid_argument,
826  ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): Input operator container must have rank 3.");
827  INTREPID2_TEST_FOR_EXCEPTION( fieldSigns.rank() != 2, std::invalid_argument,
828  ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): Input field signs container must have rank 2.");
829  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.extent(0) != fieldSigns.extent(0), std::invalid_argument,
830  ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): Zeroth dimensions (number of cells) of the operator and field signs containers must agree!");
831  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.extent(1) != fieldSigns.extent(1), std::invalid_argument,
832  ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): First dimensions (number of left fields) of the operator and field signs containers must agree!");
833 #endif
834 
835  constexpr bool are_accessible =
836  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
837  typename decltype(inoutOperator)::memory_space>::accessible &&
838  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
839  typename decltype(fieldSigns)::memory_space>::accessible;
840  static_assert(are_accessible, "FunctionSpaceTools<DeviceType>::applyLeftFieldSigns(..): input/output views' memory spaces are not compatible with DeviceType");
841 
843  <decltype(inoutOperator),decltype(fieldSigns)>;
844 
845  const ordinal_type C = inoutOperator.extent(0);
846  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, C);
847  Kokkos::parallel_for( policy, FunctorType(inoutOperator, fieldSigns) );
848  }
849 
850  // ------------------------------------------------------------------------------------
851 
852  namespace FunctorFunctionSpaceTools {
856  template<typename inoutOperatorViewType,
857  typename fieldSignViewType>
859  inoutOperatorViewType _inoutOperator;
860  const fieldSignViewType _fieldSigns;
861 
862  KOKKOS_INLINE_FUNCTION
863  F_applyRightFieldSigns( inoutOperatorViewType inoutOperator_,
864  fieldSignViewType fieldSigns_ )
865  : _inoutOperator(inoutOperator_), _fieldSigns(fieldSigns_) {}
866 
867  KOKKOS_INLINE_FUNCTION
868  void operator()(const ordinal_type cl) const {
869  const ordinal_type nlbf = _inoutOperator.extent(1);
870  const ordinal_type nrbf = _inoutOperator.extent(2);
871 
872  for (ordinal_type lbf=0;lbf<nlbf;++lbf)
873  for (ordinal_type rbf=0;rbf<nrbf;++rbf)
874  _inoutOperator(cl, lbf, rbf) *= _fieldSigns(cl, rbf);
875  }
876  };
877  }
878 
879  template<typename DeviceType>
880  template<typename inoutOperatorValueType, class ...inoutOperatorProperties,
881  typename fieldSignValueType, class ...fieldSignProperties>
882  void
884  applyRightFieldSigns( Kokkos::DynRankView<inoutOperatorValueType,inoutOperatorProperties...> inoutOperator,
885  const Kokkos::DynRankView<fieldSignValueType, fieldSignProperties...> fieldSigns ) {
886 
887 #ifdef HAVE_INTREPID2_DEBUG
888  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.rank() != 3, std::invalid_argument,
889  ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Input operator container must have rank 3.");
890  INTREPID2_TEST_FOR_EXCEPTION( fieldSigns.rank() != 2, std::invalid_argument,
891  ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Input field signs container must have rank 2.");
892  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.extent(0) != fieldSigns.extent(0), std::invalid_argument,
893  ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Zeroth dimensions (number of cells) of the operator and field signs containers must agree!");
894  INTREPID2_TEST_FOR_EXCEPTION( inoutOperator.extent(2) != fieldSigns.extent(1), std::invalid_argument,
895  ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Second dimension of the operator container and first dimension of the field signs container (number of right fields) must agree!");
896 #endif
897 
898  constexpr bool are_accessible =
899  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
900  typename decltype(inoutOperator)::memory_space>::accessible &&
901  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
902  typename decltype(fieldSigns)::memory_space>::accessible;
903  static_assert(are_accessible, "FunctionSpaceTools<DeviceType>::applyRightFieldSigns(..): input/output views' memory spaces are not compatible with DeviceType");
904 
906  <decltype(inoutOperator),decltype(fieldSigns)>;
907 
908  const ordinal_type C = inoutOperator.extent(0);
909  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, C);
910  Kokkos::parallel_for( policy, FunctorType(inoutOperator, fieldSigns) );
911  }
912 
913  // ------------------------------------------------------------------------------------
914 
915  namespace FunctorFunctionSpaceTools {
919  template<typename inoutFunctionViewType,
920  typename fieldSignViewType>
922  inoutFunctionViewType _inoutFunction;
923  const fieldSignViewType _fieldSigns;
924 
925  KOKKOS_INLINE_FUNCTION
926  F_applyFieldSigns( inoutFunctionViewType inoutFunction_,
927  fieldSignViewType fieldSigns_)
928  : _inoutFunction(inoutFunction_), _fieldSigns(fieldSigns_) {}
929 
930  KOKKOS_INLINE_FUNCTION
931  void operator()(const ordinal_type cl) const {
932  const ordinal_type nbfs = _inoutFunction.extent(1);
933  const ordinal_type npts = _inoutFunction.extent(2);
934  const ordinal_type iend = _inoutFunction.extent(3);
935  const ordinal_type jend = _inoutFunction.extent(4);
936 
937  for (ordinal_type bf=0;bf<nbfs;++bf)
938  for (ordinal_type pt=0;pt<npts;++pt)
939  for (ordinal_type i=0;i<iend;++i)
940  for (ordinal_type j=0;j<jend;++j)
941  _inoutFunction(cl, bf, pt, i, j) *= _fieldSigns(cl, bf);
942  }
943  };
944  }
945 
946  template<typename DeviceType>
947  template<typename inoutFunctionValueType, class ...inoutFunctionProperties,
948  typename fieldSignValueType, class ...fieldSignProperties>
949  void
951  applyFieldSigns( Kokkos::DynRankView<inoutFunctionValueType,inoutFunctionProperties...> inoutFunction,
952  const Kokkos::DynRankView<fieldSignValueType, fieldSignProperties...> fieldSigns ) {
953 
954 #ifdef HAVE_INTREPID2_DEBUG
955  INTREPID2_TEST_FOR_EXCEPTION( inoutFunction.rank() < 2 || inoutFunction.rank() > 5, std::invalid_argument,
956  ">>> ERROR (FunctionSpaceTools::applyFieldSigns): Input function container must have rank 2, 3, 4, or 5.");
957  INTREPID2_TEST_FOR_EXCEPTION( fieldSigns.rank() != 2, std::invalid_argument,
958  ">>> ERROR (FunctionSpaceTools::applyFieldSigns): Input field signs container must have rank 2.");
959  INTREPID2_TEST_FOR_EXCEPTION( inoutFunction.extent(0) != fieldSigns.extent(0), std::invalid_argument,
960  ">>> ERROR (FunctionSpaceTools::applyFieldSigns): Zeroth dimensions (number of integration domains) of the function and field signs containers must agree!");
961  INTREPID2_TEST_FOR_EXCEPTION( inoutFunction.extent(1) != fieldSigns.extent(1), std::invalid_argument,
962  ">>> ERROR (FunctionSpaceTools::applyFieldSigns): First dimensions (number of fields) of the function and field signs containers must agree!");
963 
964 #endif
965 
966  constexpr bool are_accessible =
967  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
968  typename decltype(inoutFunction)::memory_space>::accessible &&
969  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
970  typename decltype(fieldSigns)::memory_space>::accessible;
971  static_assert(are_accessible, "FunctionSpaceTools<DeviceType>::applyFieldSigns(..): input/output views' memory spaces are not compatible with DeviceType");
972 
974  <decltype(inoutFunction),decltype(fieldSigns)>;
975 
976  const ordinal_type C = inoutFunction.extent(0);
977  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, C);
978  Kokkos::parallel_for( policy, FunctorType(inoutFunction, fieldSigns) );
979  }
980 
981  // ------------------------------------------------------------------------------------
982 
983  namespace FunctorFunctionSpaceTools {
984 
988  template<typename outputPointViewType,
989  typename inputCoeffViewType,
990  typename inputFieldViewType>
991  struct F_evaluate {
992  outputPointViewType _outputPointVals;
993  const inputCoeffViewType _inputCoeffs;
994  const inputFieldViewType _inputFields;
995 
996  KOKKOS_INLINE_FUNCTION
997  F_evaluate( outputPointViewType outputPointVals_,
998  inputCoeffViewType inputCoeffs_,
999  inputFieldViewType inputFields_ )
1000  : _outputPointVals(outputPointVals_), _inputCoeffs(inputCoeffs_), _inputFields(inputFields_) {}
1001 
1002  KOKKOS_INLINE_FUNCTION
1003  void operator()(const ordinal_type cl) const {
1004  const ordinal_type nbfs = _inputFields.extent(1);
1005  const ordinal_type npts = _inputFields.extent(2);
1006 
1007  const ordinal_type iend = _inputFields.extent(3);
1008  const ordinal_type jend = _inputFields.extent(4);
1009 
1010  for (ordinal_type bf=0;bf<nbfs;++bf)
1011  for (ordinal_type pt=0;pt<npts;++pt)
1012  for (ordinal_type i=0;i<iend;++i)
1013  for (ordinal_type j=0;j<jend;++j)
1014  _outputPointVals(cl, pt, i, j) += _inputCoeffs(cl, bf) * _inputFields(cl, bf, pt, i, j);
1015  }
1016  };
1017  }
1018 
1019  template<typename DeviceType>
1020  template<typename outputPointValueType, class ...outputPointProperties,
1021  typename inputCoeffValueType, class ...inputCoeffProperties,
1022  typename inputFieldValueType, class ...inputFieldProperties>
1023  void
1025  evaluate( Kokkos::DynRankView<outputPointValueType,outputPointProperties...> outputPointVals,
1026  const Kokkos::DynRankView<inputCoeffValueType, inputCoeffProperties...> inputCoeffs,
1027  const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields ) {
1028 
1029 #ifdef HAVE_INTREPID2_DEBUG
1030  INTREPID2_TEST_FOR_EXCEPTION( inputFields.rank() < 3 || inputFields.rank() > 5, std::invalid_argument,
1031  ">>> ERROR (FunctionSpaceTools::evaluate): Input fields container must have rank 3, 4, or 5.");
1032  INTREPID2_TEST_FOR_EXCEPTION( inputCoeffs.rank() != 2, std::invalid_argument,
1033  ">>> ERROR (FunctionSpaceTools::evaluate): Input coefficient container must have rank 2.");
1034  INTREPID2_TEST_FOR_EXCEPTION( outputPointVals.rank() != (inputFields.rank()-1), std::invalid_argument,
1035  ">>> ERROR (FunctionSpaceTools::evaluate): Output values container must have rank one less than the rank of the input fields container.");
1036  INTREPID2_TEST_FOR_EXCEPTION( inputCoeffs.extent(0) != inputFields.extent(0), std::invalid_argument,
1037  ">>> ERROR (FunctionSpaceTools::evaluate): Zeroth dimensions (number of cells) of the coefficient and fields input containers must agree!");
1038  INTREPID2_TEST_FOR_EXCEPTION( inputCoeffs.extent(1) != inputFields.extent(1), std::invalid_argument,
1039  ">>> ERROR (FunctionSpaceTools::evaluate): First dimensions (number of fields) of the coefficient and fields input containers must agree!");
1040  INTREPID2_TEST_FOR_EXCEPTION( outputPointVals.extent(0) != inputFields.extent(0), std::invalid_argument,
1041  ">>> ERROR (FunctionSpaceTools::evaluate): Zeroth dimensions (number of cells) of the input fields container and the output values container must agree!");
1042  for (size_type i=1;i<outputPointVals.rank();++i)
1043  INTREPID2_TEST_FOR_EXCEPTION( outputPointVals.extent(i) != inputFields.extent(i+1), std::invalid_argument,
1044  ">>> ERROR (FunctionSpaceTools::evaluate): outputPointVals dimension(i) does not match to inputFields dimension(i+1).");
1045 #endif
1046 
1047  constexpr bool are_accessible =
1048  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
1049  typename decltype(outputPointVals)::memory_space>::accessible &&
1050  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
1051  typename decltype(inputCoeffs)::memory_space>::accessible &&
1052  Kokkos::Impl::MemorySpaceAccess<MemSpaceType,
1053  typename decltype(inputFields)::memory_space>::accessible;
1054  static_assert(are_accessible, "FunctionSpaceTools<DeviceType>::evaluate(..): input/output views' memory spaces are not compatible with DeviceType");
1055 
1056  using FunctorType = FunctorFunctionSpaceTools::F_evaluate
1057  <decltype(outputPointVals),decltype(inputCoeffs),decltype(inputFields)>;
1058 
1059  const ordinal_type C = inputFields.extent(0);
1060  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, C);
1061  Kokkos::parallel_for( policy, FunctorType(outputPointVals, inputCoeffs, inputFields) );
1062  }
1063 
1064  // ------------------------------------------------------------------------------------
1065 
1066 } // end namespace Intrepid2
1067 
1068 #endif
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 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 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 (...
Functor for calculation HGRADtransformGRAD, see Intrepid2::FunctionSpaceTools for more...
Functor for applyLeftFieldSigns, see Intrepid2::FunctionSpaceTools for more.
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 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 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 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 HGRADtransformCURL(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 2D curl field in the H-grad space, defined at points on a reference cell...
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 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 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 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 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 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 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 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.
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 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...
Defines TensorArgumentIterator, which allows systematic enumeration of a TensorData object...
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 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 ...
Functor for applyRightFieldSigns, see Intrepid2::FunctionSpaceTools for more.
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 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 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 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...
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 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 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 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...
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 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 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 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 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...
Defines the Intrepid2::FunctorIterator class, as well as the Intrepid2::functor_returns_ref SFINAE he...
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 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 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 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 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 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...
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 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...
Functor to evaluate functions, see Intrepid2::FunctionSpaceTools for more.
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 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 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 3D curl field in the H-curl 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)
Dot product of data and data; please read the description below.
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...
Functor for calculation of cell measure, see Intrepid2::FunctionSpaceTools for more.
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 ...
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 clone(Kokkos::DynRankView< outputValueType, outputProperties... > output, const Kokkos::DynRankView< inputValueType, inputProperties... > input)
Clone input array.
Functor for applyFieldSigns, see Intrepid2::FunctionSpaceTools for more.
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 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...