Intrepid
Intrepid_ArrayToolsDefScalar.hpp
Go to the documentation of this file.
1 #ifndef INTREPID_ARRAYTOOLSDEFSCALAR_HPP
2 #define INTREPID_ARRAYTOOLSDEFSCALAR_HPP
3 // @HEADER
4 // ************************************************************************
5 //
6 // Intrepid Package
7 // Copyright (2007) Sandia Corporation
8 //
9 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
10 // license for use of this work by or on behalf of the U.S. Government.
11 //
12 // Redistribution and use in source and binary forms, with or without
13 // modification, are permitted provided that the following conditions are
14 // met:
15 //
16 // 1. Redistributions of source code must retain the above copyright
17 // notice, this list of conditions and the following disclaimer.
18 //
19 // 2. Redistributions in binary form must reproduce the above copyright
20 // notice, this list of conditions and the following disclaimer in the
21 // documentation and/or other materials provided with the distribution.
22 //
23 // 3. Neither the name of the Corporation nor the names of the
24 // contributors may be used to endorse or promote products derived from
25 // this software without specific prior written permission.
26 //
27 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
28 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
29 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
30 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
31 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
32 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
33 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
34 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
35 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
36 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
37 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
38 //
39 // Questions? Contact Pavel Bochev (pbboche@sandia.gov)
40 // Denis Ridzal (dridzal@sandia.gov), or
41 // Kara Peterson (kjpeter@sandia.gov)
42 //
43 // ************************************************************************
44 // @HEADER
45 
52 namespace Intrepid {
53 
54 template<class Scalar, class ArrayOutFields, class ArrayInData, class ArrayInFields>
55 void ArrayTools::scalarMultiplyDataField(ArrayOutFields & outputFields,
56  const ArrayInData & inputData,
57  const ArrayInFields & inputFields,
58  const bool reciprocal) {
59 #ifdef HAVE_INTREPID_DEBUG
60  TEUCHOS_TEST_FOR_EXCEPTION( (getrank(inputData) != 2), std::invalid_argument,
61  ">>> ERROR (ArrayTools::scalarMultiplyDataField): Input data container must have rank 2.");
62  if (getrank(outputFields) <= getrank(inputFields)) {
63  TEUCHOS_TEST_FOR_EXCEPTION( ( (getrank(inputFields) < 3) || (getrank(inputFields) > 5) ), std::invalid_argument,
64  ">>> ERROR (ArrayTools::scalarMultiplyDataField): Input fields container must have rank 3, 4, or 5.");
65  TEUCHOS_TEST_FOR_EXCEPTION( (getrank(outputFields) != getrank(inputFields)), std::invalid_argument,
66  ">>> ERROR (ArrayTools::scalarMultiplyDataField): Input and output fields containers must have the same rank.");
67  TEUCHOS_TEST_FOR_EXCEPTION( (inputFields.dimension(0) != inputData.dimension(0) ), std::invalid_argument,
68  ">>> ERROR (ArrayTools::scalarMultiplyDataField): Zeroth dimensions (number of integration domains) of the fields and data input containers must agree!");
69  TEUCHOS_TEST_FOR_EXCEPTION( ( (inputFields.dimension(2) != inputData.dimension(1)) && (inputData.dimension(1) != 1) ), std::invalid_argument,
70  ">>> ERROR (ArrayTools::scalarMultiplyDataField): 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!");
71  for (size_t i=0; i<getrank(inputFields); i++) {
72  std::string errmsg = ">>> ERROR (ArrayTools::scalarMultiplyDataField): Dimension ";
73  errmsg += (char)(48+i);
74  errmsg += " of the input and output fields containers must agree!";
75  TEUCHOS_TEST_FOR_EXCEPTION( (inputFields.dimension(i) != outputFields.dimension(i)), std::invalid_argument, errmsg );
76  }
77  }
78  else {
79  TEUCHOS_TEST_FOR_EXCEPTION( ( (getrank(inputFields) < 2) || (getrank(inputFields) > 4) ), std::invalid_argument,
80  ">>> ERROR (ArrayTools::scalarMultiplyDataField): Input fields container must have rank 2, 3, or 4.");
81  TEUCHOS_TEST_FOR_EXCEPTION( (getrank(outputFields) != getrank(inputFields)+1), std::invalid_argument,
82  ">>> ERROR (ArrayTools::scalarMultiplyDataField): The rank of the input fields container must be one less than the rank of the output fields container.");
83  TEUCHOS_TEST_FOR_EXCEPTION( ( (inputFields.dimension(1) != inputData.dimension(1)) && (inputData.dimension(1) != 1) ), std::invalid_argument,
84  ">>> ERROR (ArrayTools::scalarMultiplyDataField): First dimensions of fields input container and data input container (number of integration points) must agree or first data dimension must be 1!");
85  TEUCHOS_TEST_FOR_EXCEPTION( ( inputData.dimension(0) != outputFields.dimension(0) ), std::invalid_argument,
86  ">>> ERROR (ArrayTools::scalarMultiplyDataField): Zeroth dimensions of fields output container and data input containers (number of integration domains) must agree!");
87  for (size_t i=0; i<getrank(inputFields); i++) {
88  std::string errmsg = ">>> ERROR (ArrayTools::scalarMultiplyDataField): Dimensions ";
89  errmsg += (char)(48+i);
90  errmsg += " and ";
91  errmsg += (char)(48+i+1);
92  errmsg += " of the input and output fields containers must agree!";
93  TEUCHOS_TEST_FOR_EXCEPTION( (inputFields.dimension(i) != outputFields.dimension(i+1)), std::invalid_argument, errmsg );
94  }
95  }
96 #endif
97  ArrayWrapper<Scalar,ArrayOutFields, Rank<ArrayOutFields >::value, false>outputFieldsWrap(outputFields);
98  ArrayWrapper<Scalar,ArrayInData, Rank<ArrayInData >::value, true>inputDataWrap(inputData);
99  ArrayWrapper<Scalar,ArrayInFields, Rank<ArrayInFields >::value,true>inputFieldsWrap(inputFields);
100 
101  // get sizes
102  size_t invalRank = getrank(inputFields);
103  size_t outvalRank = getrank(outputFields);
104  int numCells = outputFields.dimension(0);
105  int numFields = outputFields.dimension(1);
106  int numPoints = outputFields.dimension(2);
107  int numDataPoints = inputData.dimension(1);
108  int dim1Tens = 0;
109  int dim2Tens = 0;
110  if (outvalRank > 3) {
111  dim1Tens = outputFields.dimension(3);
112  if (outvalRank > 4) {
113  dim2Tens = outputFields.dimension(4);
114  }
115  }
116 
117  if (outvalRank == invalRank) {
118 
119  if (numDataPoints != 1) { // nonconstant data
120  switch(invalRank) {
121  case 3: {
122  if (reciprocal) {
123  for(int cl = 0; cl < numCells; cl++) {
124  for(int bf = 0; bf < numFields; bf++) {
125  for(int pt = 0; pt < numPoints; pt++) {
126  outputFieldsWrap(cl, bf, pt) = inputFieldsWrap(cl, bf, pt)/inputDataWrap(cl, pt);
127  } // P-loop
128  } // F-loop
129  } // C-loop
130  }
131  else {
132 
133 
134  for(int cl = 0; cl < numCells; cl++) {
135  for(int bf = 0; bf < numFields; bf++) {
136  for(int pt = 0; pt < numPoints; pt++) {
137  outputFieldsWrap(cl, bf, pt) = inputFieldsWrap(cl, bf, pt)*inputDataWrap(cl, pt);
138  } // P-loop
139  } // F-loop
140  } // C-loop
141 
142  }
143  }// case 3
144  break;
145 
146  case 4: {
147  if (reciprocal) {
148  for(int cl = 0; cl < numCells; cl++) {
149  for(int bf = 0; bf < numFields; bf++) {
150  for(int pt = 0; pt < numPoints; pt++) {
151  for( int iVec = 0; iVec < dim1Tens; iVec++) {
152  outputFieldsWrap(cl, bf, pt, iVec) = inputFieldsWrap(cl, bf, pt, iVec)/inputDataWrap(cl, pt);
153  } // D1-loop
154  } // P-loop
155  } // F-loop
156  } // C-loop
157  }
158  else {
159  for(int cl = 0; cl < numCells; cl++) {
160  for(int bf = 0; bf < numFields; bf++) {
161  for(int pt = 0; pt < numPoints; pt++) {
162  for( int iVec = 0; iVec < dim1Tens; iVec++) {
163  outputFieldsWrap(cl, bf, pt, iVec) = inputFieldsWrap(cl, bf, pt, iVec)*inputDataWrap(cl, pt);
164  } // D1-loop
165  } // P-loop
166  } // F-loop
167  } // C-loop
168  }
169  }// case 4
170  break;
171 
172  case 5: {
173  if (reciprocal) {
174  for(int cl = 0; cl < numCells; cl++) {
175  for(int bf = 0; bf < numFields; bf++) {
176  for(int pt = 0; pt < numPoints; pt++) {
177  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
178  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
179  outputFieldsWrap(cl, bf, pt, iTens1, iTens2) = inputFieldsWrap(cl, bf, pt, iTens1, iTens2)/inputDataWrap(cl, pt);
180  } // D2-loop
181  } // D1-loop
182  } // F-loop
183  } // P-loop
184  } // C-loop
185  }
186  else {
187  for(int cl = 0; cl < numCells; cl++) {
188  for(int bf = 0; bf < numFields; bf++) {
189  for(int pt = 0; pt < numPoints; pt++) {
190  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
191  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
192  outputFieldsWrap(cl, bf, pt, iTens1, iTens2) = inputFieldsWrap(cl, bf, pt, iTens1, iTens2)*inputDataWrap(cl, pt);
193  } // D2-loop
194  } // D1-loop
195  } // F-loop
196  } // P-loop
197  } // C-loop
198  }
199  }// case 5
200  break;
201 
202  default:
203  TEUCHOS_TEST_FOR_EXCEPTION( !( (invalRank == 3) || (invalRank == 4) || (invalRank == 5) ), std::invalid_argument,
204  ">>> ERROR (ArrayTools::scalarMultiplyDataField): This branch of the method is defined only for rank-3,4 or 5 containers.");
205  }// invalRank
206 
207  }
208  else { //constant data
209 
210  switch(invalRank) {
211  case 3: {
212  if (reciprocal) {
213  for(int cl = 0; cl < numCells; cl++) {
214  for(int bf = 0; bf < numFields; bf++) {
215  for(int pt = 0; pt < numPoints; pt++) {
216  outputFieldsWrap(cl, bf, pt) = inputFieldsWrap(cl, bf, pt)/inputDataWrap(cl, 0);
217  } // P-loop
218  } // F-loop
219  } // C-loop
220  }
221  else {
222  for(int cl = 0; cl < numCells; cl++) {
223  for(int bf = 0; bf < numFields; bf++) {
224  for(int pt = 0; pt < numPoints; pt++) {
225  outputFieldsWrap(cl, bf, pt) = inputFieldsWrap(cl, bf, pt)*inputDataWrap(cl, 0);
226  } // P-loop
227  } // F-loop
228  } // C-loop
229  }
230  }// case 3
231  break;
232 
233  case 4: {
234  if (reciprocal) {
235  for(int cl = 0; cl < numCells; cl++) {
236  for(int bf = 0; bf < numFields; bf++) {
237  for(int pt = 0; pt < numPoints; pt++) {
238  for( int iVec = 0; iVec < dim1Tens; iVec++) {
239  outputFieldsWrap(cl, bf, pt, iVec) = inputFieldsWrap(cl, bf, pt, iVec)/inputDataWrap(cl, 0);
240  } // D1-loop
241  } // P-loop
242  } // F-loop
243  } // C-loop
244  }
245  else {
246  for(int cl = 0; cl < numCells; cl++) {
247  for(int bf = 0; bf < numFields; bf++) {
248  for(int pt = 0; pt < numPoints; pt++) {
249  for( int iVec = 0; iVec < dim1Tens; iVec++) {
250  outputFieldsWrap(cl, bf, pt, iVec) = inputFieldsWrap(cl, bf, pt, iVec)*inputDataWrap(cl, 0);
251  } // D1-loop
252  } // P-loop
253  } // F-loop
254  } // C-loop
255  }
256  }// case 4
257  break;
258 
259  case 5: {
260  if (reciprocal) {
261  for(int cl = 0; cl < numCells; cl++) {
262  for(int bf = 0; bf < numFields; bf++) {
263  for(int pt = 0; pt < numPoints; pt++) {
264  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
265  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
266  outputFieldsWrap(cl, bf, pt, iTens1, iTens2) = inputFieldsWrap(cl, bf, pt, iTens1, iTens2)/inputDataWrap(cl, 0);
267  } // D2-loop
268  } // D1-loop
269  } // F-loop
270  } // P-loop
271  } // C-loop
272  }
273  else {
274  for(int cl = 0; cl < numCells; cl++) {
275  for(int bf = 0; bf < numFields; bf++) {
276  for(int pt = 0; pt < numPoints; pt++) {
277  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
278  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
279  outputFieldsWrap(cl, bf, pt, iTens1, iTens2) = inputFieldsWrap(cl, bf, pt, iTens1, iTens2)*inputDataWrap(cl, 0);
280  } // D2-loop
281  } // D1-loop
282  } // F-loop
283  } // P-loop
284  } // C-loop
285  }
286  }// case 5
287  break;
288 
289  default:
290  TEUCHOS_TEST_FOR_EXCEPTION( !( (invalRank == 3) || (invalRank == 4) || (invalRank == 5) ), std::invalid_argument,
291  ">>> ERROR (ArrayTools::scalarMultiplyDataField): This branch of the method is defined only for rank-3, 4 or 5 input containers.");
292 
293  } // invalRank
294  } // numDataPoints
295 
296  }
297  else {
298 
299  if (numDataPoints != 1) { // nonconstant data
300 
301  switch(invalRank) {
302  case 2: {
303  if (reciprocal) {
304  for(int cl = 0; cl < numCells; cl++) {
305  for(int bf = 0; bf < numFields; bf++) {
306  for(int pt = 0; pt < numPoints; pt++) {
307  outputFieldsWrap(cl, bf, pt) = inputFieldsWrap(bf, pt)/inputDataWrap(cl, pt);
308  } // P-loop
309  } // F-loop
310  } // C-loop
311  }
312  else {
313  for(int cl = 0; cl < numCells; cl++) {
314  for(int bf = 0; bf < numFields; bf++) {
315  for(int pt = 0; pt < numPoints; pt++) {
316  outputFieldsWrap(cl, bf, pt) = inputFieldsWrap(bf, pt)*inputDataWrap(cl, pt);
317  } // P-loop
318  } // F-loop
319  } // C-loop
320  }
321  }// case 2
322  break;
323 
324  case 3: {
325  if (reciprocal) {
326  for(int cl = 0; cl < numCells; cl++) {
327  for(int bf = 0; bf < numFields; bf++) {
328  for(int pt = 0; pt < numPoints; pt++) {
329  for( int iVec = 0; iVec < dim1Tens; iVec++) {
330  outputFieldsWrap(cl, bf, pt, iVec) = inputFieldsWrap(bf, pt, iVec)/inputDataWrap(cl, pt);
331  } // D1-loop
332  } // P-loop
333  } // F-loop
334  } // C-loop
335  }
336  else {
337  for(int cl = 0; cl < numCells; cl++) {
338  for(int bf = 0; bf < numFields; bf++) {
339  for(int pt = 0; pt < numPoints; pt++) {
340  for( int iVec = 0; iVec < dim1Tens; iVec++) {
341  outputFieldsWrap(cl, bf, pt, iVec) = inputFieldsWrap(bf, pt, iVec)*inputDataWrap(cl, pt);
342  } // D1-loop
343  } // P-loop
344  } // F-loop
345  } // C-loop
346  }
347  }// case 3
348  break;
349 
350  case 4: {
351  if (reciprocal) {
352  for(int cl = 0; cl < numCells; cl++) {
353  for(int bf = 0; bf < numFields; bf++) {
354  for(int pt = 0; pt < numPoints; pt++) {
355  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
356  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
357  outputFieldsWrap(cl, bf, pt, iTens1, iTens2) = inputFieldsWrap(bf, pt, iTens1, iTens2)/inputDataWrap(cl, pt);
358  } // D2-loop
359  } // D1-loop
360  } // F-loop
361  } // P-loop
362  } // C-loop
363  }
364  else {
365  for(int cl = 0; cl < numCells; cl++) {
366  for(int bf = 0; bf < numFields; bf++) {
367  for(int pt = 0; pt < numPoints; pt++) {
368  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
369  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
370  outputFieldsWrap(cl, bf, pt, iTens1, iTens2) = inputFieldsWrap(bf, pt, iTens1, iTens2)*inputDataWrap(cl, pt);
371  } // D2-loop
372  } // D1-loop
373  } // F-loop
374  } // P-loop
375  } // C-loop
376  }
377  }// case 4
378  break;
379 
380  default:
381  TEUCHOS_TEST_FOR_EXCEPTION( !( (invalRank == 2) || (invalRank == 3) || (invalRank == 4) ), std::invalid_argument,
382  ">>> ERROR (ArrayTools::scalarMultiplyDataField): This branch of the method is defined only for rank-2, 3 or 4 input containers.");
383  }// invalRank
384 
385  }
386  else { //constant data
387 
388  switch(invalRank) {
389  case 2: {
390  if (reciprocal) {
391  for(int cl = 0; cl < numCells; cl++) {
392  for(int bf = 0; bf < numFields; bf++) {
393  for(int pt = 0; pt < numPoints; pt++) {
394  outputFieldsWrap(cl, bf, pt) = inputFieldsWrap(bf, pt)/inputDataWrap(cl, 0);
395  } // P-loop
396  } // F-loop
397  } // C-loop
398  }
399  else {
400  for(int cl = 0; cl < numCells; cl++) {
401  for(int bf = 0; bf < numFields; bf++) {
402  for(int pt = 0; pt < numPoints; pt++) {
403  outputFieldsWrap(cl, bf, pt) = inputFieldsWrap(bf, pt)*inputDataWrap(cl, 0);
404  } // P-loop
405  } // F-loop
406  } // C-loop
407  }
408  }// case 2
409  break;
410 
411  case 3: {
412  if (reciprocal) {
413  for(int cl = 0; cl < numCells; cl++) {
414  for(int bf = 0; bf < numFields; bf++) {
415  for(int pt = 0; pt < numPoints; pt++) {
416  for( int iVec = 0; iVec < dim1Tens; iVec++) {
417  outputFieldsWrap(cl, bf, pt, iVec) = inputFieldsWrap(bf, pt, iVec)/inputDataWrap(cl, 0);
418  } // D1-loop
419  } // P-loop
420  } // F-loop
421  } // C-loop
422  }
423  else {
424  for(int cl = 0; cl < numCells; cl++) {
425  for(int bf = 0; bf < numFields; bf++) {
426  for(int pt = 0; pt < numPoints; pt++) {
427  for( int iVec = 0; iVec < dim1Tens; iVec++) {
428  outputFieldsWrap(cl, bf, pt, iVec) = inputFieldsWrap(bf, pt, iVec)*inputDataWrap(cl, 0);
429  } // D1-loop
430  } // P-loop
431  } // F-loop
432  } // C-loop
433  }
434  }// case 3
435  break;
436 
437  case 4: {
438  if (reciprocal) {
439  for(int cl = 0; cl < numCells; cl++) {
440  for(int bf = 0; bf < numFields; bf++) {
441  for(int pt = 0; pt < numPoints; pt++) {
442  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
443  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
444  outputFieldsWrap(cl, bf, pt, iTens1, iTens2) = inputFieldsWrap(bf, pt, iTens1, iTens2)/inputDataWrap(cl, 0);
445  } // D2-loop
446  } // D1-loop
447  } // F-loop
448  } // P-loop
449  } // C-loop
450  }
451  else {
452  for(int cl = 0; cl < numCells; cl++) {
453  for(int bf = 0; bf < numFields; bf++) {
454  for(int pt = 0; pt < numPoints; pt++) {
455  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
456  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
457  outputFieldsWrap(cl, bf, pt, iTens1, iTens2) = inputFieldsWrap(bf, pt, iTens1, iTens2)*inputDataWrap(cl, 0);
458  } // D2-loop
459  } // D1-loop
460  } // F-loop
461  } // P-loop
462  } // C-loop
463  }
464  }// case 4
465  break;
466 
467  default:
468  TEUCHOS_TEST_FOR_EXCEPTION( !( (invalRank == 2) || (invalRank == 3) || (invalRank == 3) ), std::invalid_argument,
469  ">>> ERROR (ArrayTools::scalarMultiplyDataField): This branch of the method is defined only for rank-2, 3 or 4 input containers.");
470 
471  } // invalRank
472  } // numDataPoints
473 
474  } // end if (outvalRank = invalRank)
475 
476 } // scalarMultiplyDataField
477 
478 template<class Scalar, class ArrayOutData, class ArrayInDataLeft, class ArrayInDataRight>
479 void ArrayTools::scalarMultiplyDataData(ArrayOutData & outputData,
480  const ArrayInDataLeft & inputDataLeft,
481  const ArrayInDataRight & inputDataRight,
482  const bool reciprocal) {
483 
484 #ifdef HAVE_INTREPID_DEBUG
485  TEUCHOS_TEST_FOR_EXCEPTION( (getrank(inputDataLeft) != 2), std::invalid_argument,
486  ">>> ERROR (ArrayTools::scalarMultiplyDataData): Left input data container must have rank 2.");
487  if (getrank(outputData) <= getrank(inputDataRight)) {
488  TEUCHOS_TEST_FOR_EXCEPTION( ( (getrank(inputDataRight) < 2) || (getrank(inputDataRight) > 4) ), std::invalid_argument,
489  ">>> ERROR (ArrayTools::scalarMultiplyDataData): Right input data container must have rank 2, 3, or 4.");
490  TEUCHOS_TEST_FOR_EXCEPTION( (getrank(outputData) != getrank(inputDataRight)), std::invalid_argument,
491  ">>> ERROR (ArrayTools::scalarMultiplyDataData): Right input and output data containers must have the same rank.");
492  TEUCHOS_TEST_FOR_EXCEPTION( (inputDataRight.dimension(0) != inputDataLeft.dimension(0) ), std::invalid_argument,
493  ">>> ERROR (ArrayTools::scalarMultiplyDataData): Zeroth dimensions (number of integration domains) of the left and right data input containers must agree!");
494  TEUCHOS_TEST_FOR_EXCEPTION( ( (inputDataRight.dimension(1) != inputDataLeft.dimension(1)) && (inputDataLeft.dimension(1) != 1) ), std::invalid_argument,
495  ">>> ERROR (ArrayTools::scalarMultiplyDataData): First dimensions of the left and right data input containers (number of integration points) must agree or first dimension of the left data input container must be 1!");
496  for (size_t i=0; i<getrank(inputDataRight); i++) {
497  std::string errmsg = ">>> ERROR (ArrayTools::scalarMultiplyDataData): Dimension ";
498  errmsg += (char)(48+i);
499  errmsg += " of the right input and output data containers must agree!";
500  TEUCHOS_TEST_FOR_EXCEPTION( (inputDataRight.dimension(i) != outputData.dimension(i)), std::invalid_argument, errmsg );
501  }
502  }
503  else {
504  TEUCHOS_TEST_FOR_EXCEPTION( ( (getrank(inputDataRight) < 1) || (getrank(inputDataRight) > 3) ), std::invalid_argument,
505  ">>> ERROR (ArrayTools::scalarMultiplyDataData): Right input data container must have rank 1, 2, or 3.");
506  TEUCHOS_TEST_FOR_EXCEPTION( (getrank(outputData) != getrank(inputDataRight)+1), std::invalid_argument,
507  ">>> ERROR (ArrayTools::scalarMultiplyDataData): The rank of the right input data container must be one less than the rank of the output data container.");
508  TEUCHOS_TEST_FOR_EXCEPTION( ( (inputDataRight.dimension(0) != inputDataLeft.dimension(1)) && (inputDataLeft.dimension(1) != 1) ), std::invalid_argument,
509  ">>> ERROR (ArrayTools::scalarMultiplyDataData): Zeroth dimension of the right input data container and first dimension of the left data input container (number of integration points) must agree or first dimension of the left data input container must be 1!");
510  TEUCHOS_TEST_FOR_EXCEPTION( ( inputDataLeft.dimension(0) != outputData.dimension(0) ), std::invalid_argument,
511  ">>> ERROR (ArrayTools::scalarMultiplyDataData): Zeroth dimensions of data output and left data input containers (number of integration domains) must agree!");
512  for (size_t i=0; i<getrank(inputDataRight); i++) {
513  std::string errmsg = ">>> ERROR (ArrayTools::scalarMultiplyDataData): Dimensions ";
514  errmsg += (char)(48+i);
515  errmsg += " and ";
516  errmsg += (char)(48+i+1);
517  errmsg += " of the right input and output data containers must agree!";
518  TEUCHOS_TEST_FOR_EXCEPTION( (inputDataRight.dimension(i) != outputData.dimension(i+1)), std::invalid_argument, errmsg );
519  }
520  }
521 #endif
522 
523 
524  ArrayWrapper<Scalar,ArrayOutData, Rank<ArrayOutData >::value, false>outputDataWrap(outputData);
525  ArrayWrapper<Scalar,ArrayInDataLeft, Rank<ArrayInDataLeft >::value, true>inputDataLeftWrap(inputDataLeft);
526  ArrayWrapper<Scalar,ArrayInDataRight, Rank<ArrayInDataRight >::value,true>inputDataRightWrap(inputDataRight);
527 
528 
529  // get sizes
530  size_t invalRank = getrank(inputDataRight);
531  size_t outvalRank = getrank(outputData);
532  int numCells = outputData.dimension(0);
533  int numPoints = outputData.dimension(1);
534  int numDataPoints = inputDataLeft.dimension(1);
535  int dim1Tens = 0;
536  int dim2Tens = 0;
537  if (outvalRank > 2) {
538  dim1Tens = outputData.dimension(2);
539  if (outvalRank > 3) {
540  dim2Tens = outputData.dimension(3);
541  }
542  }
543 
544  if (outvalRank == invalRank) {
545 
546  if (numDataPoints != 1) { // nonconstant data
547 
548  switch(invalRank) {
549  case 2: {
550  if (reciprocal) {
551  for(int cl = 0; cl < numCells; cl++) {
552  for(int pt = 0; pt < numPoints; pt++) {
553  outputDataWrap(cl, pt) = inputDataRightWrap(cl, pt)/inputDataLeftWrap(cl, pt);
554  } // P-loop
555  } // C-loop
556  }
557  else {
558  for(int cl = 0; cl < numCells; cl++) {
559  for(int pt = 0; pt < numPoints; pt++) {
560  outputDataWrap(cl, pt) = inputDataRightWrap(cl, pt)*inputDataLeftWrap(cl, pt);
561  } // P-loop
562  } // C-loop
563  }
564  }// case 2
565  break;
566 
567  case 3: {
568  if (reciprocal) {
569  for(int cl = 0; cl < numCells; cl++) {
570  for(int pt = 0; pt < numPoints; pt++) {
571  for( int iVec = 0; iVec < dim1Tens; iVec++) {
572  outputDataWrap(cl, pt, iVec) = inputDataRightWrap(cl, pt, iVec)/inputDataLeftWrap(cl, pt);
573  } // D1-loop
574  } // P-loop
575  } // C-loop
576  }
577  else {
578  for(int cl = 0; cl < numCells; cl++) {
579  for(int pt = 0; pt < numPoints; pt++) {
580  for( int iVec = 0; iVec < dim1Tens; iVec++) {
581  outputDataWrap(cl, pt, iVec) = inputDataRightWrap(cl, pt, iVec)*inputDataLeftWrap(cl, pt);
582  } // D1-loop
583  } // P-loop
584  } // C-loop
585  }
586  }// case 3
587  break;
588 
589  case 4: {
590  if (reciprocal) {
591  for(int cl = 0; cl < numCells; cl++) {
592  for(int pt = 0; pt < numPoints; pt++) {
593  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
594  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
595  outputDataWrap(cl, pt, iTens1, iTens2) = inputDataRightWrap(cl, pt, iTens1, iTens2)/inputDataLeftWrap(cl, pt);
596  } // D2-loop
597  } // D1-loop
598  } // P-loop
599  } // C-loop
600  }
601  else {
602  for(int cl = 0; cl < numCells; cl++) {
603  for(int pt = 0; pt < numPoints; pt++) {
604  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
605  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
606  outputDataWrap(cl, pt, iTens1, iTens2) = inputDataRightWrap(cl, pt, iTens1, iTens2)*inputDataLeftWrap(cl, pt);
607  } // D2-loop
608  } // D1-loop
609  } // P-loop
610  } // C-loop
611  }
612  }// case 4
613  break;
614 
615  default:
616  TEUCHOS_TEST_FOR_EXCEPTION( !( (invalRank == 2) || (invalRank == 3) || (invalRank == 4) ), std::invalid_argument,
617  ">>> ERROR (ArrayTools::scalarMultiplyDataData): This branch of the method is defined only for rank-2, 3 or 4 containers.");
618  }// invalRank
619 
620  }
621  else { // constant left data
622 
623  switch(invalRank) {
624  case 2: {
625  if (reciprocal) {
626  for(int cl = 0; cl < numCells; cl++) {
627  for(int pt = 0; pt < numPoints; pt++) {
628  outputDataWrap(cl, pt) = inputDataRightWrap(cl, pt)/inputDataLeftWrap(cl, 0);
629  } // P-loop
630  } // C-loop
631  }
632  else {
633  for(int cl = 0; cl < numCells; cl++) {
634  for(int pt = 0; pt < numPoints; pt++) {
635  outputDataWrap(cl, pt) = inputDataRightWrap(cl, pt)*inputDataLeftWrap(cl, 0);
636  } // P-loop
637  } // C-loop
638  }
639  }// case 2
640  break;
641 
642  case 3: {
643  if (reciprocal) {
644  for(int cl = 0; cl < numCells; cl++) {
645  for(int pt = 0; pt < numPoints; pt++) {
646  for( int iVec = 0; iVec < dim1Tens; iVec++) {
647  outputDataWrap(cl, pt, iVec) = inputDataRightWrap(cl, pt, iVec)/inputDataLeftWrap(cl, 0);
648  } // D1-loop
649  } // P-loop
650  } // C-loop
651  }
652  else {
653  for(int cl = 0; cl < numCells; cl++) {
654  for(int pt = 0; pt < numPoints; pt++) {
655  for( int iVec = 0; iVec < dim1Tens; iVec++) {
656  outputDataWrap(cl, pt, iVec) = inputDataRightWrap(cl, pt, iVec)*inputDataLeftWrap(cl, 0);
657  } // D1-loop
658  } // P-loop
659  } // C-loop
660  }
661  }// case 3
662  break;
663 
664  case 4: {
665  if (reciprocal) {
666  for(int cl = 0; cl < numCells; cl++) {
667  for(int pt = 0; pt < numPoints; pt++) {
668  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
669  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
670  outputDataWrap(cl, pt, iTens1, iTens2) = inputDataRightWrap(cl, pt, iTens1, iTens2)/inputDataLeftWrap(cl, 0);
671  } // D2-loop
672  } // D1-loop
673  } // P-loop
674  } // C-loop
675  }
676  else {
677  for(int cl = 0; cl < numCells; cl++) {
678  for(int pt = 0; pt < numPoints; pt++) {
679  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
680  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
681  outputDataWrap(cl, pt, iTens1, iTens2) = inputDataRightWrap(cl, pt, iTens1, iTens2)*inputDataLeftWrap(cl, 0);
682  } // D2-loop
683  } // D1-loop
684  } // P-loop
685  } // C-loop
686  }
687  }// case 4
688  break;
689 
690  default:
691  TEUCHOS_TEST_FOR_EXCEPTION( !( (invalRank == 2) || (invalRank == 3) || (invalRank == 4) ), std::invalid_argument,
692  ">>> ERROR (ArrayTools::scalarMultiplyDataData): This branch of the method is defined only for rank-2, 3 or 4 input containers.");
693 
694  } // invalRank
695  } // numDataPoints
696 
697  }
698  else {
699 
700  if (numDataPoints != 1) { // nonconstant data
701 
702  switch(invalRank) {
703  case 1: {
704  if (reciprocal) {
705  for(int cl = 0; cl < numCells; cl++) {
706  for(int pt = 0; pt < numPoints; pt++) {
707  outputDataWrap(cl, pt) = inputDataRightWrap(pt)/inputDataLeftWrap(cl, pt);
708  } // P-loop
709  } // C-loop
710  }
711  else {
712  for(int cl = 0; cl < numCells; cl++) {
713  for(int pt = 0; pt < numPoints; pt++) {
714  outputDataWrap(cl, pt) = inputDataRightWrap(pt)*inputDataLeftWrap(cl, pt);
715  } // P-loop
716  } // C-loop
717  }
718  }// case 1
719  break;
720 
721  case 2: {
722  if (reciprocal) {
723  for(int cl = 0; cl < numCells; cl++) {
724  for(int pt = 0; pt < numPoints; pt++) {
725  for( int iVec = 0; iVec < dim1Tens; iVec++) {
726  outputDataWrap(cl, pt, iVec) = inputDataRightWrap(pt, iVec)/inputDataLeftWrap(cl, pt);
727  } // D1-loop
728  } // P-loop
729  } // C-loop
730  }
731  else {
732  for(int cl = 0; cl < numCells; cl++) {
733  for(int pt = 0; pt < numPoints; pt++) {
734  for( int iVec = 0; iVec < dim1Tens; iVec++) {
735  outputDataWrap(cl, pt, iVec) = inputDataRightWrap(pt, iVec)*inputDataLeftWrap(cl, pt);
736  } // D1-loop
737  } // P-loop
738  } // C-loop
739  }
740  }// case 2
741  break;
742 
743  case 3: {
744  if (reciprocal) {
745  for(int cl = 0; cl < numCells; cl++) {
746  for(int pt = 0; pt < numPoints; pt++) {
747  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
748  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
749  outputDataWrap(cl, pt, iTens1, iTens2) = inputDataRightWrap(pt, iTens1, iTens2)/inputDataLeftWrap(cl, pt);
750  } // D2-loop
751  } // D1-loop
752  } // P-loop
753  } // C-loop
754  }
755  else {
756  for(int cl = 0; cl < numCells; cl++) {
757  for(int pt = 0; pt < numPoints; pt++) {
758  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
759  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
760  outputDataWrap(cl, pt, iTens1, iTens2) = inputDataRightWrap(pt, iTens1, iTens2)*inputDataLeftWrap(cl, pt);
761  } // D2-loop
762  } // D1-loop
763  } // P-loop
764  } // C-loop
765  }
766  }// case 3
767  break;
768 
769  default:
770  TEUCHOS_TEST_FOR_EXCEPTION( !( (invalRank == 1) || (invalRank == 2) || (invalRank == 3) ), std::invalid_argument,
771  ">>> ERROR (ArrayTools::scalarMultiplyDataData): This branch of the method is defined only for rank-1, 2 or 3 input containers.");
772  }// invalRank
773 
774  }
775  else { //constant data
776 
777  switch(invalRank) {
778  case 1: {
779  if (reciprocal) {
780  for(int cl = 0; cl < numCells; cl++) {
781  for(int pt = 0; pt < numPoints; pt++) {
782  outputDataWrap(cl, pt) = inputDataRightWrap(pt)/inputDataLeftWrap(cl, 0);
783  } // P-loop
784  } // C-loop
785  }
786  else {
787  for(int cl = 0; cl < numCells; cl++) {
788  for(int pt = 0; pt < numPoints; pt++) {
789  outputDataWrap(cl, pt) = inputDataRightWrap(pt)*inputDataLeftWrap(cl, 0);
790  } // P-loop
791  } // C-loop
792  }
793  }// case 1
794  break;
795 
796  case 2: {
797  if (reciprocal) {
798  for(int cl = 0; cl < numCells; cl++) {
799  for(int pt = 0; pt < numPoints; pt++) {
800  for( int iVec = 0; iVec < dim1Tens; iVec++) {
801  outputDataWrap(cl, pt, iVec) = inputDataRightWrap(pt, iVec)/inputDataLeftWrap(cl, 0);
802  } // D1-loop
803  } // P-loop
804  } // C-loop
805  }
806  else {
807  for(int cl = 0; cl < numCells; cl++) {
808  for(int pt = 0; pt < numPoints; pt++) {
809  for( int iVec = 0; iVec < dim1Tens; iVec++) {
810  outputDataWrap(cl, pt, iVec) = inputDataRightWrap(pt, iVec)*inputDataLeftWrap(cl, 0);
811  } // D1-loop
812  } // P-loop
813  } // C-loop
814  }
815  }// case 2
816  break;
817 
818  case 3: {
819  if (reciprocal) {
820  for(int cl = 0; cl < numCells; cl++) {
821  for(int pt = 0; pt < numPoints; pt++) {
822  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
823  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
824  outputDataWrap(cl, pt, iTens1, iTens2) = inputDataRightWrap(pt, iTens1, iTens2)/inputDataLeftWrap(cl, 0);
825  } // D2-loop
826  } // D1-loop
827  } // P-loop
828  } // C-loop
829  }
830  else {
831  for(int cl = 0; cl < numCells; cl++) {
832  for(int pt = 0; pt < numPoints; pt++) {
833  for( int iTens1 = 0; iTens1 < dim1Tens; iTens1++) {
834  for( int iTens2 = 0; iTens2 < dim2Tens; iTens2++) {
835  outputDataWrap(cl, pt, iTens1, iTens2) = inputDataRightWrap(pt, iTens1, iTens2)*inputDataLeftWrap(cl, 0);
836  } // D2-loop
837  } // D1-loop
838  } // P-loop
839  } // C-loop
840  }
841  }// case 3
842  break;
843 
844  default:
845  TEUCHOS_TEST_FOR_EXCEPTION( !( (invalRank == 1) || (invalRank == 2) || (invalRank == 3) ), std::invalid_argument,
846  ">>> ERROR (ArrayTools::scalarMultiplyDataData): This branch of the method is defined only for rank-1, 2 or 3 input containers.");
847 
848  } // invalRank
849  } // numDataPoints
850 
851  } // end if (outvalRank = invalRank)
852 
853 } // scalarMultiplyDataData
854 } // end namespace Intrepid
855 #endif
856 
857 
858 #if defined(Intrepid_SHOW_DEPRECATED_WARNINGS)
859 #ifdef __GNUC__
860 #warning "The Intrepid package is deprecated"
861 #endif
862 #endif
863 
Intrepid::ArrayTools::scalarMultiplyDataData
static void scalarMultiplyDataData(ArrayOutData &outputData, const ArrayInDataLeft &inputDataLeft, const ArrayInDataRight &inputDataRight, const bool reciprocal=false)
There are two use cases: (1) multiplies a rank-2, 3, or 4 container inputDataRight with dimensions (C...
Definition: Intrepid_ArrayToolsDefScalar.hpp:479
Intrepid::ArrayTools::scalarMultiplyDataField
static void scalarMultiplyDataField(ArrayOutFields &outputFields, const ArrayInData &inputData, const ArrayInFields &inputFields, const bool reciprocal=false)
There are two use cases: (1) multiplies a rank-3, 4, or 5 container inputFields with dimensions (C...
Definition: Intrepid_ArrayToolsDefScalar.hpp:55
Generated by   doxygen 1.8.5