Intrepid2
Intrepid2_ArrayToolsDefContractions.hpp
Go to the documentation of this file.
1// @HEADER
2// *****************************************************************************
3// Intrepid2 Package
4//
5// Copyright 2007 NTESS and the Intrepid2 contributors.
6// SPDX-License-Identifier: BSD-3-Clause
7// *****************************************************************************
8// @HEADER
9
16#ifndef __INTREPID2_ARRAYTOOLS_DEF_CONTRACTIONS_HPP__
17#define __INTREPID2_ARRAYTOOLS_DEF_CONTRACTIONS_HPP__
18
19namespace Intrepid2 {
20
21
22 namespace FunctorArrayTools {
26 template < typename outFieldViewType , typename leftFieldViewType , typename rightFieldViewType >
28 outFieldViewType _outputFields;
29 leftFieldViewType _leftFields;
30 rightFieldViewType _rightFields;
31 const bool _sumInto;
32
33 KOKKOS_INLINE_FUNCTION
34 F_contractFieldFieldScalar(outFieldViewType outputFields_,
35 leftFieldViewType leftFields_,
36 rightFieldViewType rightFields_,
37 const bool sumInto_)
38 : _outputFields(outputFields_), _leftFields(leftFields_), _rightFields(rightFields_), _sumInto(sumInto_) {}
39
40 KOKKOS_INLINE_FUNCTION
41 void operator()(const size_type cl, const size_type lbf, const size_type rbf) const {
42 const ordinal_type npts = _leftFields.extent(2);
43
44 _outputFields( cl, lbf, rbf ) *= (_sumInto ? 1.0 : 0.0);
45 for (ordinal_type qp = 0; qp < npts; ++qp)
46 _outputFields( cl, lbf, rbf ) += _leftFields(cl, lbf, qp)*_rightFields(cl, rbf, qp);
47 }
48 };
49
50 template < typename outFieldViewType , typename leftFieldViewType , typename rightFieldViewType >
52 outFieldViewType _outputFields;
53 leftFieldViewType _leftFields;
54 rightFieldViewType _rightFields;
55 const bool _sumInto;
56
57 KOKKOS_INLINE_FUNCTION
58 F_contractFieldFieldVector(outFieldViewType outputFields_,
59 leftFieldViewType leftFields_,
60 rightFieldViewType rightFields_,
61 const bool sumInto_)
62 : _outputFields(outputFields_), _leftFields(leftFields_), _rightFields(rightFields_), _sumInto(sumInto_) {}
63
64 KOKKOS_INLINE_FUNCTION
65 void operator()(const size_type cl, const size_type lbf, const size_type rbf) const {
66 const ordinal_type npts = _leftFields.extent(2);
67 const ordinal_type iend = _leftFields.extent(3);
68
69 _outputFields( cl, lbf, rbf ) *= (_sumInto ? 1.0 : 0.0);
70 for (ordinal_type qp = 0; qp < npts; ++qp)
71 for (ordinal_type i = 0; i < iend; ++i)
72 _outputFields( cl, lbf, rbf ) += _leftFields(cl, lbf, qp, i)*_rightFields(cl, rbf, qp, i);
73 }
74 };
75
76 template < typename outFieldViewType , typename leftFieldViewType , typename rightFieldViewType >
78 outFieldViewType _outputFields;
79 leftFieldViewType _leftFields;
80 rightFieldViewType _rightFields;
81 const bool _sumInto;
82
83 KOKKOS_INLINE_FUNCTION
84 F_contractFieldFieldTensor(outFieldViewType outputFields_,
85 leftFieldViewType leftFields_,
86 rightFieldViewType rightFields_,
87 const bool sumInto_)
88 : _outputFields(outputFields_), _leftFields(leftFields_), _rightFields(rightFields_), _sumInto(sumInto_) {}
89
90 KOKKOS_INLINE_FUNCTION
91 void operator()(const size_type cl, const size_type lbf, const size_type rbf) const {
92 const ordinal_type npts = _leftFields.extent(2);
93 const ordinal_type iend = _leftFields.extent(3);
94 const ordinal_type jend = _leftFields.extent(4);
95
96 _outputFields( cl, lbf, rbf ) *= (_sumInto ? 1.0 : 0.0);
97 for (ordinal_type qp = 0; qp < npts; ++qp)
98 for (ordinal_type i = 0; i < iend; ++i)
99 for (ordinal_type j = 0; j < jend; ++j)
100 _outputFields( cl, lbf, rbf ) += _leftFields(cl, lbf, qp, i, j)*_rightFields(cl, rbf, qp, i, j);
101 }
102 };
103
104
105 } //end namespace
106
107 template<typename DeviceType>
108 template<typename outputFieldValueType, class ...outputFieldProperties,
109 typename leftFieldValueType, class ...leftFieldProperties,
110 typename rightFieldValueType, class ...rightFieldProperties>
111 void
113 contractFieldField( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
114 const Kokkos::DynRankView<leftFieldValueType, leftFieldProperties...> leftFields,
115 const Kokkos::DynRankView<rightFieldValueType, rightFieldProperties...> rightFields,
116 const bool sumInto ) {
117
118 using outFieldViewType = Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...>;
119 using leftFieldViewType = Kokkos::DynRankView<leftFieldValueType,leftFieldProperties...>;
120 using rightFieldViewType = Kokkos::DynRankView<rightFieldValueType,rightFieldProperties...>;
121
122 using range_policy_type = Kokkos::MDRangePolicy< ExecSpaceType, Kokkos::Rank<3>, Kokkos::IndexType<ordinal_type> >;
123 const range_policy_type policy( { 0, 0, 0 },
124 { /*C*/ leftFields.extent(0), /*F*/ leftFields.extent(1), /*F*/ rightFields.extent(1) } );
125 if (rightFields.rank() == 3) {
126 using FunctorType = FunctorArrayTools::F_contractFieldFieldScalar<outFieldViewType, leftFieldViewType, rightFieldViewType>;
127 Kokkos::parallel_for( policy, FunctorType(outputFields, leftFields, rightFields, sumInto) );
128 } else
129 if (rightFields.rank() == 4) {
130 using FunctorType = FunctorArrayTools::F_contractFieldFieldVector<outFieldViewType, leftFieldViewType, rightFieldViewType>;
131 Kokkos::parallel_for( policy, FunctorType(outputFields, leftFields, rightFields, sumInto) );
132 } else {
133 using FunctorType = FunctorArrayTools::F_contractFieldFieldTensor<outFieldViewType, leftFieldViewType, rightFieldViewType>;
134 Kokkos::parallel_for( policy, FunctorType(outputFields, leftFields, rightFields, sumInto) );
135 }
136 }
137
138
139 namespace FunctorArrayTools {
143 template < typename outputFieldsViewType , typename inputDataViewType , typename inputFieldsViewType >
145 outputFieldsViewType _outputFields;
146 inputDataViewType _inputData;
147 inputFieldsViewType _inputFields;
148 const bool _sumInto;
149
150 KOKKOS_INLINE_FUNCTION
151 F_contractDataFieldScalar(outputFieldsViewType outputFields_,
152 inputDataViewType inputData_,
153 inputFieldsViewType inputFields_,
154 const bool sumInto_)
155 : _outputFields(outputFields_), _inputData(inputData_), _inputFields(inputFields_), _sumInto(sumInto_) {}
156
157 KOKKOS_DEFAULTED_FUNCTION
158 ~F_contractDataFieldScalar() = default;
159
160 KOKKOS_INLINE_FUNCTION
161 void operator()(const size_type cl, const size_type bf) const {
162 const size_type npts = _inputFields.extent(2);
163 _outputFields(cl, bf) *= (_sumInto ? 1 : 0);
164
165 if(_inputData.extent(1) != 1)
166 for (size_type qp = 0; qp < npts; ++qp)
167 _outputFields(cl, bf) += _inputFields(cl, bf, qp) * _inputData(cl, qp);
168 else
169 for (size_type qp = 0; qp < npts; ++qp)
170 _outputFields(cl, bf) += _inputFields(cl, bf, qp) * _inputData(cl, 0);
171 }
172 };
173
174 template < typename outputFieldsViewType , typename inputDataViewType , typename inputFieldsViewType >
176 outputFieldsViewType _outputFields;
177 inputDataViewType _inputData;
178 inputFieldsViewType _inputFields;
179 const bool _sumInto;
180
181 KOKKOS_INLINE_FUNCTION
182 F_contractDataFieldVector(outputFieldsViewType outputFields_,
183 inputDataViewType inputData_,
184 inputFieldsViewType inputFields_,
185 const bool sumInto_)
186 : _outputFields(outputFields_), _inputData(inputData_), _inputFields(inputFields_), _sumInto(sumInto_) {}
187
188 KOKKOS_DEFAULTED_FUNCTION
189 ~F_contractDataFieldVector() = default;
190
191 KOKKOS_INLINE_FUNCTION
192 void operator()(const size_type cl, const size_type bf) const {
193 const size_type npts = _inputFields.extent(2);
194 const ordinal_type iend = _inputFields.extent(3);
195
196 _outputFields(cl, bf) *= (_sumInto ? 1 : 0);
197
198 if(_inputData.extent(1) != 1)
199 for (size_type qp = 0; qp < npts; ++qp)
200 for (ordinal_type i = 0; i < iend; ++i)
201 _outputFields(cl, bf) += _inputFields(cl, bf, qp, i) * _inputData(cl, qp, i);
202 else
203 for (size_type qp = 0; qp < npts; ++qp)
204 for (ordinal_type i = 0; i < iend; ++i)
205 _outputFields(cl, bf) += _inputFields(cl, bf, qp, i) * _inputData(cl, 0, i);
206 }
207 };
208
209
210 template < typename outputFieldsViewType , typename inputDataViewType , typename inputFieldsViewType >
212 outputFieldsViewType _outputFields;
213 inputDataViewType _inputData;
214 inputFieldsViewType _inputFields;
215 const bool _sumInto;
216
217 KOKKOS_INLINE_FUNCTION
218 F_contractDataFieldTensor(outputFieldsViewType outputFields_,
219 inputDataViewType inputData_,
220 inputFieldsViewType inputFields_,
221 const bool sumInto_)
222 : _outputFields(outputFields_), _inputData(inputData_), _inputFields(inputFields_), _sumInto(sumInto_) {}
223
224 KOKKOS_DEFAULTED_FUNCTION
225 ~F_contractDataFieldTensor() = default;
226
227 KOKKOS_INLINE_FUNCTION
228 void operator()(const size_type cl, const size_type bf) const {
229 const size_type npts = _inputFields.extent(2);
230 const ordinal_type iend = _inputFields.extent(3);
231 const ordinal_type jend = _inputFields.extent(4);
232
233 _outputFields(cl, bf) *= (_sumInto ? 1 : 0);
234
235 if(_inputData.extent(1) != 1)
236 for (size_type qp = 0; qp < npts; ++qp)
237 for (ordinal_type i = 0; i < iend; ++i)
238 for (ordinal_type j = 0; j < jend; ++j)
239 _outputFields(cl, bf) += _inputFields(cl, bf, qp, i, j) * _inputData(cl, qp, i, j);
240 else
241 for (size_type qp = 0; qp < npts; ++qp)
242 for (ordinal_type i = 0; i < iend; ++i)
243 for (ordinal_type j = 0; j < jend; ++j)
244 _outputFields(cl, bf) += _inputFields(cl, bf, qp, i, j) * _inputData(cl, 0, i, j);
245 }
246 };
247
248 } //namespace
249
250 template<typename DeviceType>
251 template<typename outputFieldValueType, class ...outputFieldProperties,
252 typename inputDataValueType, class ...inputDataProperties,
253 typename inputFieldValueType, class ...inputFieldProperties>
254 void
256 contractDataField( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
257 const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
258 const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields,
259 const bool sumInto ) {
260
261 using outputFieldsViewType = Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...>;
262 using inputDataViewType = Kokkos::DynRankView<inputDataValueType, inputDataProperties...>;
263 using inputFieldsViewType = Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...>;
264
265 using range_policy_type = Kokkos::MDRangePolicy< ExecSpaceType, Kokkos::Rank<2>, Kokkos::IndexType<ordinal_type> >;
266 const range_policy_type policy( { 0, 0 }, { /*C*/ inputFields.extent(0), /*F*/ inputFields.extent(1)} );
267
268 if (inputFields.rank() == 3) {
269 using FunctorType = FunctorArrayTools::F_contractDataFieldScalar<outputFieldsViewType, inputDataViewType, inputFieldsViewType>;
270 Kokkos::parallel_for( policy, FunctorType(outputFields, inputData, inputFields, sumInto) );
271 }
272 else if (inputFields.rank() == 4) {
273 using FunctorType = FunctorArrayTools::F_contractDataFieldVector<outputFieldsViewType, inputDataViewType, inputFieldsViewType>;
274 Kokkos::parallel_for( policy, FunctorType(outputFields, inputData, inputFields, sumInto) );
275 }
276 else {
277 using FunctorType = FunctorArrayTools::F_contractDataFieldTensor<outputFieldsViewType, inputDataViewType, inputFieldsViewType>;
278 Kokkos::parallel_for( policy, FunctorType(outputFields, inputData, inputFields, sumInto) );
279 }
280 }
281
282
283 namespace FunctorArrayTools {
287 template < typename outputDataViewType , typename inputDataLeftViewType , typename inputDataRightViewType >
289 outputDataViewType _outputData;
290 inputDataLeftViewType _inputDataLeft;
291 inputDataRightViewType _inputDataRight;
292 const bool _sumInto;
293
294 KOKKOS_INLINE_FUNCTION
295 F_contractDataDataScalar(outputDataViewType outputData_,
296 inputDataLeftViewType inputDataLeft_,
297 inputDataRightViewType inputDataRight_,
298 const bool sumInto_)
299 : _outputData(outputData_), _inputDataLeft(inputDataLeft_), _inputDataRight(inputDataRight_), _sumInto(sumInto_) {}
300
301 KOKKOS_DEFAULTED_FUNCTION
302 ~F_contractDataDataScalar() = default;
303
304 KOKKOS_INLINE_FUNCTION
305 void operator()(const size_type cl) const {
306 size_type npts = _inputDataLeft.extent(1);
307 _outputData(cl) *= (_sumInto ? 1 : 0);
308 for (size_type qp = 0; qp < npts; ++qp)
309 _outputData(cl) += _inputDataLeft(cl, qp)*_inputDataRight(cl, qp);
310 }
311 };
312
313 template < typename outputDataViewType , typename inputDataLeftViewType , typename inputDataRightViewType >
315 outputDataViewType _outputData;
316 inputDataLeftViewType _inputDataLeft;
317 inputDataRightViewType _inputDataRight;
318 const bool _sumInto;
319
320 KOKKOS_INLINE_FUNCTION
321 F_contractDataDataVector(outputDataViewType outputData_,
322 inputDataLeftViewType inputDataLeft_,
323 inputDataRightViewType inputDataRight_,
324 const bool sumInto_)
325 : _outputData(outputData_), _inputDataLeft(inputDataLeft_), _inputDataRight(inputDataRight_), _sumInto(sumInto_) {}
326
327 KOKKOS_DEFAULTED_FUNCTION
328 ~F_contractDataDataVector() = default;
329
330 KOKKOS_INLINE_FUNCTION
331 void operator()(const size_type cl) const {
332 size_type npts = _inputDataLeft.extent(1);
333 ordinal_type iend = _inputDataLeft.extent(2);
334
335 _outputData(cl) *= (_sumInto ? 1 : 0);
336 for (size_type qp = 0; qp < npts; ++qp)
337 for (ordinal_type i = 0; i < iend; ++i)
338 _outputData(cl) += _inputDataLeft(cl, qp, i)*_inputDataRight(cl, qp, i);
339 }
340 };
341
342 template < typename outputDataViewType , typename inputDataLeftViewType , typename inputDataRightViewType >
344 outputDataViewType _outputData;
345 inputDataLeftViewType _inputDataLeft;
346 inputDataRightViewType _inputDataRight;
347 const bool _sumInto;
348
349 KOKKOS_INLINE_FUNCTION
350 F_contractDataDataTensor(outputDataViewType outputData_,
351 inputDataLeftViewType inputDataLeft_,
352 inputDataRightViewType inputDataRight_,
353 const bool sumInto_)
354 : _outputData(outputData_), _inputDataLeft(inputDataLeft_), _inputDataRight(inputDataRight_), _sumInto(sumInto_) {}
355
356 KOKKOS_DEFAULTED_FUNCTION
357 ~F_contractDataDataTensor() = default;
358
359 KOKKOS_INLINE_FUNCTION
360 void operator()(const size_type cl) const {
361 size_type npts = _inputDataLeft.extent(1);
362 ordinal_type iend = _inputDataLeft.extent(2);
363 ordinal_type jend = _inputDataLeft.extent(3);
364
365 _outputData(cl) *= (_sumInto ? 1 : 0);
366 for (size_type qp = 0; qp < npts; ++qp)
367 for (ordinal_type i = 0; i < iend; ++i)
368 for (ordinal_type j = 0; j < jend; ++j)
369 _outputData(cl) += _inputDataLeft(cl, qp, i, j)*_inputDataRight(cl, qp, i, j);
370 }
371 };
372 } //namespace
373
374 template<typename DeviceType>
375 template<typename outputDataValueType, class ...outputDataProperties,
376 typename inputDataLeftValueType, class ...inputDataLeftProperties,
377 typename inputDataRightValueType, class ...inputDataRightProperties>
378 void
380 contractDataData( Kokkos::DynRankView<outputDataValueType, outputDataProperties...> outputData,
381 const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
382 const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight,
383 const bool sumInto ) {
384 using outputDataViewType = Kokkos::DynRankView<outputDataValueType, outputDataProperties...>;
385 using inputDataLeftViewType = Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...>;
386 using inputDataRightViewType = Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...>;
387
388 const size_type loopSize = inputDataLeft.extent(0);
389 Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, loopSize);
390
391 if (inputDataLeft.rank() == 2) {
393 Kokkos::parallel_for( policy, FunctorType(outputData, inputDataLeft, inputDataRight, sumInto) );
394 }
395 else if (inputDataLeft.rank() == 3) {
396 using FunctorType = FunctorArrayTools::F_contractDataDataVector<outputDataViewType, inputDataLeftViewType, inputDataRightViewType>;
397 Kokkos::parallel_for( policy, FunctorType(outputData, inputDataLeft, inputDataRight, sumInto) );
398 }
399 else {
400 using FunctorType = FunctorArrayTools::F_contractDataDataTensor<outputDataViewType, inputDataLeftViewType, inputDataRightViewType>;
401 Kokkos::parallel_for( policy, FunctorType(outputData, inputDataLeft, inputDataRight, sumInto) );
402 }
403 }
404
405
406
407 template<typename DeviceType>
408 template<typename outputFieldValueType, class ...outputFieldProperties,
409 typename leftFieldValueType, class ...leftFieldProperties,
410 typename rightFieldValueType, class ...rightFieldProperties>
411 void
413 contractFieldFieldScalar( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
414 const Kokkos::DynRankView<leftFieldValueType, leftFieldProperties...> leftFields,
415 const Kokkos::DynRankView<rightFieldValueType, rightFieldProperties...> rightFields,
416 const bool sumInto ) {
417
418#ifdef HAVE_INTREPID2_DEBUG
419 {
420 INTREPID2_TEST_FOR_EXCEPTION( leftFields.rank() != 3, std::invalid_argument,
421 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Rank of the left input argument must equal 3!");
422 INTREPID2_TEST_FOR_EXCEPTION( rightFields.rank() != 3, std::invalid_argument,
423 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Rank of right input argument must equal 3!");
424 INTREPID2_TEST_FOR_EXCEPTION( outputFields.rank() != 3, std::invalid_argument,
425 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Rank of output argument must equal 3!");
426 INTREPID2_TEST_FOR_EXCEPTION( leftFields.extent(0) != rightFields.extent(0), std::invalid_argument,
427 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
428 INTREPID2_TEST_FOR_EXCEPTION( leftFields.extent(2) != rightFields.extent(2), std::invalid_argument,
429 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Second dimensions (numbers of integration points) of the left and right input containers must agree!");
430 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(0) != rightFields.extent(0), std::invalid_argument,
431 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
432 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(1) != leftFields.extent(1), std::invalid_argument,
433 ">>> ERROR (ArrayTools::contractFieldFieldScalar): First dimension of output container and first dimension of left input container must agree!");
434 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(2) != rightFields.extent(1), std::invalid_argument,
435 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Second dimension of output container and first dimension of right input container must agree!");
436
437 }
438#endif
439
441 leftFields,
442 rightFields,
443 sumInto );
444 }
445
446
447 template<typename DeviceType>
448 template<typename outputFieldValueType, class ...outputFieldProperties,
449 typename leftFieldValueType, class ...leftFieldProperties,
450 typename rightFieldValueType, class ...rightFieldProperties>
451 void
453 contractFieldFieldVector( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
454 const Kokkos::DynRankView<leftFieldValueType, leftFieldProperties...> leftFields,
455 const Kokkos::DynRankView<rightFieldValueType, rightFieldProperties...> rightFields,
456 const bool sumInto ) {
457
458#ifdef HAVE_INTREPID2_DEBUG
459 {
460 INTREPID2_TEST_FOR_EXCEPTION( leftFields.rank() != 4, std::invalid_argument,
461 ">>> ERROR (ArrayTools::contractFieldFieldVector): Rank of the left input argument must equal 4!");
462 INTREPID2_TEST_FOR_EXCEPTION( rightFields.rank() != 4, std::invalid_argument,
463 ">>> ERROR (ArrayTools::contractFieldFieldVector): Rank of right input argument must equal 4!");
464 INTREPID2_TEST_FOR_EXCEPTION( outputFields.rank() != 3, std::invalid_argument,
465 ">>> ERROR (ArrayTools::contractFieldFieldVector): Rank of output argument must equal 3!");
466 INTREPID2_TEST_FOR_EXCEPTION( leftFields.extent(0) != rightFields.extent(0), std::invalid_argument,
467 ">>> ERROR (ArrayTools::contractFieldFieldVector): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
468 INTREPID2_TEST_FOR_EXCEPTION( leftFields.extent(2) != rightFields.extent(2), std::invalid_argument,
469 ">>> ERROR (ArrayTools::contractFieldFieldVector): Second dimensions (numbers of integration points) of the left and right input containers must agree!");
470 INTREPID2_TEST_FOR_EXCEPTION( leftFields.extent(3) != rightFields.extent(3), std::invalid_argument,
471 ">>> ERROR (ArrayTools::contractFieldFieldVector): Third dimensions (numbers of vector components) of the left and right input containers must agree!");
472 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(0) != rightFields.extent(0), std::invalid_argument,
473 ">>> ERROR (ArrayTools::contractFieldFieldVector): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
474 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(1) != leftFields.extent(1), std::invalid_argument,
475 ">>> ERROR (ArrayTools::contractFieldFieldVector): First dimension of output container and first dimension of left input container must agree!");
476 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(2) != rightFields.extent(1), std::invalid_argument,
477 ">>> ERROR (ArrayTools::contractFieldFieldVector): Second dimension of output container and first dimension of right input container must agree!");
478 }
479#endif
480
482 leftFields,
483 rightFields,
484 sumInto );
485 }
486
487
488 template<typename DeviceType>
489 template<typename outputFieldValueType, class ...outputFieldProperties,
490 typename leftFieldValueType, class ...leftFieldProperties,
491 typename rightFieldValueType, class ...rightFieldProperties>
492 void
494 contractFieldFieldTensor( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
495 const Kokkos::DynRankView<leftFieldValueType, leftFieldProperties...> leftFields,
496 const Kokkos::DynRankView<rightFieldValueType, rightFieldProperties...> rightFields,
497 const bool sumInto ) {
498
499#ifdef HAVE_INTREPID2_DEBUG
500 {
501 INTREPID2_TEST_FOR_EXCEPTION( leftFields.rank() != 5, std::invalid_argument,
502 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Rank of the left input argument must equal 5!");
503 INTREPID2_TEST_FOR_EXCEPTION( rightFields.rank() != 5, std::invalid_argument,
504 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Rank of right input argument must equal 5!");
505 INTREPID2_TEST_FOR_EXCEPTION( outputFields.rank() != 3, std::invalid_argument,
506 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Rank of output argument must equal 3!");
507 INTREPID2_TEST_FOR_EXCEPTION( leftFields.extent(0) != rightFields.extent(0), std::invalid_argument,
508 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
509 INTREPID2_TEST_FOR_EXCEPTION( leftFields.extent(2) != rightFields.extent(2), std::invalid_argument,
510 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Second dimensions (numbers of integration points) of the left and right input containers must agree!");
511 INTREPID2_TEST_FOR_EXCEPTION( leftFields.extent(3) != rightFields.extent(3), std::invalid_argument,
512 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Third dimensions (first tensor dimensions) of the left and right input containers must agree!");
513 INTREPID2_TEST_FOR_EXCEPTION( leftFields.extent(4) != rightFields.extent(4), std::invalid_argument,
514 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Fourth dimensions (second tensor dimensions) of the left and right input containers must agree!");
515 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(0) != rightFields.extent(0), std::invalid_argument,
516 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
517 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(1) != leftFields.extent(1), std::invalid_argument,
518 ">>> ERROR (ArrayTools::contractFieldFieldTensor): First dimension of output container and first dimension of left input container must agree!");
519 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(2) != rightFields.extent(1), std::invalid_argument,
520 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Second dimension of output container and first dimension of right input container must agree!");
521 }
522#endif
523
525 leftFields,
526 rightFields,
527 sumInto );
528 }
529
530
531 template<typename DeviceType>
532 template<typename outputFieldValueType, class ...outputFieldProperties,
533 typename inputDataValueType, class ...inputDataProperties,
534 typename inputFieldValueType, class ...inputFieldProperties>
535 void
537 contractDataFieldScalar( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
538 const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
539 const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields,
540 const bool sumInto ) {
541
542#ifdef HAVE_INTREPID2_DEBUG
543 {
544 INTREPID2_TEST_FOR_EXCEPTION( inputFields.rank() != 3, std::invalid_argument,
545 ">>> ERROR (ArrayTools::contractDataFieldScalar): Rank of the fields input argument must equal 3!");
546 INTREPID2_TEST_FOR_EXCEPTION( inputData.rank() != 2, std::invalid_argument,
547 ">>> ERROR (ArrayTools::contractDataFieldScalar): Rank of the data input argument must equal 2!");
548 INTREPID2_TEST_FOR_EXCEPTION( outputFields.rank() != 2, std::invalid_argument,
549 ">>> ERROR (ArrayTools::contractDataFieldScalar): Rank of output argument must equal 2!");
550 INTREPID2_TEST_FOR_EXCEPTION( inputFields.extent(0) != inputData.extent(0), std::invalid_argument,
551 ">>> ERROR (ArrayTools::contractDataFieldScalar): Zeroth dimensions (number of integration domains) of the fields and data input containers must agree!");
552
553 INTREPID2_TEST_FOR_EXCEPTION( inputData.extent(1) != inputFields.extent(2) &&
554 inputData.extent(1) != 1, std::invalid_argument,
555 ">>> ERROR (ArrayTools::contractDataFieldScalar): Second dimension of fields input container and first dimension of data input container (number of integration points) must agree or first data dimension must be 1!");
556 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(0) != inputFields.extent(0), std::invalid_argument,
557 ">>> ERROR (ArrayTools::contractDataFieldScalar): Zeroth dimensions (numbers of integration domains) of the fields input and output containers must agree!");
558 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(1) != inputFields.extent(1), std::invalid_argument,
559 ">>> ERROR (ArrayTools::contractDataFieldScalar): First dimensions (number of fields) of the fields input and output containers must agree!");
560 }
561#endif
562
564 inputData,
565 inputFields,
566 sumInto );
567 }
568
569
570 template<typename DeviceType>
571 template<typename outputFieldValueType, class ...outputFieldProperties,
572 typename inputDataValueType, class ...inputDataProperties,
573 typename inputFieldValueType, class ...inputFieldProperties>
574 void
576 contractDataFieldVector( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
577 const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
578 const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields,
579 const bool sumInto ) {
580
581#ifdef HAVE_INTREPID2_DEBUG
582 {
583 INTREPID2_TEST_FOR_EXCEPTION( inputFields.rank() != 4, std::invalid_argument,
584 ">>> ERROR (ArrayTools::contractDataFieldVector): Rank of the fields input argument must equal 4!");
585 INTREPID2_TEST_FOR_EXCEPTION( inputData.rank() != 3, std::invalid_argument,
586 ">>> ERROR (ArrayTools::contractDataFieldVector): Rank of the data input argument must equal 3!");
587 INTREPID2_TEST_FOR_EXCEPTION( outputFields.rank() != 2, std::invalid_argument,
588 ">>> ERROR (ArrayTools::contractDataFieldVector): Rank of output argument must equal 2!");
589 INTREPID2_TEST_FOR_EXCEPTION( inputFields.extent(0) != inputData.extent(0), std::invalid_argument,
590 ">>> ERROR (ArrayTools::contractDataFieldVector): Zeroth dimensions (number of integration domains) of the fields and data input containers must agree!");
591 INTREPID2_TEST_FOR_EXCEPTION( inputData.extent(1) != inputFields.extent(2) &&
592 inputData.extent(1) != 1, std::invalid_argument,
593 ">>> ERROR (ArrayTools::contractDataFieldVector): Second dimension of the fields input container and first dimension of data input container (number of integration points) must agree or first data dimension must be 1!");
594 INTREPID2_TEST_FOR_EXCEPTION( inputFields.extent(3) != inputData.extent(2), std::invalid_argument,
595 ">>> ERROR (ArrayTools::contractDataFieldVector): Third dimension of the fields input container and second dimension of data input container (vector index) must agree!");
596 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(0) != inputFields.extent(0), std::invalid_argument,
597 ">>> ERROR (ArrayTools::contractDataFieldVector): Zeroth dimensions (numbers of integration domains) of the fields input and output containers must agree!");
598 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(1) != inputFields.extent(1), std::invalid_argument,
599 ">>> ERROR (ArrayTools::contractDataFieldVector): First dimensions of output container and fields input container (number of fields) must agree!");
600 }
601#endif
602
604 inputData,
605 inputFields,
606 sumInto );
607 }
608
609
610
611 template<typename DeviceType>
612 template<typename outputFieldValueType, class ...outputFieldProperties,
613 typename inputDataValueType, class ...inputDataProperties,
614 typename inputFieldValueType, class ...inputFieldProperties>
615 void
617 contractDataFieldTensor( Kokkos::DynRankView<outputFieldValueType,outputFieldProperties...> outputFields,
618 const Kokkos::DynRankView<inputDataValueType, inputDataProperties...> inputData,
619 const Kokkos::DynRankView<inputFieldValueType, inputFieldProperties...> inputFields,
620 const bool sumInto ) {
621
622#ifdef HAVE_INTREPID2_DEBUG
623 {
624 INTREPID2_TEST_FOR_EXCEPTION( inputFields.rank() != 5, std::invalid_argument,
625 ">>> ERROR (ArrayTools::contractDataFieldTensor): Rank of the fields input argument must equal 5!");
626 INTREPID2_TEST_FOR_EXCEPTION( inputData.rank() != 4, std::invalid_argument,
627 ">>> ERROR (ArrayTools::contractDataFieldTensor): Rank of the data input argument must equal 4!");
628 INTREPID2_TEST_FOR_EXCEPTION( outputFields.rank() != 2, std::invalid_argument,
629 ">>> ERROR (ArrayTools::contractDataFieldTensor): Rank of output argument must equal 2!");
630 INTREPID2_TEST_FOR_EXCEPTION( inputFields.extent(0) != inputData.extent(0), std::invalid_argument,
631 ">>> ERROR (ArrayTools::contractDataFieldTensor): Zeroth dimensions (number of integration domains) of the fields and data input containers must agree!");
632 INTREPID2_TEST_FOR_EXCEPTION( inputData.extent(1) != inputFields.extent(2) && inputData.extent(1) != 1, std::invalid_argument,
633 ">>> ERROR (ArrayTools::contractDataFieldTensor): Second dimension of the fields input container and first dimension of data input container (number of integration points) must agree or first data dimension must be 1!");
634 INTREPID2_TEST_FOR_EXCEPTION( inputFields.extent(3) != inputData.extent(2), std::invalid_argument,
635 ">>> ERROR (ArrayTools::contractDataFieldTensor): Third dimension of the fields input container and second dimension of data input container (first tensor dimension) must agree!");
636 INTREPID2_TEST_FOR_EXCEPTION( inputFields.extent(4) != inputData.extent(3), std::invalid_argument,
637 ">>> ERROR (ArrayTools::contractDataFieldTensor): Fourth dimension of the fields input container and third dimension of data input container (second tensor dimension) must agree!");
638 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(0) != inputFields.extent(0), std::invalid_argument,
639 ">>> ERROR (ArrayTools::contractDataFieldTensor): Zeroth dimensions (numbers of integration domains) of the fields input and output containers must agree!");
640 INTREPID2_TEST_FOR_EXCEPTION( outputFields.extent(1) != inputFields.extent(1), std::invalid_argument,
641 ">>> ERROR (ArrayTools::contractDataFieldTensor): First dimensions (number of fields) of output container and fields input container must agree!");
642 }
643#endif
644
646 inputData,
647 inputFields,
648 sumInto );
649 }
650
651
652
653 template<typename DeviceType>
654 template<typename outputDataValueType, class ...outputDataProperties,
655 typename inputDataLeftValueType, class ...inputDataLeftProperties,
656 typename inputDataRightValueType, class ...inputDataRightProperties>
657 void
659 contractDataDataScalar( Kokkos::DynRankView<outputDataValueType, outputDataProperties...> outputData,
660 const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
661 const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight,
662 const bool sumInto ) {
663
664#ifdef HAVE_INTREPID2_DEBUG
665 {
666 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.rank() != 2, std::invalid_argument,
667 ">>> ERROR (ArrayTools::contractDataDataScalar): Rank of the left input argument must equal 2!");
668 INTREPID2_TEST_FOR_EXCEPTION( inputDataRight.rank() != 2, std::invalid_argument,
669 ">>> ERROR (ArrayTools::contractDataDataScalar): Rank of right input argument must equal 2!");
670 INTREPID2_TEST_FOR_EXCEPTION( outputData.rank() != 1, std::invalid_argument,
671 ">>> ERROR (ArrayTools::contractDataDataScalar): Rank of output argument must equal 1!");
672 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.extent(0) != inputDataRight.extent(0), std::invalid_argument,
673 ">>> ERROR (ArrayTools::contractDataDataScalar): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
674 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.extent(1) != inputDataRight.extent(1), std::invalid_argument,
675 ">>> ERROR (ArrayTools::contractDataDataScalar): First dimensions (numbers of integration points) of the left and right input containers must agree!");
676 INTREPID2_TEST_FOR_EXCEPTION( outputData.extent(0) != inputDataRight.extent(0), std::invalid_argument,
677 ">>> ERROR (ArrayTools::contractDataDataScalar): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
678 }
679#endif
680
682 inputDataLeft,
683 inputDataRight,
684 sumInto );
685 }
686
687
688 template<typename DeviceType>
689 template<typename outputDataValueType, class ...outputDataProperties,
690 typename inputDataLeftValueType, class ...inputDataLeftProperties,
691 typename inputDataRightValueType, class ...inputDataRightProperties>
692 void
694 contractDataDataVector( Kokkos::DynRankView<outputDataValueType, outputDataProperties...> outputData,
695 const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
696 const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight,
697 const bool sumInto ) {
698
699#ifdef HAVE_INTREPID2_DEBUG
700 {
701 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.rank() != 3, std::invalid_argument,
702 ">>> ERROR (ArrayTools::contractDataDataVector): Rank of the left input argument must equal 3!");
703 INTREPID2_TEST_FOR_EXCEPTION( inputDataRight.rank() != 3, std::invalid_argument,
704 ">>> ERROR (ArrayTools::contractDataDataVector): Rank of right input argument must equal 3!");
705 INTREPID2_TEST_FOR_EXCEPTION( outputData.rank() != 1, std::invalid_argument,
706 ">>> ERROR (ArrayTools::contractDataDataVector): Rank of output argument must equal 1!");
707 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.extent(0) != inputDataRight.extent(0), std::invalid_argument,
708 ">>> ERROR (ArrayTools::contractDataDataVector): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
709 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.extent(1) != inputDataRight.extent(1), std::invalid_argument,
710 ">>> ERROR (ArrayTools::contractDataDataVector): First dimensions (numbers of integration points) of the left and right input containers must agree!");
711 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.extent(2) != inputDataRight.extent(2), std::invalid_argument,
712 ">>> ERROR (ArrayTools::contractDataDataVector): Second dimensions (numbers of vector components) of the left and right input containers must agree!");
713 INTREPID2_TEST_FOR_EXCEPTION( outputData.extent(0) != inputDataRight.extent(0), std::invalid_argument,
714 ">>> ERROR (ArrayTools::contractDataDataVector): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
715 }
716#endif
717
719 inputDataLeft,
720 inputDataRight,
721 sumInto );
722 }
723
724
725 template<typename DeviceType>
726 template<typename outputDataValueType, class ...outputDataProperties,
727 typename inputDataLeftValueType, class ...inputDataLeftProperties,
728 typename inputDataRightValueType, class ...inputDataRightProperties>
729 void
731 contractDataDataTensor( Kokkos::DynRankView<outputDataValueType, outputDataProperties...> outputData,
732 const Kokkos::DynRankView<inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft,
733 const Kokkos::DynRankView<inputDataRightValueType,inputDataRightProperties...> inputDataRight,
734 const bool sumInto ) {
735
736#ifdef HAVE_INTREPID2_DEBUG
737 {
738 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.rank() != 4, std::invalid_argument,
739 ">>> ERROR (ArrayTools::contractDataDataTensor): Rank of the left input argument must equal 4");
740 INTREPID2_TEST_FOR_EXCEPTION( inputDataRight.rank() != 4, std::invalid_argument,
741 ">>> ERROR (ArrayTools::contractDataDataTensor): Rank of right input argument must equal 4!");
742 INTREPID2_TEST_FOR_EXCEPTION( outputData.rank() != 1, std::invalid_argument,
743 ">>> ERROR (ArrayTools::contractDataDataTensor): Rank of output argument must equal 1!");
744 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.extent(0) != inputDataRight.extent(0), std::invalid_argument,
745 ">>> ERROR (ArrayTools::contractDataDataTensor): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
746 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.extent(1) != inputDataRight.extent(1), std::invalid_argument,
747 ">>> ERROR (ArrayTools::contractDataDataTensor): First dimensions (numbers of integration points) of the left and right input containers must agree!");
748 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.extent(2) != inputDataRight.extent(2), std::invalid_argument,
749 ">>> ERROR (ArrayTools::contractDataDataTensor): Second dimensions (first tensor dimensions) of the left and right input containers must agree!");
750 INTREPID2_TEST_FOR_EXCEPTION( inputDataLeft.extent(3) != inputDataRight.extent(3), std::invalid_argument,
751 ">>> ERROR (ArrayTools::contractDataDataTensor): Third dimensions (second tensor dimensions) of the left and right input containers must agree!");
752 INTREPID2_TEST_FOR_EXCEPTION( outputData.extent(0) != inputDataRight.extent(0), std::invalid_argument,
753 ">>> ERROR (ArrayTools::contractDataDataTensor): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
754 }
755#endif
756
758 inputDataLeft,
759 inputDataRight,
760 sumInto );
761 }
762
763}
764#endif
Utility class that provides methods for higher-order algebraic manipulation of user-defined arrays,...
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 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 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 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 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 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 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 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 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...
Functor to contractDataData see Intrepid2::ArrayTools for more.
Functor to contractDataField see Intrepid2::ArrayTools for more.
Functor to contractFieldField see Intrepid2::ArrayTools for more.