Intrepid
test_01.cpp
1 // @HEADER
2 // ************************************************************************
3 //
4 // Intrepid Package
5 // Copyright (2007) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Pavel Bochev (pbboche@sandia.gov)
38 // Denis Ridzal (dridzal@sandia.gov), or
39 // Kara Peterson (kjpeter@sandia.gov)
40 //
41 // ************************************************************************
42 // @HEADER
43 
49 #include "Intrepid_HGRAD_QUAD_C1_FEM.hpp"
50 #include "Teuchos_oblackholestream.hpp"
51 #include "Teuchos_RCP.hpp"
52 #include "Teuchos_GlobalMPISession.hpp"
53 
54 using namespace std;
55 using namespace Intrepid;
56 
57 #define INTREPID_TEST_COMMAND( S , throwCounter, nException ) \
58 { \
59  ++nException; \
60  try { \
61  S ; \
62  } \
63  catch (const std::logic_error & err) { \
64  ++throwCounter; \
65  *outStream << "Expected Error " << nException << " -------------------------------------------------------------\n"; \
66  *outStream << err.what() << '\n'; \
67  *outStream << "-------------------------------------------------------------------------------" << "\n\n"; \
68  }; \
69 }
70 
71 int main(int argc, char *argv[]) {
72 
73  Teuchos::GlobalMPISession mpiSession(&argc, &argv);
74 
75  // This little trick lets us print to std::cout only if
76  // a (dummy) command-line argument is provided.
77  int iprint = argc - 1;
78  Teuchos::RCP<std::ostream> outStream;
79  Teuchos::oblackholestream bhs; // outputs nothing
80  if (iprint > 0)
81  outStream = Teuchos::rcp(&std::cout, false);
82  else
83  outStream = Teuchos::rcp(&bhs, false);
84 
85  // Save the format state of the original std::cout.
86  Teuchos::oblackholestream oldFormatState;
87  oldFormatState.copyfmt(std::cout);
88 
89  *outStream \
90  << "===============================================================================\n" \
91  << "| |\n" \
92  << "| Unit Test (Basis_HGRAD_QUAD_C1_FEM) |\n" \
93  << "| |\n" \
94  << "| 1) Conversion of Dof tags into Dof ordinals and back |\n" \
95  << "| 2) Basis values for VALUE, GRAD, CURL, and Dk operators |\n" \
96  << "| |\n" \
97  << "| Questions? Contact Pavel Bochev (pbboche@sandia.gov), |\n" \
98  << "| Denis Ridzal (dridzal@sandia.gov), |\n" \
99  << "| Kara Peterson (kjpeter@sandia.gov). |\n" \
100  << "| |\n" \
101  << "| Intrepid's website: http://trilinos.sandia.gov/packages/intrepid |\n" \
102  << "| Trilinos website: http://trilinos.sandia.gov |\n" \
103  << "| |\n" \
104  << "===============================================================================\n"\
105  << "| TEST 1: Basis creation, exception testing |\n"\
106  << "===============================================================================\n";
107 
108  // Define basis and error flag
110  int errorFlag = 0;
111 
112  // Initialize throw counter for exception testing
113  int nException = 0;
114  int throwCounter = 0;
115 
116  // Define array containing the 4 vertices of the reference QUAD and its center.
117  FieldContainer<double> quadNodes(5, 2);
118  quadNodes(0,0) = -1.0; quadNodes(0,1) = -1.0;
119  quadNodes(1,0) = 1.0; quadNodes(1,1) = -1.0;
120  quadNodes(2,0) = 1.0; quadNodes(2,1) = 1.0;
121  quadNodes(3,0) = -1.0; quadNodes(3,1) = 1.0;
122  quadNodes(4,0) = 0.0; quadNodes(4,1) = 0.0;
123 
124  // Generic array for the output values; needs to be properly resized depending on the operator type
126 
127  try{
128  // exception #1: DIV cannot be applied to scalar functions
129  // resize vals to rank-2 container with dimensions (num. points, num. basis functions)
130  vals.resize(quadBasis.getCardinality(), quadNodes.dimension(0));
131  INTREPID_TEST_COMMAND( quadBasis.getValues(vals, quadNodes, OPERATOR_DIV), throwCounter, nException );
132 
133  // Exceptions 2-6: all bf tags/bf Ids below are wrong and should cause getDofOrdinal() and
134  // getDofTag() to access invalid array elements thereby causing bounds check exception
135  // exception #2
136  INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(3,0,0), throwCounter, nException );
137  // exception #3
138  INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(1,1,1), throwCounter, nException );
139  // exception #4
140  INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(0,4,0), throwCounter, nException );
141  // exception #5
142  INTREPID_TEST_COMMAND( quadBasis.getDofTag(5), throwCounter, nException );
143  // exception #6
144  INTREPID_TEST_COMMAND( quadBasis.getDofTag(-1), throwCounter, nException );
145 
146 #ifdef HAVE_INTREPID_DEBUG
147  // Exceptions 7- test exception handling with incorrectly dimensioned input/output arrays
148  // exception #7: input points array must be of rank-2
149  FieldContainer<double> badPoints1(4, 5, 3);
150  INTREPID_TEST_COMMAND( quadBasis.getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException );
151 
152  // exception #8 dimension 1 in the input point array must equal space dimension of the cell
153  FieldContainer<double> badPoints2(4, 3);
154  INTREPID_TEST_COMMAND( quadBasis.getValues(vals, badPoints2, OPERATOR_VALUE), throwCounter, nException );
155 
156  // exception #9 output values must be of rank-2 for OPERATOR_VALUE
157  FieldContainer<double> badVals1(4, 3, 1);
158  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals1, quadNodes, OPERATOR_VALUE), throwCounter, nException );
159 
160  // exception #10 output values must be of rank-3 for OPERATOR_GRAD
161  FieldContainer<double> badVals2(4, 3);
162  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_GRAD), throwCounter, nException );
163 
164  // exception #11 output values must be of rank-3 for OPERATOR_CURL
165  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_CURL), throwCounter, nException );
166 
167  // exception #12 output values must be of rank-3 for OPERATOR_D2
168  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_D2), throwCounter, nException );
169 
170  // exception #13 incorrect 0th dimension of output array (must equal number of basis functions)
171  FieldContainer<double> badVals3(quadBasis.getCardinality() + 1, quadNodes.dimension(0));
172  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals3, quadNodes, OPERATOR_VALUE), throwCounter, nException );
173 
174  // exception #14 incorrect 1st dimension of output array (must equal number of points)
175  FieldContainer<double> badVals4(quadBasis.getCardinality(), quadNodes.dimension(0) + 1);
176  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals4, quadNodes, OPERATOR_VALUE), throwCounter, nException );
177 
178  // exception #15: incorrect 2nd dimension of output array (must equal the space dimension)
179  FieldContainer<double> badVals5(quadBasis.getCardinality(), quadNodes.dimension(0), 4);
180  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals5, quadNodes, OPERATOR_GRAD), throwCounter, nException );
181 
182  // exception #16: incorrect 2nd dimension of output array (must equal D2 cardinality in 2D)
183  FieldContainer<double> badVals6(quadBasis.getCardinality(), quadNodes.dimension(0), 40);
184  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals6, quadNodes, OPERATOR_D2), throwCounter, nException );
185 
186  // exception #17: incorrect 2nd dimension of output array (must equal D3 cardinality in 2D)
187  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals6, quadNodes, OPERATOR_D3), throwCounter, nException );
188 #endif
189 
190  }
191  catch (const std::logic_error & err) {
192  *outStream << "UNEXPECTED ERROR !!! ----------------------------------------------------------\n";
193  *outStream << err.what() << '\n';
194  *outStream << "-------------------------------------------------------------------------------" << "\n\n";
195  errorFlag = -1000;
196  };
197 
198  // Check if number of thrown exceptions matches the one we expect
199  if (throwCounter != nException) {
200  errorFlag++;
201  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
202  }
203 
204  *outStream \
205  << "\n"
206  << "===============================================================================\n"\
207  << "| TEST 2: correctness of tag to enum and enum to tag lookups |\n"\
208  << "===============================================================================\n";
209 
210  try{
211  std::vector<std::vector<int> > allTags = quadBasis.getAllDofTags();
212 
213  // Loop over all tags, lookup the associated dof enumeration and then lookup the tag again
214  for (unsigned i = 0; i < allTags.size(); i++) {
215  int bfOrd = quadBasis.getDofOrdinal(allTags[i][0], allTags[i][1], allTags[i][2]);
216 
217  std::vector<int> myTag = quadBasis.getDofTag(bfOrd);
218  if( !( (myTag[0] == allTags[i][0]) &&
219  (myTag[1] == allTags[i][1]) &&
220  (myTag[2] == allTags[i][2]) &&
221  (myTag[3] == allTags[i][3]) ) ) {
222  errorFlag++;
223  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
224  *outStream << " getDofOrdinal( {"
225  << allTags[i][0] << ", "
226  << allTags[i][1] << ", "
227  << allTags[i][2] << ", "
228  << allTags[i][3] << "}) = " << bfOrd <<" but \n";
229  *outStream << " getDofTag(" << bfOrd << ") = { "
230  << myTag[0] << ", "
231  << myTag[1] << ", "
232  << myTag[2] << ", "
233  << myTag[3] << "}\n";
234  }
235  }
236 
237  // Now do the same but loop over basis functions
238  for( int bfOrd = 0; bfOrd < quadBasis.getCardinality(); bfOrd++) {
239  std::vector<int> myTag = quadBasis.getDofTag(bfOrd);
240  int myBfOrd = quadBasis.getDofOrdinal(myTag[0], myTag[1], myTag[2]);
241  if( bfOrd != myBfOrd) {
242  errorFlag++;
243  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
244  *outStream << " getDofTag(" << bfOrd << ") = { "
245  << myTag[0] << ", "
246  << myTag[1] << ", "
247  << myTag[2] << ", "
248  << myTag[3] << "} but getDofOrdinal({"
249  << myTag[0] << ", "
250  << myTag[1] << ", "
251  << myTag[2] << ", "
252  << myTag[3] << "} ) = " << myBfOrd << "\n";
253  }
254  }
255  }
256  catch (const std::logic_error & err){
257  *outStream << err.what() << "\n\n";
258  errorFlag = -1000;
259  };
260 
261  *outStream \
262  << "\n"
263  << "===============================================================================\n"\
264  << "| TEST 3: correctness of basis function values |\n"\
265  << "===============================================================================\n";
266 
267  outStream -> precision(20);
268 
269  // VALUE: Each row gives the 4 correct basis set values at an evaluation point
270  double basisValues[] = {
271  1.0, 0.0, 0.0, 0.0,
272  0.0, 1.0, 0.0, 0.0,
273  0.0, 0.0, 1.0, 0.0,
274  0.0, 0.0, 0.0, 1.0,
275  0.25,0.25,0.25,0.25
276  };
277 
278  // GRAD and D1: each row gives the 8 correct values of the gradients of the 4 basis functions
279  double basisGrads[] = {
280  -0.5, -0.5, 0.5, 0.0, 0.0, 0.0, 0.0, 0.5,
281  -0.5, 0.0, 0.5, -0.5, 0.0, 0.5, 0.0, 0.0,
282  0.0, 0.0, 0.0, -0.5, 0.5, 0.5, -0.5, 0.0,
283  0.0, -0.5, 0.0, 0.0, 0.5, 0.0, -0.5, 0.5,
284  -0.25,-0.25, 0.25,-0.25, 0.25, 0.25, -0.25, 0.25
285  };
286 
287  // CURL: each row gives the 8 correct values of the curls of the 4 basis functions
288  double basisCurls[] = {
289  -0.5, 0.5, 0.0, -0.5, 0.0, 0.0, 0.5, 0.0,
290  0.0, 0.5, -0.5, -0.5, 0.5, 0.0, 0.0, 0.0,
291  0.0, 0.0, -0.5, 0.0, 0.5, -0.5, 0.0, 0.5,
292  -0.5, 0.0, 0.0, 0.0, 0.0, -0.5, 0.5, 0.5,
293  -0.25, 0.25, -0.25,-0.25, 0.25,-0.25, 0.25, 0.25
294  };
295 
296  //D2: each row gives the 12 correct values of all 2nd derivatives of the 4 basis functions
297  double basisD2[] = {
298  0.0, 0.25, 0.0, 0.0,-0.25, 0.0, 0.0, 0.25, 0.0, 0.0,-0.25, 0.0,
299  0.0, 0.25, 0.0, 0.0,-0.25, 0.0, 0.0, 0.25, 0.0, 0.0,-0.25, 0.0,
300  0.0, 0.25, 0.0, 0.0,-0.25, 0.0, 0.0, 0.25, 0.0, 0.0,-0.25, 0.0,
301  0.0, 0.25, 0.0, 0.0,-0.25, 0.0, 0.0, 0.25, 0.0, 0.0,-0.25, 0.0,
302  0.0, 0.25, 0.0, 0.0,-0.25, 0.0, 0.0, 0.25, 0.0, 0.0,-0.25, 0.0,
303  };
304 
305  try{
306 
307  // Dimensions for the output arrays:
308  int numFields = quadBasis.getCardinality();
309  int numPoints = quadNodes.dimension(0);
310  int spaceDim = quadBasis.getBaseCellTopology().getDimension();
311  int D2Cardin = Intrepid::getDkCardinality(OPERATOR_D2, spaceDim);
312 
313  // Generic array for values, grads, curls, etc. that will be properly sized before each call
315 
316  // Check VALUE of basis functions: resize vals to rank-2 container:
317  vals.resize(numFields, numPoints);
318  quadBasis.getValues(vals, quadNodes, OPERATOR_VALUE);
319  for (int i = 0; i < numFields; i++) {
320  for (int j = 0; j < numPoints; j++) {
321  int l = i + j * numFields;
322  if (std::abs(vals(i,j) - basisValues[l]) > INTREPID_TOL) {
323  errorFlag++;
324  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
325 
326  // Output the multi-index of the value where the error is:
327  *outStream << " At multi-index { ";
328  *outStream << i << " ";*outStream << j << " ";
329  *outStream << "} computed value: " << vals(i,j)
330  << " but reference value: " << basisValues[l] << "\n";
331  }
332  }
333  }
334 
335  // Check GRAD of basis function: resize vals to rank-3 container
336  vals.resize(numFields, numPoints, spaceDim);
337  quadBasis.getValues(vals, quadNodes, OPERATOR_GRAD);
338  for (int i = 0; i < numFields; i++) {
339  for (int j = 0; j < numPoints; j++) {
340  for (int k = 0; k < spaceDim; k++) {
341  int l = k + i * spaceDim + j * spaceDim * numFields;
342  if (std::abs(vals(i,j,k) - basisGrads[l]) > INTREPID_TOL) {
343  errorFlag++;
344  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
345 
346  // Output the multi-index of the value where the error is:
347  *outStream << " At multi-index { ";
348  *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
349  *outStream << "} computed grad component: " << vals(i,j,k)
350  << " but reference grad component: " << basisGrads[l] << "\n";
351  }
352  }
353  }
354  }
355 
356 
357  // Check D1 of basis function (do not resize vals because it has the correct size: D1 = GRAD)
358  quadBasis.getValues(vals, quadNodes, OPERATOR_D1);
359  for (int i = 0; i < numFields; i++) {
360  for (int j = 0; j < numPoints; j++) {
361  for (int k = 0; k < spaceDim; k++) {
362  int l = k + i * spaceDim + j * spaceDim * numFields;
363  if (std::abs(vals(i,j,k) - basisGrads[l]) > INTREPID_TOL) {
364  errorFlag++;
365  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
366 
367  // Output the multi-index of the value where the error is:
368  *outStream << " At multi-index { ";
369  *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
370  *outStream << "} computed D1 component: " << vals(i,j,k)
371  << " but reference D1 component: " << basisGrads[l] << "\n";
372  }
373  }
374  }
375  }
376 
377 
378  // Check CURL of basis function: resize vals just for illustration!
379  vals.resize(numFields, numPoints, spaceDim);
380  quadBasis.getValues(vals, quadNodes, OPERATOR_CURL);
381  for (int i = 0; i < numFields; i++) {
382  for (int j = 0; j < numPoints; j++) {
383  for (int k = 0; k < spaceDim; k++) {
384  int l = k + i * spaceDim + j * spaceDim * numFields;
385  if (std::abs(vals(i,j,k) - basisCurls[l]) > INTREPID_TOL) {
386  errorFlag++;
387  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
388 
389  // Output the multi-index of the value where the error is:
390  *outStream << " At multi-index { ";
391  *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
392  *outStream << "} computed curl component: " << vals(i,j,k)
393  << " but reference curl component: " << basisCurls[l] << "\n";
394  }
395  }
396  }
397  }
398 
399  // Check D2 of basis function
400  vals.resize(numFields, numPoints, D2Cardin);
401  quadBasis.getValues(vals, quadNodes, OPERATOR_D2);
402  for (int i = 0; i < numFields; i++) {
403  for (int j = 0; j < numPoints; j++) {
404  for (int k = 0; k < D2Cardin; k++) {
405  int l = k + i * D2Cardin + j * D2Cardin * numFields;
406  if (std::abs(vals(i,j,k) - basisD2[l]) > INTREPID_TOL) {
407  errorFlag++;
408  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
409 
410  // Output the multi-index of the value where the error is:
411  *outStream << " At multi-index { ";
412  *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
413  *outStream << "} computed D2 component: " << vals(i,j,k)
414  << " but reference D2 component: " << basisD2[l] << "\n";
415  }
416  }
417  }
418  }
419 
420 
421  // Check all higher derivatives - must be zero.
422  for(EOperator op = OPERATOR_D3; op < OPERATOR_MAX; op++) {
423 
424  // The last dimension is the number of kth derivatives and needs to be resized for every Dk
425  int DkCardin = Intrepid::getDkCardinality(op, spaceDim);
426  vals.resize(numFields, numPoints, DkCardin);
427 
428  quadBasis.getValues(vals, quadNodes, op);
429  for (int i = 0; i < vals.size(); i++) {
430  if (std::abs(vals[i]) > INTREPID_TOL) {
431  errorFlag++;
432  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
433 
434  // Get the multi-index of the value where the error is and the operator order
435  std::vector<int> myIndex;
436  vals.getMultiIndex(myIndex,i);
437  int ord = Intrepid::getOperatorOrder(op);
438  *outStream << " At multi-index { ";
439  for(int j = 0; j < vals.rank(); j++) {
440  *outStream << myIndex[j] << " ";
441  }
442  *outStream << "} computed D"<< ord <<" component: " << vals[i]
443  << " but reference D" << ord << " component: 0 \n";
444  }
445  }
446  }
447  }
448  // Catch unexpected errors
449  catch (const std::logic_error & err) {
450  *outStream << err.what() << "\n\n";
451  errorFlag = -1000;
452  };
453 
454 
455  *outStream \
456  << "\n"
457  << "===============================================================================\n"\
458  << "| TEST 4: correctness of DoF locations |\n"\
459  << "===============================================================================\n";
460 
461  try{
462  Teuchos::RCP<Basis<double, FieldContainer<double> > > basis =
463  Teuchos::rcp(new Basis_HGRAD_QUAD_C1_FEM<double, FieldContainer<double> >);
464  Teuchos::RCP<DofCoordsInterface<FieldContainer<double> > > coord_iface =
465  Teuchos::rcp_dynamic_cast<DofCoordsInterface<FieldContainer<double> > >(basis);
466 
468  FieldContainer<double> bvals(basis->getCardinality(), basis->getCardinality());
469 
470  // Check exceptions.
471 #ifdef HAVE_INTREPID_DEBUG
472  cvals.resize(1,2,3);
473  INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
474  cvals.resize(3,2);
475  INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
476  cvals.resize(4,3);
477  INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
478 #endif
479  cvals.resize(4,2);
480  INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException ); nException--;
481  // Check if number of thrown exceptions matches the one we expect
482  if (throwCounter != nException) {
483  errorFlag++;
484  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
485  }
486 
487  // Check mathematical correctness.
488  basis->getValues(bvals, cvals, OPERATOR_VALUE);
489  char buffer[120];
490  for (int i=0; i<bvals.dimension(0); i++) {
491  for (int j=0; j<bvals.dimension(1); j++) {
492  if ((i != j) && (std::abs(bvals(i,j) - 0.0) > INTREPID_TOL)) {
493  errorFlag++;
494  sprintf(buffer, "\nValue of basis function %d at (%6.4e, %6.4e) is %6.4e but should be %6.4e!\n", i, cvals(i,0), cvals(i,1), bvals(i,j), 0.0);
495  *outStream << buffer;
496  }
497  else if ((i == j) && (std::abs(bvals(i,j) - 1.0) > INTREPID_TOL)) {
498  errorFlag++;
499  sprintf(buffer, "\nValue of basis function %d at (%6.4e, %6.4e) is %6.4e but should be %6.4e!\n", i, cvals(i,0), cvals(i,1), bvals(i,j), 1.0);
500  *outStream << buffer;
501  }
502  }
503  }
504 
505  }
506  catch (const std::logic_error & err){
507  *outStream << err.what() << "\n\n";
508  errorFlag = -1000;
509  };
510 
511  if (errorFlag != 0)
512  std::cout << "End Result: TEST FAILED\n";
513  else
514  std::cout << "End Result: TEST PASSED\n";
515 
516  // reset format state of std::cout
517  std::cout.copyfmt(oldFormatState);
518 
519  return errorFlag;
520 }
int size() const
Returns size of the FieldContainer defined as the product of its dimensions.
virtual const std::vector< int > & getDofTag(const int dofOrd)
DoF ordinal to DoF tag lookup.
virtual int getCardinality() const
Returns cardinality of the basis.
This is an interface class for bases whose degrees of freedom can be associated with spatial location...
virtual const std::vector< std::vector< int > > & getAllDofTags()
Retrieves all DoF tags.
Header file for utility class to provide multidimensional containers.
void getMultiIndex(int &i0, const int valueEnum) const
Returns the multi-index of a value, based on its enumeration, as a list, for rank-1 containers...
int rank() const
Return rank of the FieldContainer = number of indices used to tag the multi-indexed value...
Implementation of the default H(grad)-compatible FEM basis of degree 1 on Quadrilateral cell...
void resize(const int dim0)
Resizes FieldContainer to a rank-1 container with the specified dimension, initialized by 0...
void getValues(ArrayScalar &outputValues, const ArrayScalar &inputPoints, const EOperator operatorType) const
FEM basis evaluation on a reference Quadrilateral cell.
virtual int getDofOrdinal(const int subcDim, const int subcOrd, const int subcDofOrd)
DoF tag to ordinal lookup.
virtual const shards::CellTopology getBaseCellTopology() const
Returns the base cell topology for which the basis is defined. See Shards documentation http://trilin...