52 template<
class Scalar,
class ArrayTypeOut,
class ArrayTypeIn>
54 const ArrayTypeIn & inVals) {
56 ArrayTools::cloneFields<Scalar>(outVals, inVals);
60 template<
class Scalar,
class ArrayTypeOut,
class ArrayTypeJac,
class ArrayTypeIn>
62 const ArrayTypeJac & jacobianInverse,
63 const ArrayTypeIn & inVals,
64 const char transpose) {
66 ArrayTools::matvecProductDataField<Scalar>(outVals, jacobianInverse, inVals, transpose);
70 template<
class Scalar,
class ArrayTypeOut,
class ArrayTypeJac,
class ArrayTypeIn>
72 const ArrayTypeJac & jacobianInverse,
73 const ArrayTypeIn & inVals,
74 const char transpose) {
76 ArrayTools::matvecProductDataField<Scalar>(outVals, jacobianInverse, inVals, transpose);
80 template<
class Scalar,
class ArrayTypeOut,
class ArrayTypeJac,
class ArrayTypeDet,
class ArrayTypeIn>
82 const ArrayTypeJac & jacobian,
83 const ArrayTypeDet & jacobianDet,
84 const ArrayTypeIn & inVals,
85 const char transpose) {
87 ArrayTools::matvecProductDataField<Scalar>(outVals, jacobian, inVals, transpose);
88 ArrayTools::scalarMultiplyDataField<Scalar>(outVals, jacobianDet, outVals,
true);
92 template<
class Scalar,
class ArrayTypeOut,
class ArrayTypeJac,
class ArrayTypeDet,
class ArrayTypeIn>
94 const ArrayTypeJac & jacobian,
95 const ArrayTypeDet & jacobianDet,
96 const ArrayTypeIn & inVals,
97 const char transpose) {
99 ArrayTools::matvecProductDataField<Scalar>(outVals, jacobian, inVals, transpose);
100 ArrayTools::scalarMultiplyDataField<Scalar>(outVals, jacobianDet, outVals,
true);
104 template<
class Scalar,
class ArrayTypeOut,
class ArrayTypeDet,
class ArrayTypeIn>
106 const ArrayTypeDet & jacobianDet,
107 const ArrayTypeIn & inVals) {
109 ArrayTools::scalarMultiplyDataField<Scalar>(outVals, jacobianDet, inVals,
true);
113 template<
class Scalar,
class ArrayTypeOut,
class ArrayTypeDet,
class ArrayTypeIn>
115 const ArrayTypeDet & jacobianDet,
116 const ArrayTypeIn & inVals) {
118 ArrayTools::scalarMultiplyDataField<Scalar>(outVals, jacobianDet, inVals,
true);
121 template<
class Scalar>
125 const ECompEngine compEngine,
126 const bool sumInto) {
127 int outRank = getrank(outputValues);
128 int lRank = getrank(leftValues);
133 #ifdef HAVE_INTREPID_DEBUG
134 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
rank() != 2 ), std::invalid_argument,
135 ">>> ERROR (ArrayTools::contractDataDataScalar): Rank of the left input argument must equal 2!");
136 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
rank() != 2 ), std::invalid_argument,
137 ">>> ERROR (ArrayTools::contractDataDataScalar): Rank of right input argument must equal 2!");
138 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
rank() != 1 ), std::invalid_argument,
139 ">>> ERROR (ArrayTools::contractDataDataScalar): Rank of output argument must equal 1!");
140 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
141 ">>> ERROR (ArrayTools::contractDataDataScalar): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
142 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(1) != rightValues.
dimension(1) ), std::invalid_argument,
143 ">>> ERROR (ArrayTools::contractDataDataScalar): First dimensions (numbers of integration points) of the left and right input containers must agree!");
144 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
145 ">>> ERROR (ArrayTools::contractDataDataScalar): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
155 for (
int cl = 0; cl < numCells; cl++) {
157 for (
int qp = 0; qp < numPoints; qp++) {
158 tmpVal += leftValues(cl, qp)*rightValues(cl, qp);
160 outputValues(cl) += tmpVal;
164 for (
int cl = 0; cl < numCells; cl++) {
166 for (
int qp = 0; qp < numPoints; qp++) {
167 tmpVal += leftValues(cl, qp)*rightValues(cl, qp);
169 outputValues(cl) = tmpVal;
178 for (
int cl=0; cl < numCells; cl++) {
179 Teuchos::BLAS<int, Scalar> myblas;
180 outputValues(cl) += myblas.DOT(numPoints, &leftValues[cl*numPoints], incr, &rightValues[cl*numPoints], incr);
184 for (
int cl=0; cl < numCells; cl++) {
185 Teuchos::BLAS<int, Scalar> myblas;
186 outputValues(cl) = myblas.DOT(numPoints, &leftValues[cl*numPoints], incr, &rightValues[cl*numPoints], incr);
193 TEUCHOS_TEST_FOR_EXCEPTION( ( ~isValidCompEngine(compEngine) ), std::invalid_argument,
194 ">>> ERROR (ArrayTools::contractDataDataScalar): Computational engine not defined!");
199 #ifdef HAVE_INTREPID_DEBUG
200 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
rank() != 3 ), std::invalid_argument,
201 ">>> ERROR (ArrayTools::contractDataDataVector): Rank of the left input argument must equal 3!");
202 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
rank() != 3 ), std::invalid_argument,
203 ">>> ERROR (ArrayTools::contractDataDataVector): Rank of right input argument must equal 3!");
204 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
rank() != 1 ), std::invalid_argument,
205 ">>> ERROR (ArrayTools::contractDataDataVector): Rank of output argument must equal 1!");
206 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
207 ">>> ERROR (ArrayTools::contractDataDataVector): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
208 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(1) != rightValues.
dimension(1) ), std::invalid_argument,
209 ">>> ERROR (ArrayTools::contractDataDataVector): First dimensions (numbers of integration points) of the left and right input containers must agree!");
210 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(2) != rightValues.
dimension(2) ), std::invalid_argument,
211 ">>> ERROR (ArrayTools::contractDataDataVector): Second dimensions (numbers of vector components) of the left and right input containers must agree!");
212 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
213 ">>> ERROR (ArrayTools::contractDataDataVector): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
224 for (
int cl = 0; cl < numCells; cl++) {
226 for (
int qp = 0; qp < numPoints; qp++) {
227 for (
int iVec = 0; iVec < dimVec; iVec++) {
228 tmpVal += leftValues(cl, qp, iVec)*rightValues(cl, qp, iVec);
231 outputValues(cl) += tmpVal;
235 for (
int cl = 0; cl < numCells; cl++) {
237 for (
int qp = 0; qp < numPoints; qp++) {
238 for (
int iVec = 0; iVec < dimVec; iVec++) {
239 tmpVal += leftValues(cl, qp, iVec)*rightValues(cl, qp, iVec);
242 outputValues(cl) = tmpVal;
249 int skip = numPoints*dimVec;
252 for (
int cl=0; cl < numCells; cl++) {
253 Teuchos::BLAS<int, Scalar> myblas;
254 outputValues(cl) += myblas.DOT(skip, &leftValues[cl*skip], incr, &rightValues[cl*skip], incr);
258 for (
int cl=0; cl < numCells; cl++) {
259 Teuchos::BLAS<int, Scalar> myblas;
260 outputValues(cl) = myblas.DOT(skip, &leftValues[cl*skip], incr, &rightValues[cl*skip], incr);
267 TEUCHOS_TEST_FOR_EXCEPTION( ( ~isValidCompEngine(compEngine) ), std::invalid_argument,
268 ">>> ERROR (ArrayTools::contractDataDataVector): Computational engine not defined!");
276 #ifdef HAVE_INTREPID_DEBUG
277 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
rank() != 4 ), std::invalid_argument,
278 ">>> ERROR (ArrayTools::contractDataDataTensor): Rank of the left input argument must equal 4");
279 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
rank() != 4 ), std::invalid_argument,
280 ">>> ERROR (ArrayTools::contractDataDataTensor): Rank of right input argument must equal 4!");
281 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
rank() != 1 ), std::invalid_argument,
282 ">>> ERROR (ArrayTools::contractDataDataTensor): Rank of output argument must equal 1!");
283 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
284 ">>> ERROR (ArrayTools::contractDataDataTensor): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
285 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(1) != rightValues.
dimension(1) ), std::invalid_argument,
286 ">>> ERROR (ArrayTools::contractDataDataTensor): First dimensions (numbers of integration points) of the left and right input containers must agree!");
287 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(2) != rightValues.
dimension(2) ), std::invalid_argument,
288 ">>> ERROR (ArrayTools::contractDataDataTensor): Second dimensions (first tensor dimensions) of the left and right input containers must agree!");
289 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(3) != rightValues.
dimension(3) ), std::invalid_argument,
290 ">>> ERROR (ArrayTools::contractDataDataTensor): Third dimensions (second tensor dimensions) of the left and right input containers must agree!");
291 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
292 ">>> ERROR (ArrayTools::contractDataDataTensor): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
298 int dim1Tensor = leftValues.
dimension(2);
299 int dim2Tensor = leftValues.
dimension(3);
304 for (
int cl = 0; cl < numCells; cl++) {
306 for (
int qp = 0; qp < numPoints; qp++) {
307 for (
int iTens1 = 0; iTens1 < dim1Tensor; iTens1++) {
308 for (
int iTens2 = 0; iTens2 < dim2Tensor; iTens2++) {
309 tmpVal += leftValues(cl, qp, iTens1, iTens2)*rightValues(cl, qp, iTens1, iTens2);
313 outputValues(cl) += tmpVal;
317 for (
int cl = 0; cl < numCells; cl++) {
319 for (
int qp = 0; qp < numPoints; qp++) {
320 for (
int iTens1 = 0; iTens1 < dim1Tensor; iTens1++) {
321 for (
int iTens2 = 0; iTens2 < dim2Tensor; iTens2++) {
322 tmpVal += leftValues(cl, qp, iTens1, iTens2)*rightValues(cl, qp, iTens1, iTens2);
326 outputValues(cl) = tmpVal;
333 int skip = numPoints*dim1Tensor*dim2Tensor;
336 for (
int cl=0; cl < numCells; cl++) {
337 Teuchos::BLAS<int, Scalar> myblas;
338 outputValues(cl) += myblas.DOT(skip, &leftValues[cl*skip], incr, &rightValues[cl*skip], incr);
342 for (
int cl=0; cl < numCells; cl++) {
343 Teuchos::BLAS<int, Scalar> myblas;
344 outputValues(cl) = myblas.DOT(skip, &leftValues[cl*skip], incr, &rightValues[cl*skip], incr);
351 TEUCHOS_TEST_FOR_EXCEPTION( ( ~isValidCompEngine(compEngine) ), std::invalid_argument,
352 ">>> ERROR (ArrayTools::contractDataDataTensor): Computational engine not defined!");
358 #ifdef HAVE_INTREPID_DEBUG
360 TEUCHOS_TEST_FOR_EXCEPTION( ((lRank != 2) && (lRank != 3) && (lRank != 4)), std::invalid_argument,
361 ">>> ERROR (FunctionSpaceTools::dataIntegral): Left data input container must have rank 2, 3 or 4.");
373 #ifdef HAVE_INTREPID_DEBUG
374 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
rank() != 3 ), std::invalid_argument,
375 ">>> ERROR (ArrayTools::contractDataFieldScalar): Rank of the fields input argument must equal 3!");
376 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
rank() != 2 ), std::invalid_argument,
377 ">>> ERROR (ArrayTools::contractDataFieldScalar): Rank of the data input argument must equal 2!");
378 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
rank() != 2 ), std::invalid_argument,
379 ">>> ERROR (ArrayTools::contractDataFieldScalar): Rank of output argument must equal 2!");
380 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
dimension(0) != leftValues.
dimension(0) ), std::invalid_argument,
381 ">>> ERROR (ArrayTools::contractDataFieldScalar): Zeroth dimensions (number of integration domains) of the fields and data input containers must agree!");
382 TEUCHOS_TEST_FOR_EXCEPTION( ( (rightValues.
dimension(2) != leftValues.
dimension(1)) && (leftValues.
dimension(1) != 1) ), std::invalid_argument,
383 ">>> 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!");
384 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
385 ">>> ERROR (ArrayTools::contractDataFieldScalar): Zeroth dimensions (numbers of integration domains) of the fields input and output containers must agree!");
386 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(1) != rightValues.
dimension(1) ), std::invalid_argument,
387 ">>> ERROR (ArrayTools::contractDataFieldScalar): First dimensions (number of fields) of the fields input and output containers must agree!");
392 int numFields = rightValues.
dimension(1);
393 int numPoints = rightValues.
dimension(2);
394 int numDataPoints = leftValues.
dimension(1);
396 ECompEngine myCompEngine = (numDataPoints == 1 ? COMP_CPP : compEngine);
398 switch(myCompEngine) {
401 if (numDataPoints != 1) {
402 for (
int cl = 0; cl < numCells; cl++) {
403 for (
int lbf = 0; lbf < numFields; lbf++) {
405 for (
int qp = 0; qp < numPoints; qp++) {
406 tmpVal += rightValues(cl, lbf, qp)*leftValues(cl, qp);
408 outputValues(cl, lbf) += tmpVal;
413 for (
int cl = 0; cl < numCells; cl++) {
414 for (
int lbf = 0; lbf < numFields; lbf++) {
416 for (
int qp = 0; qp < numPoints; qp++) {
417 tmpVal += rightValues(cl, lbf, qp)*leftValues(cl, 0);
419 outputValues(cl, lbf) += tmpVal;
425 if (numDataPoints != 1) {
426 for (
int cl = 0; cl < numCells; cl++) {
427 for (
int lbf = 0; lbf < numFields; lbf++) {
429 for (
int qp = 0; qp < numPoints; qp++) {
430 tmpVal += rightValues(cl, lbf, qp)*leftValues(cl, qp);
432 outputValues(cl, lbf) = tmpVal;
437 for (
int cl = 0; cl < numCells; cl++) {
438 for (
int lbf = 0; lbf < numFields; lbf++) {
440 for (
int qp = 0; qp < numPoints; qp++) {
441 tmpVal += rightValues(cl, lbf, qp)*leftValues(cl, 0);
443 outputValues(cl, lbf) = tmpVal;
472 int numData = numPoints;
473 int skipL = numFields*numPoints;
474 int skipR = numPoints;
475 int skipOp = numFields;
482 for (
int cl=0; cl < numCells; cl++) {
484 Teuchos::BLAS<int, Scalar> myblas;
485 myblas.GEMM(Teuchos::TRANS, Teuchos::NO_TRANS,
486 1, numFields, numData,
487 alpha, &leftValues[cl*skipR], numData,
488 &rightValues[cl*skipL], numData,
489 beta, &outputValues[cl*skipOp], 1);
503 TEUCHOS_TEST_FOR_EXCEPTION( ( ~isValidCompEngine(compEngine) ), std::invalid_argument,
504 ">>> ERROR (ArrayTools::contractDataFieldScalar): Computational engine not defined!");
509 #ifdef HAVE_INTREPID_DEBUG
510 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
rank() != 4 ), std::invalid_argument,
511 ">>> ERROR (ArrayTools::contractDataFieldVector): Rank of the fields input argument must equal 4!");
512 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
rank() != 3 ), std::invalid_argument,
513 ">>> ERROR (ArrayTools::contractDataFieldVector): Rank of the data input argument must equal 3!");
514 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
rank() != 2 ), std::invalid_argument,
515 ">>> ERROR (ArrayTools::contractDataFieldVector): Rank of output argument must equal 2!");
516 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
dimension(0) != leftValues.
dimension(0) ), std::invalid_argument,
517 ">>> ERROR (ArrayTools::contractDataFieldVector): Zeroth dimensions (number of integration domains) of the fields and data input containers must agree!");
518 TEUCHOS_TEST_FOR_EXCEPTION( ( (rightValues.
dimension(2) != leftValues.
dimension(1)) && (leftValues.
dimension(1) != 1) ), std::invalid_argument,
519 ">>> 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!");
520 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
dimension(3) != leftValues.
dimension(2) ), std::invalid_argument,
521 ">>> ERROR (ArrayTools::contractDataFieldVector): Third dimension of the fields input container and second dimension of data input container (vector index) must agree!");
522 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
523 ">>> ERROR (ArrayTools::contractDataFieldVector): Zeroth dimensions (numbers of integration domains) of the fields input and output containers must agree!");
524 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(1) != rightValues.
dimension(1) ), std::invalid_argument,
525 ">>> ERROR (ArrayTools::contractDataFieldVector): First dimensions of output container and fields input container (number of fields) must agree!");
530 int numFields = rightValues.
dimension(1);
531 int numPoints = rightValues.
dimension(2);
533 int numDataPoints = leftValues.
dimension(1);
535 ECompEngine myCompEngine = (numDataPoints == 1 ? COMP_CPP : compEngine);
537 switch(myCompEngine) {
540 if (numDataPoints != 1) {
541 for (
int cl = 0; cl < numCells; cl++) {
542 for (
int lbf = 0; lbf < numFields; lbf++) {
544 for (
int qp = 0; qp < numPoints; qp++) {
545 for (
int iVec = 0; iVec < dimVec; iVec++) {
546 tmpVal += rightValues(cl, lbf, qp, iVec)*leftValues(cl, qp, iVec);
549 outputValues(cl, lbf) += tmpVal;
554 for (
int cl = 0; cl < numCells; cl++) {
555 for (
int lbf = 0; lbf < numFields; lbf++) {
557 for (
int qp = 0; qp < numPoints; qp++) {
558 for (
int iVec = 0; iVec < dimVec; iVec++) {
559 tmpVal += rightValues(cl, lbf, qp, iVec)*leftValues(cl, 0, iVec);
562 outputValues(cl, lbf) += tmpVal;
568 if (numDataPoints != 1) {
569 for (
int cl = 0; cl < numCells; cl++) {
570 for (
int lbf = 0; lbf < numFields; lbf++) {
572 for (
int qp = 0; qp < numPoints; qp++) {
573 for (
int iVec = 0; iVec < dimVec; iVec++) {
574 tmpVal += rightValues(cl, lbf, qp, iVec)*leftValues(cl, qp, iVec);
577 outputValues(cl, lbf) = tmpVal;
582 for (
int cl = 0; cl < numCells; cl++) {
583 for (
int lbf = 0; lbf < numFields; lbf++) {
585 for (
int qp = 0; qp < numPoints; qp++) {
586 for (
int iVec = 0; iVec < dimVec; iVec++) {
587 tmpVal += rightValues(cl, lbf, qp, iVec)*leftValues(cl, 0, iVec);
590 outputValues(cl, lbf) = tmpVal;
619 int numData = numPoints*dimVec;
620 int skipL = numFields*numData;
622 int skipOp = numFields;
629 for (
int cl=0; cl < numCells; cl++) {
631 Teuchos::BLAS<int, Scalar> myblas;
632 myblas.GEMM(Teuchos::TRANS, Teuchos::NO_TRANS,
633 1, numFields, numData,
634 alpha, &leftValues[cl*skipR], numData,
635 &rightValues[cl*skipL], numData,
636 beta, &outputValues[cl*skipOp], 1);
650 TEUCHOS_TEST_FOR_EXCEPTION( ( ~isValidCompEngine(compEngine) ), std::invalid_argument,
651 ">>> ERROR (ArrayTools::contractDataFieldVector): Computational engine not defined!");
656 #ifdef HAVE_INTREPID_DEBUG
657 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
rank() != 5 ), std::invalid_argument,
658 ">>> ERROR (ArrayTools::contractDataFieldTensor): Rank of the fields input argument must equal 5!");
659 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
rank() != 4 ), std::invalid_argument,
660 ">>> ERROR (ArrayTools::contractDataFieldTensor): Rank of the data input argument must equal 4!");
661 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
rank() != 2 ), std::invalid_argument,
662 ">>> ERROR (ArrayTools::contractDataFieldTensor): Rank of output argument must equal 2!");
663 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
dimension(0) != leftValues.
dimension(0) ), std::invalid_argument,
664 ">>> ERROR (ArrayTools::contractDataFieldTensor): Zeroth dimensions (number of integration domains) of the fields and data input containers must agree!");
665 TEUCHOS_TEST_FOR_EXCEPTION( ( (rightValues.
dimension(2) != leftValues.
dimension(1)) && (leftValues.
dimension(1) != 1) ), std::invalid_argument,
666 ">>> 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!");
667 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
dimension(3) != leftValues.
dimension(2) ), std::invalid_argument,
668 ">>> ERROR (ArrayTools::contractDataFieldTensor): Third dimension of the fields input container and second dimension of data input container (first tensor dimension) must agree!");
669 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
dimension(4) != leftValues.
dimension(3) ), std::invalid_argument,
670 ">>> ERROR (ArrayTools::contractDataFieldTensor): Fourth dimension of the fields input container and third dimension of data input container (second tensor dimension) must agree!");
671 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
672 ">>> ERROR (ArrayTools::contractDataFieldTensor): Zeroth dimensions (numbers of integration domains) of the fields input and output containers must agree!");
673 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(1) != rightValues.
dimension(1) ), std::invalid_argument,
674 ">>> ERROR (ArrayTools::contractDataFieldTensor): First dimensions (number of fields) of output container and fields input container must agree!");
679 int numFields = rightValues.
dimension(1);
680 int numPoints = rightValues.
dimension(2);
683 int numDataPoints = leftValues.
dimension(1);
685 ECompEngine myCompEngine = (numDataPoints == 1 ? COMP_CPP : compEngine);
687 switch(myCompEngine) {
690 if (numDataPoints != 1) {
691 for (
int cl = 0; cl < numCells; cl++) {
692 for (
int lbf = 0; lbf < numFields; lbf++) {
694 for (
int qp = 0; qp < numPoints; qp++) {
695 for (
int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
696 for (
int iTens2 =0; iTens2 < dim2Tens; iTens2++) {
697 tmpVal += rightValues(cl, lbf, qp, iTens1, iTens2)*leftValues(cl, qp, iTens1, iTens2);
701 outputValues(cl, lbf) += tmpVal;
706 for (
int cl = 0; cl < numCells; cl++) {
707 for (
int lbf = 0; lbf < numFields; lbf++) {
709 for (
int qp = 0; qp < numPoints; qp++) {
710 for (
int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
711 for (
int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
712 tmpVal += rightValues(cl, lbf, qp, iTens1, iTens2)*leftValues(cl, 0, iTens1, iTens2);
716 outputValues(cl, lbf) += tmpVal;
722 if (numDataPoints != 1) {
723 for (
int cl = 0; cl < numCells; cl++) {
724 for (
int lbf = 0; lbf < numFields; lbf++) {
726 for (
int qp = 0; qp < numPoints; qp++) {
727 for (
int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
728 for (
int iTens2 =0; iTens2 < dim2Tens; iTens2++) {
729 tmpVal += rightValues(cl, lbf, qp, iTens1, iTens2)*leftValues(cl, qp, iTens1, iTens2);
733 outputValues(cl, lbf) = tmpVal;
738 for (
int cl = 0; cl < numCells; cl++) {
739 for (
int lbf = 0; lbf < numFields; lbf++) {
741 for (
int qp = 0; qp < numPoints; qp++) {
742 for (
int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
743 for (
int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
744 tmpVal += rightValues(cl, lbf, qp, iTens1, iTens2)*leftValues(cl, 0, iTens1, iTens2);
748 outputValues(cl, lbf) = tmpVal;
777 int numData = numPoints*dim1Tens*dim2Tens;
778 int skipL = numFields*numData;
780 int skipOp = numFields;
787 for (
int cl=0; cl < numCells; cl++) {
789 Teuchos::BLAS<int, Scalar> myblas;
790 myblas.GEMM(Teuchos::TRANS, Teuchos::NO_TRANS,
791 1, numFields, numData,
792 alpha, &leftValues[cl*skipR], numData,
793 &rightValues[cl*skipL], numData,
794 beta, &outputValues[cl*skipOp], 1);
808 TEUCHOS_TEST_FOR_EXCEPTION( ( ~isValidCompEngine(compEngine) ), std::invalid_argument,
809 ">>> ERROR (ArrayTools::contractDataFieldTensor): Computational engine not defined!");
814 #ifdef HAVE_INTREPID_DEBUG
816 TEUCHOS_TEST_FOR_EXCEPTION( ((lRank != 2) && (lRank != 3) && (lRank != 4)), std::invalid_argument,
817 ">>> ERROR (FunctionSpaceTools::functionalIntegral): Data input container must have rank 2, 3 or 4.");
828 #ifdef HAVE_INTREPID_DEBUG
829 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
rank() != 3 ), std::invalid_argument,
830 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Rank of the left input argument must equal 3!");
831 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
rank() != 3 ), std::invalid_argument,
832 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Rank of right input argument must equal 3!");
833 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
rank() != 3 ), std::invalid_argument,
834 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Rank of output argument must equal 3!");
835 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
836 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
837 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(2) != rightValues.
dimension(2) ), std::invalid_argument,
838 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Second dimensions (numbers of integration points) of the left and right input containers must agree!");
839 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
840 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
841 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(1) != leftValues.
dimension(1) ), std::invalid_argument,
842 ">>> ERROR (ArrayTools::contractFieldFieldScalar): First dimension of output container and first dimension of left input container must agree!");
843 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(2) != rightValues.
dimension(1) ), std::invalid_argument,
844 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Second dimension of output container and first dimension of right input container must agree!");
849 int numLeftFields = leftValues.
dimension(1);
850 int numRightFields = rightValues.
dimension(1);
856 for (
int cl = 0; cl < numCells; cl++) {
857 for (
int lbf = 0; lbf < numLeftFields; lbf++) {
858 for (
int rbf = 0; rbf < numRightFields; rbf++) {
860 for (
int qp = 0; qp < numPoints; qp++) {
861 tmpVal += leftValues(cl, lbf, qp)*rightValues(cl, rbf, qp);
863 outputValues(cl, lbf, rbf) += tmpVal;
869 for (
int cl = 0; cl < numCells; cl++) {
870 for (
int lbf = 0; lbf < numLeftFields; lbf++) {
871 for (
int rbf = 0; rbf < numRightFields; rbf++) {
873 for (
int qp = 0; qp < numPoints; qp++) {
874 tmpVal += leftValues(cl, lbf, qp)*rightValues(cl, rbf, qp);
876 outputValues(cl, lbf, rbf) = tmpVal;
905 int skipL = numLeftFields*numPoints;
906 int skipR = numRightFields*numPoints;
907 int skipOp = numLeftFields*numRightFields;
914 for (
int cl=0; cl < numCells; cl++) {
916 Teuchos::BLAS<int, Scalar> myblas;
917 myblas.GEMM(Teuchos::TRANS, Teuchos::NO_TRANS,
918 numRightFields, numLeftFields, numPoints,
919 alpha, &rightValues[cl*skipR], numPoints,
920 &leftValues[cl*skipL], numPoints,
921 beta, &outputValues[cl*skipOp], numRightFields);
935 TEUCHOS_TEST_FOR_EXCEPTION( ( ~isValidCompEngine(compEngine) ), std::invalid_argument,
936 ">>> ERROR (ArrayTools::contractFieldFieldScalar): Computational engine not defined!");
942 #ifdef HAVE_INTREPID_DEBUG
943 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
rank() != 4 ), std::invalid_argument,
944 ">>> ERROR (ArrayTools::contractFieldFieldVector): Rank of the left input argument must equal 4!");
945 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
rank() != 4 ), std::invalid_argument,
946 ">>> ERROR (ArrayTools::contractFieldFieldVector): Rank of right input argument must equal 4!");
947 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
rank() != 3 ), std::invalid_argument,
948 ">>> ERROR (ArrayTools::contractFieldFieldVector): Rank of output argument must equal 3!");
949 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
950 ">>> ERROR (ArrayTools::contractFieldFieldVector): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
951 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(2) != rightValues.
dimension(2) ), std::invalid_argument,
952 ">>> ERROR (ArrayTools::contractFieldFieldVector): Second dimensions (numbers of integration points) of the left and right input containers must agree!");
953 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(3) != rightValues.
dimension(3) ), std::invalid_argument,
954 ">>> ERROR (ArrayTools::contractFieldFieldVector): Third dimensions (numbers of vector components) of the left and right input containers must agree!");
955 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
956 ">>> ERROR (ArrayTools::contractFieldFieldVector): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
957 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(1) != leftValues.
dimension(1) ), std::invalid_argument,
958 ">>> ERROR (ArrayTools::contractFieldFieldVector): First dimension of output container and first dimension of left input container must agree!");
959 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(2) != rightValues.
dimension(1) ), std::invalid_argument,
960 ">>> ERROR (ArrayTools::contractFieldFieldVector): Second dimension of output container and first dimension of right input container must agree!");
965 int numLeftFields = leftValues.
dimension(1);
966 int numRightFields = rightValues.
dimension(1);
973 for (
int cl = 0; cl < numCells; cl++) {
974 for (
int lbf = 0; lbf < numLeftFields; lbf++) {
975 for (
int rbf = 0; rbf < numRightFields; rbf++) {
977 for (
int qp = 0; qp < numPoints; qp++) {
978 for (
int iVec = 0; iVec < dimVec; iVec++) {
979 tmpVal += leftValues(cl, lbf, qp, iVec)*rightValues(cl, rbf, qp, iVec);
982 outputValues(cl, lbf, rbf) += tmpVal;
988 for (
int cl = 0; cl < numCells; cl++) {
989 for (
int lbf = 0; lbf < numLeftFields; lbf++) {
990 for (
int rbf = 0; rbf < numRightFields; rbf++) {
992 for (
int qp = 0; qp < numPoints; qp++) {
993 for (
int iVec = 0; iVec < dimVec; iVec++) {
994 tmpVal += leftValues(cl, lbf, qp, iVec)*rightValues(cl, rbf, qp, iVec);
997 outputValues(cl, lbf, rbf) = tmpVal;
1026 int numData = numPoints*dimVec;
1027 int skipL = numLeftFields*numData;
1028 int skipR = numRightFields*numData;
1029 int skipOp = numLeftFields*numRightFields;
1036 for (
int cl=0; cl < numCells; cl++) {
1038 Teuchos::BLAS<int, Scalar> myblas;
1039 myblas.GEMM(Teuchos::TRANS, Teuchos::NO_TRANS,
1040 numRightFields, numLeftFields, numData,
1041 alpha, &rightValues[cl*skipR], numData,
1042 &leftValues[cl*skipL], numData,
1043 beta, &outputValues[cl*skipOp], numRightFields);
1057 TEUCHOS_TEST_FOR_EXCEPTION( ( ~isValidCompEngine(compEngine) ), std::invalid_argument,
1058 ">>> ERROR (ArrayTools::contractFieldFieldVector): Computational engine not defined!");
1063 #ifdef HAVE_INTREPID_DEBUG
1064 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
rank() != 5 ), std::invalid_argument,
1065 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Rank of the left input argument must equal 5!");
1066 TEUCHOS_TEST_FOR_EXCEPTION( (rightValues.
rank() != 5 ), std::invalid_argument,
1067 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Rank of right input argument must equal 5!");
1068 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
rank() != 3 ), std::invalid_argument,
1069 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Rank of output argument must equal 3!");
1070 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
1071 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Zeroth dimensions (number of integration domains) of the left and right input containers must agree!");
1072 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(2) != rightValues.
dimension(2) ), std::invalid_argument,
1073 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Second dimensions (numbers of integration points) of the left and right input containers must agree!");
1074 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(3) != rightValues.
dimension(3) ), std::invalid_argument,
1075 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Third dimensions (first tensor dimensions) of the left and right input containers must agree!");
1076 TEUCHOS_TEST_FOR_EXCEPTION( (leftValues.
dimension(4) != rightValues.
dimension(4) ), std::invalid_argument,
1077 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Fourth dimensions (second tensor dimensions) of the left and right input containers must agree!");
1078 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(0) != rightValues.
dimension(0) ), std::invalid_argument,
1079 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Zeroth dimensions (numbers of integration domains) of the input and output containers must agree!");
1080 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(1) != leftValues.
dimension(1) ), std::invalid_argument,
1081 ">>> ERROR (ArrayTools::contractFieldFieldTensor): First dimension of output container and first dimension of left input container must agree!");
1082 TEUCHOS_TEST_FOR_EXCEPTION( (outputValues.
dimension(2) != rightValues.
dimension(1) ), std::invalid_argument,
1083 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Second dimension of output container and first dimension of right input container must agree!");
1088 int numLeftFields = leftValues.
dimension(1);
1089 int numRightFields = rightValues.
dimension(1);
1090 int numPoints = leftValues.
dimension(2);
1091 int dim1Tensor = leftValues.
dimension(3);
1092 int dim2Tensor = leftValues.
dimension(4);
1094 switch(compEngine) {
1097 for (
int cl = 0; cl < numCells; cl++) {
1098 for (
int lbf = 0; lbf < numLeftFields; lbf++) {
1099 for (
int rbf = 0; rbf < numRightFields; rbf++) {
1101 for (
int qp = 0; qp < numPoints; qp++) {
1102 for (
int iTens1 = 0; iTens1 < dim1Tensor; iTens1++) {
1103 for (
int iTens2 = 0; iTens2 < dim2Tensor; iTens2++) {
1104 tmpVal += leftValues(cl, lbf, qp, iTens1, iTens2)*rightValues(cl, rbf, qp, iTens1, iTens2);
1108 outputValues(cl, lbf, rbf) += tmpVal;
1114 for (
int cl = 0; cl < numCells; cl++) {
1115 for (
int lbf = 0; lbf < numLeftFields; lbf++) {
1116 for (
int rbf = 0; rbf < numRightFields; rbf++) {
1118 for (
int qp = 0; qp < numPoints; qp++) {
1119 for (
int iTens1 = 0; iTens1 < dim1Tensor; iTens1++) {
1120 for (
int iTens2 = 0; iTens2 < dim2Tensor; iTens2++) {
1121 tmpVal += leftValues(cl, lbf, qp, iTens1, iTens2)*rightValues(cl, rbf, qp, iTens1, iTens2);
1125 outputValues(cl, lbf, rbf) = tmpVal;
1154 int numData = numPoints*dim1Tensor*dim2Tensor;
1155 int skipL = numLeftFields*numData;
1156 int skipR = numRightFields*numData;
1157 int skipOp = numLeftFields*numRightFields;
1164 for (
int cl=0; cl < numCells; cl++) {
1166 Teuchos::BLAS<int, Scalar> myblas;
1167 myblas.GEMM(Teuchos::TRANS, Teuchos::NO_TRANS,
1168 numRightFields, numLeftFields, numData,
1169 alpha, &rightValues[cl*skipR], numData,
1170 &leftValues[cl*skipL], numData,
1171 beta, &outputValues[cl*skipOp], numRightFields);
1185 TEUCHOS_TEST_FOR_EXCEPTION( ( ~isValidCompEngine(compEngine) ), std::invalid_argument,
1186 ">>> ERROR (ArrayTools::contractFieldFieldTensor): Computational engine not defined!");
1191 #ifdef HAVE_INTREPID_DEBUG
1193 TEUCHOS_TEST_FOR_EXCEPTION( ((lRank != 3) && (lRank != 4) && (lRank != 5)), std::invalid_argument,
1194 ">>> ERROR (FunctionSpaceTools::operatorIntegral): Left fields input container must have rank 3, 4 or 5.");
1205 #ifdef HAVE_INTREPID_DEBUG
1207 TEUCHOS_TEST_FOR_EXCEPTION( ((outRank != 1) && (outRank != 2) && (outRank != 3)), std::invalid_argument,
1208 ">>> ERROR (FunctionSpaceTools::integrate): Output container must have rank 1, 2 or 3.");
1214 template<
class Scalar,
class ArrayOut,
class ArrayInLeft,
class ArrayInRight>
1216 const ArrayInLeft & leftValues,
1217 const ArrayInRight & rightValues,
1218 const ECompEngine compEngine,
1219 const bool sumInto) {
1220 ArrayWrapper<Scalar,ArrayOut, Rank<ArrayOut >::value,
false>outputValuesWrap(outputValues);
1221 ArrayWrapper<Scalar,ArrayInLeft, Rank<ArrayInLeft >::value,
true>leftValuesWrap(leftValues);
1222 ArrayWrapper<Scalar,ArrayInRight, Rank<ArrayInRight >::value,
true>rightValuesWrap(rightValues);
1223 int outRank = getrank(outputValues);
1226 dataIntegral<Scalar>(outputValuesWrap, leftValuesWrap, rightValuesWrap, compEngine, sumInto);
1229 functionalIntegral<Scalar>(outputValuesWrap, leftValuesWrap, rightValuesWrap, compEngine, sumInto);
1232 operatorIntegral<Scalar>(outputValuesWrap, leftValuesWrap, rightValuesWrap, compEngine, sumInto);
1235 #ifdef HAVE_INTREPID_DEBUG
1237 TEUCHOS_TEST_FOR_EXCEPTION( ((outRank != 1) && (outRank != 2) && (outRank != 3)), std::invalid_argument,
1238 ">>> ERROR (FunctionSpaceTools::integrate): Output container must have rank 1, 2 or 3.");
1245 template<
class Scalar,
class ArrayOutFields,
class ArrayInFieldsLeft,
class ArrayInFieldsRight>
1247 const ArrayInFieldsLeft & leftFields,
1248 const ArrayInFieldsRight & rightFields,
1249 const ECompEngine compEngine,
1250 const bool sumInto) {
1251 int lRank = getrank(leftFields);
1255 ArrayTools::contractFieldFieldScalar<Scalar>(outputFields, leftFields, rightFields, compEngine, sumInto);
1258 ArrayTools::contractFieldFieldVector<Scalar>(outputFields, leftFields, rightFields, compEngine, sumInto);
1261 ArrayTools::contractFieldFieldTensor<Scalar>(outputFields, leftFields, rightFields, compEngine, sumInto);
1264 #ifdef HAVE_INTREPID_DEBUG
1266 TEUCHOS_TEST_FOR_EXCEPTION( ((lRank != 3) && (lRank != 4) && (lRank != 5)), std::invalid_argument,
1267 ">>> ERROR (FunctionSpaceTools::operatorIntegral): Left fields input container must have rank 3, 4 or 5.");
1276 template<
class Scalar,
class ArrayOutFields,
class ArrayInData,
class ArrayInFields>
1278 const ArrayInData & inputData,
1279 const ArrayInFields & inputFields,
1280 const ECompEngine compEngine,
1281 const bool sumInto) {
1282 int dRank = getrank(inputData);
1286 ArrayTools::contractDataFieldScalar<Scalar>(outputFields, inputData, inputFields, compEngine, sumInto);
1289 ArrayTools::contractDataFieldVector<Scalar>(outputFields, inputData, inputFields, compEngine, sumInto);
1292 ArrayTools::contractDataFieldTensor<Scalar>(outputFields, inputData, inputFields, compEngine, sumInto);
1295 #ifdef HAVE_INTREPID_DEBUG
1297 TEUCHOS_TEST_FOR_EXCEPTION( ((dRank != 2) && (dRank != 3) && (dRank != 4)), std::invalid_argument,
1298 ">>> ERROR (FunctionSpaceTools::functionalIntegral): Data input container must have rank 2, 3 or 4.");
1306 template<
class Scalar,
class ArrayOutData,
class ArrayInDataLeft,
class ArrayInDataRight>
1308 const ArrayInDataLeft & inputDataLeft,
1309 const ArrayInDataRight & inputDataRight,
1310 const ECompEngine compEngine,
1311 const bool sumInto) {
1312 int lRank = getrank(inputDataLeft);
1316 ArrayTools::contractDataDataScalar<Scalar>(outputData, inputDataLeft, inputDataRight, compEngine, sumInto);
1319 ArrayTools::contractDataDataVector<Scalar>(outputData, inputDataLeft, inputDataRight, compEngine, sumInto);
1322 ArrayTools::contractDataDataTensor<Scalar>(outputData, inputDataLeft, inputDataRight, compEngine, sumInto);
1325 #ifdef HAVE_INTREPID_DEBUG
1327 TEUCHOS_TEST_FOR_EXCEPTION( ((lRank != 2) && (lRank != 3) && (lRank != 4)), std::invalid_argument,
1328 ">>> ERROR (FunctionSpaceTools::dataIntegral): Left data input container must have rank 2, 3 or 4.");
1336 template<
class Scalar,
class ArrayOut,
class ArrayDet,
class ArrayWeights>
1338 const ArrayDet & inDet,
1339 const ArrayWeights & inWeights) {
1340 #ifdef HAVE_INTREPID_DEBUG
1342 TEUCHOS_TEST_FOR_EXCEPTION( (inDet.rank() != 2), std::invalid_argument,
1343 ">>> ERROR (FunctionSpaceTools::computeCellMeasure): Input determinants container must have rank 2.");
1347 ArrayTools::scalarMultiplyDataData<Scalar>(outVals, inDet, inWeights);
1349 for (
int cell=0; cell<outVals.dimension(0); cell++) {
1350 if (inDet(cell,0) < 0.0) {
1351 for (
int point=0; point<outVals.dimension(1); point++) {
1352 outVals(cell, point) *= -1.0;
1359 template<
class Scalar,
class ArrayOut,
class ArrayJac,
class ArrayWeights>
1361 const ArrayJac & inJac,
1362 const ArrayWeights & inWeights,
1363 const int whichFace,
1364 const shards::CellTopology & parentCell) {
1366 #ifdef HAVE_INTREPID_DEBUG
1368 TEUCHOS_TEST_FOR_EXCEPTION( (inJac.rank() != 4), std::invalid_argument,
1369 ">>> ERROR (FunctionSpaceTools::computeFaceMeasure): Input Jacobian container must have rank 4.");
1383 ArrayTools::scalarMultiplyDataData<Scalar>(outVals, outVals, inWeights);
1388 template<
class Scalar,
class ArrayOut,
class ArrayJac,
class ArrayWeights>
1390 const ArrayJac & inJac,
1391 const ArrayWeights & inWeights,
1392 const int whichEdge,
1393 const shards::CellTopology & parentCell) {
1395 #ifdef HAVE_INTREPID_DEBUG
1397 TEUCHOS_TEST_FOR_EXCEPTION( (inJac.rank() != 4), std::invalid_argument,
1398 ">>> ERROR (FunctionSpaceTools::computeEdgeMeasure): Input Jacobian container must have rank 4.");
1412 ArrayTools::scalarMultiplyDataData<Scalar>(outVals, outVals, inWeights);
1418 template<
class Scalar,
class ArrayTypeOut,
class ArrayTypeMeasure,
class ArrayTypeIn>
1420 const ArrayTypeMeasure & inMeasure,
1421 const ArrayTypeIn & inVals) {
1423 ArrayTools::scalarMultiplyDataField<Scalar>(outVals, inMeasure, inVals);
1428 template<
class Scalar,
class ArrayOutFields,
class ArrayInData,
class ArrayInFields>
1430 ArrayInData & inputData,
1431 ArrayInFields & inputFields,
1432 const bool reciprocal) {
1434 ArrayTools::scalarMultiplyDataField<Scalar>(outputFields, inputData, inputFields, reciprocal);
1439 template<
class Scalar,
class ArrayOutData,
class ArrayInDataLeft,
class ArrayInDataRight>
1441 ArrayInDataLeft & inputDataLeft,
1442 ArrayInDataRight & inputDataRight,
1443 const bool reciprocal) {
1445 ArrayTools::scalarMultiplyDataData<Scalar>(outputData, inputDataLeft, inputDataRight, reciprocal);
1450 template<
class Scalar,
class ArrayOutFields,
class ArrayInData,
class ArrayInFields>
1452 const ArrayInData & inputData,
1453 const ArrayInFields & inputFields) {
1455 ArrayTools::dotMultiplyDataField<Scalar>(outputFields, inputData, inputFields);
1460 template<
class Scalar,
class ArrayOutData,
class ArrayInDataLeft,
class ArrayInDataRight>
1462 const ArrayInDataLeft & inputDataLeft,
1463 const ArrayInDataRight & inputDataRight) {
1465 ArrayTools::dotMultiplyDataData<Scalar>(outputData, inputDataLeft, inputDataRight);
1470 template<
class Scalar,
class ArrayOutFields,
class ArrayInData,
class ArrayInFields>
1472 const ArrayInData & inputData,
1473 const ArrayInFields & inputFields) {
1475 int outRank = getrank(outputFields);
1480 ArrayTools::crossProductDataField<Scalar>(outputFields, inputData, inputFields);
1483 ArrayTools::outerProductDataField<Scalar>(outputFields, inputData, inputFields);
1486 TEUCHOS_TEST_FOR_EXCEPTION( ((outRank != 3) && (outRank != 4) && (outRank != 5)), std::invalid_argument,
1487 ">>> ERROR (FunctionSpaceTools::vectorMultiplyDataField): Output container must have rank 3, 4 or 5.");
1493 template<
class Scalar,
class ArrayOutData,
class ArrayInDataLeft,
class ArrayInDataRight>
1495 const ArrayInDataLeft & inputDataLeft,
1496 const ArrayInDataRight & inputDataRight) {
1498 int outRank = getrank(outputData);
1503 ArrayTools::crossProductDataData<Scalar>(outputData, inputDataLeft, inputDataRight);
1506 ArrayTools::outerProductDataData<Scalar>(outputData, inputDataLeft, inputDataRight);
1509 TEUCHOS_TEST_FOR_EXCEPTION( ((outRank != 2) && (outRank != 3) && (outRank != 4)), std::invalid_argument,
1510 ">>> ERROR (FunctionSpaceTools::vectorMultiplyDataData): Output container must have rank 2, 3 or 4.");
1516 template<
class Scalar,
class ArrayOutFields,
class ArrayInData,
class ArrayInFields>
1518 const ArrayInData & inputData,
1519 const ArrayInFields & inputFields,
1520 const char transpose) {
1522 int outRank = outputFields.rank();
1526 ArrayTools::matvecProductDataField<Scalar>(outputFields, inputData, inputFields, transpose);
1529 ArrayTools::matmatProductDataField<Scalar>(outputFields, inputData, inputFields, transpose);
1532 TEUCHOS_TEST_FOR_EXCEPTION( ((outRank != 4) && (outRank != 5)), std::invalid_argument,
1533 ">>> ERROR (FunctionSpaceTools::tensorMultiplyDataField): Output container must have rank 4 or 5.");
1538 template<
class Scalar,
class ArrayOutData,
class ArrayInDataLeft,
class ArrayInDataRight>
1541 const ArrayInDataLeft & inputDataLeft,
1542 const ArrayInDataRight & inputDataRight,
1543 const char transpose) {
1544 int outRank = getrank(outputData);
1548 ArrayTools::matvecProductDataData<Scalar>(outputData, inputDataLeft, inputDataRight, transpose);
1551 ArrayTools::matmatProductDataData<Scalar>(outputData, inputDataLeft, inputDataRight, transpose);
1556 template<
class Scalar,
class ArrayOutData,
class ArrayInDataLeft,
class ArrayInDataRight>
1559 const ArrayInDataLeft & inputDataLeft,
1560 const ArrayInDataRight & inputDataRight,
1561 const char transpose) {
1562 ArrayTools::matvecProductDataData<Scalar>(outputData, inputDataLeft, inputDataRight, transpose);
1566 template<
class Scalar,
class ArrayOutData,
class ArrayInDataLeft,
class ArrayInDataRight>
1569 const ArrayInDataLeft & inputDataLeft,
1570 const ArrayInDataRight & inputDataRight,
1571 const char transpose) {
1572 ArrayTools::matmatProductDataData<Scalar>(outputData, inputDataLeft, inputDataRight, transpose);
1576 template<
class Scalar,
class ArrayOutData,
class ArrayInDataLeft,
class ArrayInDataRight>
1578 const ArrayInDataLeft & inputDataLeft,
1579 const ArrayInDataRight & inputDataRight,
1580 const char transpose){
1584 template<
class Scalar,
class ArrayTypeInOut,
class ArrayTypeSign>
1586 const ArrayTypeSign & fieldSigns) {
1587 #ifdef HAVE_INTREPID_DEBUG
1588 TEUCHOS_TEST_FOR_EXCEPTION( (inoutOperator.rank() != 3), std::invalid_argument,
1589 ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): Input operator container must have rank 3.");
1590 TEUCHOS_TEST_FOR_EXCEPTION( (fieldSigns.rank() != 2), std::invalid_argument,
1591 ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): Input field signs container must have rank 2.");
1592 TEUCHOS_TEST_FOR_EXCEPTION( (inoutOperator.dimension(0) != fieldSigns.dimension(0) ), std::invalid_argument,
1593 ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): Zeroth dimensions (number of cells) of the operator and field signs containers must agree!");
1594 TEUCHOS_TEST_FOR_EXCEPTION( (inoutOperator.dimension(1) != fieldSigns.dimension(1) ), std::invalid_argument,
1595 ">>> ERROR (FunctionSpaceTools::applyLeftFieldSigns): First dimensions (number of left fields) of the operator and field signs containers must agree!");
1598 for (
int cell=0; cell<inoutOperator.dimension(0); cell++) {
1599 for (
int lbf=0; lbf<inoutOperator.dimension(1); lbf++) {
1600 for (
int rbf=0; rbf<inoutOperator.dimension(2); rbf++) {
1601 inoutOperator(cell, lbf, rbf) *= fieldSigns(cell, lbf);
1609 template<
class Scalar,
class ArrayTypeInOut,
class ArrayTypeSign>
1611 const ArrayTypeSign & fieldSigns) {
1612 #ifdef HAVE_INTREPID_DEBUG
1613 TEUCHOS_TEST_FOR_EXCEPTION( (getrank(inoutOperator) != 3), std::invalid_argument,
1614 ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Input operator container must have rank 3.");
1615 TEUCHOS_TEST_FOR_EXCEPTION( (getrank(fieldSigns) != 2), std::invalid_argument,
1616 ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Input field signs container must have rank 2.");
1617 TEUCHOS_TEST_FOR_EXCEPTION( (inoutOperator.dimension(0) != fieldSigns.dimension(0) ), std::invalid_argument,
1618 ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Zeroth dimensions (number of cells) of the operator and field signs containers must agree!");
1619 TEUCHOS_TEST_FOR_EXCEPTION( (inoutOperator.dimension(2) != fieldSigns.dimension(1) ), std::invalid_argument,
1620 ">>> ERROR (FunctionSpaceTools::applyRightFieldSigns): Second dimension of the operator container and first dimension of the field signs container (number of right fields) must agree!");
1623 for (
int cell=0; cell<inoutOperator.dimension(0); cell++) {
1624 for (
int lbf=0; lbf<inoutOperator.dimension(1); lbf++) {
1625 for (
int rbf=0; rbf<inoutOperator.dimension(2); rbf++) {
1626 inoutOperator(cell, lbf, rbf) *= fieldSigns(cell, rbf);
1635 template<
class Scalar,
class ArrayTypeInOut,
class ArrayTypeSign>
1637 const ArrayTypeSign & fieldSigns) {
1639 #ifdef HAVE_INTREPID_DEBUG
1641 TEUCHOS_TEST_FOR_EXCEPTION( ((inoutFunction.rank() < 2) || (inoutFunction.rank() > 5)), std::invalid_argument,
1642 ">>> ERROR (FunctionSpaceTools::applyFieldSigns): Input function container must have rank 2, 3, 4, or 5.");
1643 TEUCHOS_TEST_FOR_EXCEPTION( (fieldSigns.rank() != 2), std::invalid_argument,
1644 ">>> ERROR (FunctionSpaceTools::applyFieldSigns): Input field signs container must have rank 2.");
1645 TEUCHOS_TEST_FOR_EXCEPTION( (inoutFunction.dimension(0) != fieldSigns.dimension(0) ), std::invalid_argument,
1646 ">>> ERROR (FunctionSpaceTools::applyFieldSigns): Zeroth dimensions (number of integration domains) of the function and field signs containers must agree!");
1647 TEUCHOS_TEST_FOR_EXCEPTION( (inoutFunction.dimension(1) != fieldSigns.dimension(1) ), std::invalid_argument,
1648 ">>> ERROR (FunctionSpaceTools::applyFieldSigns): First dimensions (number of fields) of the function and field signs containers must agree!");
1652 ArrayWrapper<Scalar,ArrayTypeInOut, Rank<ArrayTypeInOut >::value,
false>inoutFunctionWrap(inoutFunction);
1653 ArrayWrapper<Scalar,ArrayTypeSign, Rank<ArrayTypeSign >::value,
true>fieldSignsWrap(fieldSigns);
1656 int numCells = inoutFunction.dimension(0);
1657 int numFields = inoutFunction.dimension(1);
1658 int fRank = getrank(inoutFunction);
1662 for (
int cell=0; cell<numCells; cell++) {
1663 for (
int bf=0; bf<numFields; bf++) {
1664 inoutFunctionWrap(cell, bf) *= fieldSignsWrap(cell, bf);
1671 int numPoints = inoutFunction.dimension(2);
1672 for (
int cell=0; cell<numCells; cell++) {
1673 for (
int bf=0; bf<numFields; bf++) {
1674 for (
int pt=0; pt<numPoints; pt++) {
1675 inoutFunctionWrap(cell, bf, pt) *= fieldSignsWrap(cell, bf);
1683 int numPoints = inoutFunction.dimension(2);
1684 int spaceDim1 = inoutFunction.dimension(3);
1685 for (
int cell=0; cell<numCells; cell++) {
1686 for (
int bf=0; bf<numFields; bf++) {
1687 for (
int pt=0; pt<numPoints; pt++) {
1688 for (
int d1=0; d1<spaceDim1; d1++) {
1689 inoutFunctionWrap(cell, bf, pt, d1) *= fieldSignsWrap(cell, bf);
1698 int numPoints = inoutFunction.dimension(2);
1699 int spaceDim1 = inoutFunction.dimension(3);
1700 int spaceDim2 = inoutFunction.dimension(4);
1701 for (
int cell=0; cell<numCells; cell++) {
1702 for (
int bf=0; bf<numFields; bf++) {
1703 for (
int pt=0; pt<numPoints; pt++) {
1704 for (
int d1=0; d1<spaceDim1; d1++) {
1705 for (
int d2=0; d2<spaceDim2; d2++) {
1706 inoutFunctionWrap(cell, bf, pt, d1, d2) *= fieldSignsWrap(cell, bf);
1716 #ifdef HAVE_INTREPID_DEBUG
1718 TEUCHOS_TEST_FOR_EXCEPTION( !( (inoutFunction.rank() == 2) || (inoutFunction.rank() == 3) || (inoutFunction.rank() == 4) || (inoutFunction.rank() == 5)), std::invalid_argument,
1719 ">>> ERROR (FunctionSpaceTools::applyFieldSigns): Method defined only for rank-2, 3, 4, or 5 input function containers.");
1728 template<
class Scalar,
class ArrayOutPo
intVals,
class ArrayInCoeffs,
class ArrayInFields>
1730 const ArrayInCoeffs & inCoeffs,
1731 const ArrayInFields & inFields) {
1733 #ifdef HAVE_INTREPID_DEBUG
1735 TEUCHOS_TEST_FOR_EXCEPTION( ((inFields.rank() < 3) || (inFields.rank() > 5)), std::invalid_argument,
1736 ">>> ERROR (FunctionSpaceTools::evaluate): Input fields container must have rank 3, 4, or 5.");
1737 TEUCHOS_TEST_FOR_EXCEPTION( (inCoeffs.rank() != 2), std::invalid_argument,
1738 ">>> ERROR (FunctionSpaceTools::evaluate): Input coefficient container must have rank 2.");
1739 TEUCHOS_TEST_FOR_EXCEPTION( (outPointVals.rank() != inFields.rank()-1), std::invalid_argument,
1740 ">>> ERROR (FunctionSpaceTools::evaluate): Output values container must have rank one less than the rank of the input fields container.");
1741 TEUCHOS_TEST_FOR_EXCEPTION( (inCoeffs.dimension(0) != inFields.dimension(0) ), std::invalid_argument,
1742 ">>> ERROR (FunctionSpaceTools::evaluate): Zeroth dimensions (number of cells) of the coefficient and fields input containers must agree!");
1743 TEUCHOS_TEST_FOR_EXCEPTION( (inCoeffs.dimension(1) != inFields.dimension(1) ), std::invalid_argument,
1744 ">>> ERROR (FunctionSpaceTools::evaluate): First dimensions (number of fields) of the coefficient and fields input containers must agree!");
1745 TEUCHOS_TEST_FOR_EXCEPTION( (outPointVals.dimension(0) != inFields.dimension(0) ), std::invalid_argument,
1746 ">>> ERROR (FunctionSpaceTools::evaluate): Zeroth dimensions (number of cells) of the input fields container and the output values container must agree!");
1747 for (
int i=1; i<outPointVals.rank(); i++) {
1748 std::string errmsg =
">>> ERROR (FunctionSpaceTools::evaluate): Dimensions ";
1749 errmsg += (char)(48+i);
1751 errmsg += (char)(48+i+1);
1752 errmsg +=
" of the output values and input fields containers must agree!";
1753 TEUCHOS_TEST_FOR_EXCEPTION( (outPointVals.dimension(i) != inFields.dimension(i+1)), std::invalid_argument, errmsg );
1757 ArrayWrapper<Scalar,ArrayOutPointVals, Rank<ArrayOutPointVals >::value,
false>outPointValsWrap(outPointVals);
1758 ArrayWrapper<Scalar,ArrayInCoeffs, Rank<ArrayInCoeffs>::value,
true>inCoeffsWrap(inCoeffs);
1759 ArrayWrapper<Scalar,ArrayInFields, Rank<ArrayInFields>::value,
true>inFieldsWrap(inFields);
1761 int numCells = inFields.dimension(0);
1762 int numFields = inFields.dimension(1);
1763 int numPoints = inFields.dimension(2);
1764 int fRank = getrank(inFields);
1768 for (
int cell=0; cell<numCells; cell++) {
1769 for (
int pt=0; pt<numPoints; pt++) {
1770 for (
int bf=0; bf<numFields; bf++) {
1771 outPointValsWrap(cell, pt) += inCoeffsWrap(cell, bf) * inFieldsWrap(cell, bf, pt);
1779 int spaceDim1 = inFields.dimension(3);
1780 for (
int cell=0; cell<numCells; cell++) {
1781 for (
int pt=0; pt<numPoints; pt++) {
1782 for (
int d1=0; d1<spaceDim1; d1++) {
1783 for (
int bf=0; bf<numFields; bf++) {
1784 outPointValsWrap(cell, pt, d1) += inCoeffsWrap(cell, bf) * inFieldsWrap(cell, bf, pt, d1);
1793 int spaceDim1 = inFields.dimension(3);
1794 int spaceDim2 = inFields.dimension(4);
1795 for (
int cell=0; cell<numCells; cell++) {
1796 for (
int pt=0; pt<numPoints; pt++) {
1797 for (
int d1=0; d1<spaceDim1; d1++) {
1798 for (
int d2=0; d2<spaceDim2; d2++) {
1799 for (
int bf=0; bf<numFields; bf++) {
1800 outPointValsWrap(cell, pt, d1, d2) += inCoeffsWrap(cell, bf) * inFieldsWrap(cell, bf, pt, d1, d2);
1810 #ifdef HAVE_INTREPID_DEBUG
1812 TEUCHOS_TEST_FOR_EXCEPTION( !( (fRank == 3) || (fRank == 4) || (fRank == 5)), std::invalid_argument,
1813 ">>> ERROR (FunctionSpaceTools::evaluate): Method defined only for rank-3, 4, or 5 input fields containers.");
int dimension(const int whichDim) const
Returns the specified dimension.
int rank() const
Return rank of the FieldContainer = number of indices used to tag the multi-indexed value...