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 
125  try{
126  // Generic array for the output values; needs to be properly resized depending on the operator type
128 
129  // exception #1: DIV cannot be applied to scalar functions
130  // resize vals to rank-2 container with dimensions (num. points, num. basis functions)
131  vals.resize(quadBasis.getCardinality(), quadNodes.dimension(0));
132  INTREPID_TEST_COMMAND( quadBasis.getValues(vals, quadNodes, OPERATOR_DIV), throwCounter, nException );
133 
134  // Exceptions 2-6: all bf tags/bf Ids below are wrong and should cause getDofOrdinal() and
135  // getDofTag() to access invalid array elements thereby causing bounds check exception
136  // exception #2
137  INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(3,0,0), throwCounter, nException );
138  // exception #3
139  INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(1,1,1), throwCounter, nException );
140  // exception #4
141  INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(0,4,0), throwCounter, nException );
142  // exception #5
143  INTREPID_TEST_COMMAND( quadBasis.getDofTag(5), throwCounter, nException );
144  // exception #6
145  INTREPID_TEST_COMMAND( quadBasis.getDofTag(-1), throwCounter, nException );
146 
147 #ifdef HAVE_INTREPID_DEBUG
148  // Exceptions 7- test exception handling with incorrectly dimensioned input/output arrays
149  // exception #7: input points array must be of rank-2
150  FieldContainer<double> badPoints1(4, 5, 3);
151  INTREPID_TEST_COMMAND( quadBasis.getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException );
152 
153  // exception #8 dimension 1 in the input point array must equal space dimension of the cell
154  FieldContainer<double> badPoints2(4, 3);
155  INTREPID_TEST_COMMAND( quadBasis.getValues(vals, badPoints2, OPERATOR_VALUE), throwCounter, nException );
156 
157  // exception #9 output values must be of rank-2 for OPERATOR_VALUE
158  FieldContainer<double> badVals1(4, 3, 1);
159  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals1, quadNodes, OPERATOR_VALUE), throwCounter, nException );
160 
161  // exception #10 output values must be of rank-3 for OPERATOR_GRAD
162  FieldContainer<double> badVals2(4, 3);
163  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_GRAD), throwCounter, nException );
164 
165  // exception #11 output values must be of rank-3 for OPERATOR_CURL
166  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_CURL), throwCounter, nException );
167 
168  // exception #12 output values must be of rank-3 for OPERATOR_D2
169  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_D2), throwCounter, nException );
170 
171  // exception #13 incorrect 0th dimension of output array (must equal number of basis functions)
172  FieldContainer<double> badVals3(quadBasis.getCardinality() + 1, quadNodes.dimension(0));
173  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals3, quadNodes, OPERATOR_VALUE), throwCounter, nException );
174 
175  // exception #14 incorrect 1st dimension of output array (must equal number of points)
176  FieldContainer<double> badVals4(quadBasis.getCardinality(), quadNodes.dimension(0) + 1);
177  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals4, quadNodes, OPERATOR_VALUE), throwCounter, nException );
178 
179  // exception #15: incorrect 2nd dimension of output array (must equal the space dimension)
180  FieldContainer<double> badVals5(quadBasis.getCardinality(), quadNodes.dimension(0), 4);
181  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals5, quadNodes, OPERATOR_GRAD), throwCounter, nException );
182 
183  // exception #16: incorrect 2nd dimension of output array (must equal D2 cardinality in 2D)
184  FieldContainer<double> badVals6(quadBasis.getCardinality(), quadNodes.dimension(0), 40);
185  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals6, quadNodes, OPERATOR_D2), throwCounter, nException );
186 
187  // exception #17: incorrect 2nd dimension of output array (must equal D3 cardinality in 2D)
188  INTREPID_TEST_COMMAND( quadBasis.getValues(badVals6, quadNodes, OPERATOR_D3), throwCounter, nException );
189 #endif
190 
191  }
192  catch (const std::logic_error & err) {
193  *outStream << "UNEXPECTED ERROR !!! ----------------------------------------------------------\n";
194  *outStream << err.what() << '\n';
195  *outStream << "-------------------------------------------------------------------------------" << "\n\n";
196  errorFlag = -1000;
197  };
198 
199  // Check if number of thrown exceptions matches the one we expect
200  if (throwCounter != nException) {
201  errorFlag++;
202  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
203  }
204 
205  *outStream \
206  << "\n"
207  << "===============================================================================\n"\
208  << "| TEST 2: correctness of tag to enum and enum to tag lookups |\n"\
209  << "===============================================================================\n";
210 
211  try{
212  std::vector<std::vector<int> > allTags = quadBasis.getAllDofTags();
213 
214  // Loop over all tags, lookup the associated dof enumeration and then lookup the tag again
215  for (unsigned i = 0; i < allTags.size(); i++) {
216  int bfOrd = quadBasis.getDofOrdinal(allTags[i][0], allTags[i][1], allTags[i][2]);
217 
218  std::vector<int> myTag = quadBasis.getDofTag(bfOrd);
219  if( !( (myTag[0] == allTags[i][0]) &&
220  (myTag[1] == allTags[i][1]) &&
221  (myTag[2] == allTags[i][2]) &&
222  (myTag[3] == allTags[i][3]) ) ) {
223  errorFlag++;
224  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
225  *outStream << " getDofOrdinal( {"
226  << allTags[i][0] << ", "
227  << allTags[i][1] << ", "
228  << allTags[i][2] << ", "
229  << allTags[i][3] << "}) = " << bfOrd <<" but \n";
230  *outStream << " getDofTag(" << bfOrd << ") = { "
231  << myTag[0] << ", "
232  << myTag[1] << ", "
233  << myTag[2] << ", "
234  << myTag[3] << "}\n";
235  }
236  }
237 
238  // Now do the same but loop over basis functions
239  for( int bfOrd = 0; bfOrd < quadBasis.getCardinality(); bfOrd++) {
240  std::vector<int> myTag = quadBasis.getDofTag(bfOrd);
241  int myBfOrd = quadBasis.getDofOrdinal(myTag[0], myTag[1], myTag[2]);
242  if( bfOrd != myBfOrd) {
243  errorFlag++;
244  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
245  *outStream << " getDofTag(" << bfOrd << ") = { "
246  << myTag[0] << ", "
247  << myTag[1] << ", "
248  << myTag[2] << ", "
249  << myTag[3] << "} but getDofOrdinal({"
250  << myTag[0] << ", "
251  << myTag[1] << ", "
252  << myTag[2] << ", "
253  << myTag[3] << "} ) = " << myBfOrd << "\n";
254  }
255  }
256  }
257  catch (const std::logic_error & err){
258  *outStream << err.what() << "\n\n";
259  errorFlag = -1000;
260  };
261 
262  *outStream \
263  << "\n"
264  << "===============================================================================\n"\
265  << "| TEST 3: correctness of basis function values |\n"\
266  << "===============================================================================\n";
267 
268  outStream -> precision(20);
269 
270  // VALUE: Each row gives the 4 correct basis set values at an evaluation point
271  double basisValues[] = {
272  1.0, 0.0, 0.0, 0.0,
273  0.0, 1.0, 0.0, 0.0,
274  0.0, 0.0, 1.0, 0.0,
275  0.0, 0.0, 0.0, 1.0,
276  0.25,0.25,0.25,0.25
277  };
278 
279  // GRAD and D1: each row gives the 8 correct values of the gradients of the 4 basis functions
280  double basisGrads[] = {
281  -0.5, -0.5, 0.5, 0.0, 0.0, 0.0, 0.0, 0.5,
282  -0.5, 0.0, 0.5, -0.5, 0.0, 0.5, 0.0, 0.0,
283  0.0, 0.0, 0.0, -0.5, 0.5, 0.5, -0.5, 0.0,
284  0.0, -0.5, 0.0, 0.0, 0.5, 0.0, -0.5, 0.5,
285  -0.25,-0.25, 0.25,-0.25, 0.25, 0.25, -0.25, 0.25
286  };
287 
288  // CURL: each row gives the 8 correct values of the curls of the 4 basis functions
289  double basisCurls[] = {
290  -0.5, 0.5, 0.0, -0.5, 0.0, 0.0, 0.5, 0.0,
291  0.0, 0.5, -0.5, -0.5, 0.5, 0.0, 0.0, 0.0,
292  0.0, 0.0, -0.5, 0.0, 0.5, -0.5, 0.0, 0.5,
293  -0.5, 0.0, 0.0, 0.0, 0.0, -0.5, 0.5, 0.5,
294  -0.25, 0.25, -0.25,-0.25, 0.25,-0.25, 0.25, 0.25
295  };
296 
297  //D2: each row gives the 12 correct values of all 2nd derivatives of the 4 basis functions
298  double basisD2[] = {
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  0.0, 0.25, 0.0, 0.0,-0.25, 0.0, 0.0, 0.25, 0.0, 0.0,-0.25, 0.0,
304  };
305 
306  try{
307 
308  // Dimensions for the output arrays:
309  int numFields = quadBasis.getCardinality();
310  int numPoints = quadNodes.dimension(0);
311  int spaceDim = quadBasis.getBaseCellTopology().getDimension();
312  int D2Cardin = Intrepid::getDkCardinality(OPERATOR_D2, spaceDim);
313 
314  // Generic array for values, grads, curls, etc. that will be properly sized before each call
316 
317  // Check VALUE of basis functions: resize vals to rank-2 container:
318  vals.resize(numFields, numPoints);
319  quadBasis.getValues(vals, quadNodes, OPERATOR_VALUE);
320  for (int i = 0; i < numFields; i++) {
321  for (int j = 0; j < numPoints; j++) {
322  int l = i + j * numFields;
323  if (std::abs(vals(i,j) - basisValues[l]) > INTREPID_TOL) {
324  errorFlag++;
325  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
326 
327  // Output the multi-index of the value where the error is:
328  *outStream << " At multi-index { ";
329  *outStream << i << " ";*outStream << j << " ";
330  *outStream << "} computed value: " << vals(i,j)
331  << " but reference value: " << basisValues[l] << "\n";
332  }
333  }
334  }
335 
336  // Check GRAD of basis function: resize vals to rank-3 container
337  vals.resize(numFields, numPoints, spaceDim);
338  quadBasis.getValues(vals, quadNodes, OPERATOR_GRAD);
339  for (int i = 0; i < numFields; i++) {
340  for (int j = 0; j < numPoints; j++) {
341  for (int k = 0; k < spaceDim; k++) {
342  int l = k + i * spaceDim + j * spaceDim * numFields;
343  if (std::abs(vals(i,j,k) - basisGrads[l]) > INTREPID_TOL) {
344  errorFlag++;
345  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
346 
347  // Output the multi-index of the value where the error is:
348  *outStream << " At multi-index { ";
349  *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
350  *outStream << "} computed grad component: " << vals(i,j,k)
351  << " but reference grad component: " << basisGrads[l] << "\n";
352  }
353  }
354  }
355  }
356 
357 
358  // Check D1 of basis function (do not resize vals because it has the correct size: D1 = GRAD)
359  quadBasis.getValues(vals, quadNodes, OPERATOR_D1);
360  for (int i = 0; i < numFields; i++) {
361  for (int j = 0; j < numPoints; j++) {
362  for (int k = 0; k < spaceDim; k++) {
363  int l = k + i * spaceDim + j * spaceDim * numFields;
364  if (std::abs(vals(i,j,k) - basisGrads[l]) > INTREPID_TOL) {
365  errorFlag++;
366  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
367 
368  // Output the multi-index of the value where the error is:
369  *outStream << " At multi-index { ";
370  *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
371  *outStream << "} computed D1 component: " << vals(i,j,k)
372  << " but reference D1 component: " << basisGrads[l] << "\n";
373  }
374  }
375  }
376  }
377 
378 
379  // Check CURL of basis function: resize vals just for illustration!
380  vals.resize(numFields, numPoints, spaceDim);
381  quadBasis.getValues(vals, quadNodes, OPERATOR_CURL);
382  for (int i = 0; i < numFields; i++) {
383  for (int j = 0; j < numPoints; j++) {
384  for (int k = 0; k < spaceDim; k++) {
385  int l = k + i * spaceDim + j * spaceDim * numFields;
386  if (std::abs(vals(i,j,k) - basisCurls[l]) > INTREPID_TOL) {
387  errorFlag++;
388  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
389 
390  // Output the multi-index of the value where the error is:
391  *outStream << " At multi-index { ";
392  *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
393  *outStream << "} computed curl component: " << vals(i,j,k)
394  << " but reference curl component: " << basisCurls[l] << "\n";
395  }
396  }
397  }
398  }
399 
400  // Check D2 of basis function
401  vals.resize(numFields, numPoints, D2Cardin);
402  quadBasis.getValues(vals, quadNodes, OPERATOR_D2);
403  for (int i = 0; i < numFields; i++) {
404  for (int j = 0; j < numPoints; j++) {
405  for (int k = 0; k < D2Cardin; k++) {
406  int l = k + i * D2Cardin + j * D2Cardin * numFields;
407  if (std::abs(vals(i,j,k) - basisD2[l]) > INTREPID_TOL) {
408  errorFlag++;
409  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
410 
411  // Output the multi-index of the value where the error is:
412  *outStream << " At multi-index { ";
413  *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
414  *outStream << "} computed D2 component: " << vals(i,j,k)
415  << " but reference D2 component: " << basisD2[l] << "\n";
416  }
417  }
418  }
419  }
420 
421 
422  // Check all higher derivatives - must be zero.
423  for(EOperator op = OPERATOR_D3; op < OPERATOR_MAX; op++) {
424 
425  // The last dimension is the number of kth derivatives and needs to be resized for every Dk
426  int DkCardin = Intrepid::getDkCardinality(op, spaceDim);
427  vals.resize(numFields, numPoints, DkCardin);
428 
429  quadBasis.getValues(vals, quadNodes, op);
430  for (int i = 0; i < vals.size(); i++) {
431  if (std::abs(vals[i]) > INTREPID_TOL) {
432  errorFlag++;
433  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
434 
435  // Get the multi-index of the value where the error is and the operator order
436  std::vector<int> myIndex;
437  vals.getMultiIndex(myIndex,i);
438  int ord = Intrepid::getOperatorOrder(op);
439  *outStream << " At multi-index { ";
440  for(int j = 0; j < vals.rank(); j++) {
441  *outStream << myIndex[j] << " ";
442  }
443  *outStream << "} computed D"<< ord <<" component: " << vals[i]
444  << " but reference D" << ord << " component: 0 \n";
445  }
446  }
447  }
448  }
449  // Catch unexpected errors
450  catch (const std::logic_error & err) {
451  *outStream << err.what() << "\n\n";
452  errorFlag = -1000;
453  };
454 
455 
456  *outStream \
457  << "\n"
458  << "===============================================================================\n"\
459  << "| TEST 4: correctness of DoF locations |\n"\
460  << "===============================================================================\n";
461 
462  try{
463  Teuchos::RCP<Basis<double, FieldContainer<double> > > basis =
464  Teuchos::rcp(new Basis_HGRAD_QUAD_C1_FEM<double, FieldContainer<double> >);
465  Teuchos::RCP<DofCoordsInterface<FieldContainer<double> > > coord_iface =
466  Teuchos::rcp_dynamic_cast<DofCoordsInterface<FieldContainer<double> > >(basis);
467 
469  FieldContainer<double> bvals(basis->getCardinality(), basis->getCardinality());
470 
471  // Check exceptions.
472 #ifdef HAVE_INTREPID_DEBUG
473  cvals.resize(1,2,3);
474  INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
475  cvals.resize(3,2);
476  INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
477  cvals.resize(4,3);
478  INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
479 #endif
480  cvals.resize(4,2);
481  INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException ); nException--;
482  // Check if number of thrown exceptions matches the one we expect
483  if (throwCounter != nException) {
484  errorFlag++;
485  *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
486  }
487 
488  // Check mathematical correctness.
489  basis->getValues(bvals, cvals, OPERATOR_VALUE);
490  char buffer[120];
491  for (int i=0; i<bvals.dimension(0); i++) {
492  for (int j=0; j<bvals.dimension(1); j++) {
493  if ((i != j) && (std::abs(bvals(i,j) - 0.0) > INTREPID_TOL)) {
494  errorFlag++;
495  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);
496  *outStream << buffer;
497  }
498  else if ((i == j) && (std::abs(bvals(i,j) - 1.0) > INTREPID_TOL)) {
499  errorFlag++;
500  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);
501  *outStream << buffer;
502  }
503  }
504  }
505 
506  }
507  catch (const std::logic_error & err){
508  *outStream << err.what() << "\n\n";
509  errorFlag = -1000;
510  };
511 
512  if (errorFlag != 0)
513  std::cout << "End Result: TEST FAILED\n";
514  else
515  std::cout << "End Result: TEST PASSED\n";
516 
517  // reset format state of std::cout
518  std::cout.copyfmt(oldFormatState);
519 
520  return errorFlag;
521 }
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...