51 #include "Teuchos_oblackholestream.hpp"
52 #include "Teuchos_RCP.hpp"
53 #include "Teuchos_GlobalMPISession.hpp"
56 using namespace Intrepid;
58 #define INTREPID_TEST_COMMAND( S , throwCounter, nException ) \
64 catch (const std::logic_error & err) { \
66 *outStream << "Expected Error " << nException << " -------------------------------------------------------------\n"; \
67 *outStream << err.what() << '\n'; \
68 *outStream << "-------------------------------------------------------------------------------" << "\n\n"; \
72 int main(
int argc,
char *argv[]) {
74 Teuchos::GlobalMPISession mpiSession(&argc, &argv);
78 int iprint = argc - 1;
79 Teuchos::RCP<std::ostream> outStream;
80 Teuchos::oblackholestream bhs;
82 outStream = Teuchos::rcp(&std::cout,
false);
84 outStream = Teuchos::rcp(&bhs,
false);
87 Teuchos::oblackholestream oldFormatState;
88 oldFormatState.copyfmt(std::cout);
91 <<
"===============================================================================\n" \
93 <<
"| Unit Test (Basis_HGRAD_HEX_C2_FEM) |\n" \
95 <<
"| 1) Conversion of Dof tags into Dof ordinals and back |\n" \
96 <<
"| 2) Basis values for VALUE, GRAD, and Dk operators |\n" \
98 <<
"| Questions? Contact Pavel Bochev (pbboche@sandia.gov), |\n" \
99 <<
"| Robert Kirby (robert.c.kirby@ttu.edu), |\n" \
100 <<
"| Denis Ridzal (dridzal@sandia.gov), |\n" \
101 <<
"| Kara Peterson (kjpeter@sandia.gov). |\n" \
103 <<
"| Intrepid's website: http://trilinos.sandia.gov/packages/intrepid |\n" \
104 <<
"| Trilinos website: http://trilinos.sandia.gov |\n" \
106 <<
"===============================================================================\n"\
107 <<
"| TEST 1: Basis creation, exception testing |\n"\
108 <<
"===============================================================================\n";
113 shards::CellTopology line(shards::getCellTopologyData< shards::Line<> >());
115 PointTools::getLattice<double,FieldContainer<double> >(pts,line,deg);
123 int throwCounter = 0;
128 hexNodes(0, 0) = -1.0; hexNodes(0, 1) = -1.0; hexNodes(0, 2) = -1.0;
129 hexNodes(1, 0) = 0.0; hexNodes(1, 1) = -1.0; hexNodes(1, 2) = -1.0;
130 hexNodes(2, 0) = 1.0; hexNodes(2, 1) = -1.0; hexNodes(2, 2) = -1.0;
131 hexNodes(3, 0) = -1.0; hexNodes(3, 1) = 0.0; hexNodes(3, 2) = -1.0;
132 hexNodes(4, 0) = 0.0; hexNodes(4, 1) = 0.0; hexNodes(4, 2) = -1.0;
133 hexNodes(5, 0) = 1.0; hexNodes(5, 1) = 0.0; hexNodes(5, 2) = -1.0;
134 hexNodes(6, 0) = -1.0; hexNodes(6, 1) = 1.0; hexNodes(6, 2) = -1.0;
135 hexNodes(7, 0) = 0.0; hexNodes(7, 1) = 1.0; hexNodes(7, 2) = -1.0;
136 hexNodes(8, 0) = 1.0; hexNodes(8, 1) = 1.0; hexNodes(8, 2) = -1.0;
137 hexNodes(9, 0) = -1.0; hexNodes(9, 1) = -1.0; hexNodes(9, 2) = 0.0;
138 hexNodes(10, 0) = 0.0; hexNodes(10, 1) = -1.0; hexNodes(10, 2) = 0.0;
139 hexNodes(11, 0) = 1.0; hexNodes(11, 1) = -1.0; hexNodes(11, 2) = 0.0;
140 hexNodes(12, 0) = -1.0; hexNodes(12, 1) = 0.0; hexNodes(12, 2) = 0.0;
141 hexNodes(13, 0) = 0.0; hexNodes(13, 1) = 0.0; hexNodes(13, 2) = 0.0;
142 hexNodes(14, 0) = 1.0; hexNodes(14, 1) = 0.0; hexNodes(14, 2) = 0.0;
143 hexNodes(15, 0) = -1.0; hexNodes(15, 1) = 1.0; hexNodes(15, 2) = 0.0;
144 hexNodes(16, 0) = 0.0; hexNodes(16, 1) = 1.0; hexNodes(16, 2) = 0.0;
145 hexNodes(17, 0) = 1.0; hexNodes(17, 1) = 1.0; hexNodes(17, 2) = 0.0;
146 hexNodes(18, 0) = -1.0; hexNodes(18, 1) = -1.0; hexNodes(18, 2) = 1.0;
147 hexNodes(19, 0) = 0.0; hexNodes(19, 1) = -1.0; hexNodes(19, 2) = 1.0;
148 hexNodes(20, 0) = 1.0; hexNodes(20, 1) = -1.0; hexNodes(20, 2) = 1.0;
149 hexNodes(21, 0) = -1.0; hexNodes(21, 1) = 0.0; hexNodes(21, 2) = 1.0;
150 hexNodes(22, 0) = 0.0; hexNodes(22, 1) = 0.0; hexNodes(22, 2) = 1.0;
151 hexNodes(23, 0) = 1.0; hexNodes(23, 1) = 0.0; hexNodes(23, 2) = 1.0;
152 hexNodes(24, 0) = -1.0; hexNodes(24, 1) = 1.0; hexNodes(24, 2) = 1.0;
153 hexNodes(25, 0) = 0.0; hexNodes(25, 1) = 1.0; hexNodes(25, 2) = 1.0;
154 hexNodes(26, 0) = 1.0; hexNodes(26, 1) = 1.0; hexNodes(26, 2) = 1.0;
163 vals.
resize(hexBasis.getCardinality(), hexNodes.dimension(0), 4 );
164 INTREPID_TEST_COMMAND( hexBasis.getValues(vals, hexNodes, OPERATOR_CURL), throwCounter, nException );
168 vals.
resize(hexBasis.getCardinality(), hexNodes.dimension(0) );
169 INTREPID_TEST_COMMAND( hexBasis.getValues(vals, hexNodes, OPERATOR_DIV), throwCounter, nException );
174 INTREPID_TEST_COMMAND( hexBasis.getDofOrdinal(3,10,0), throwCounter, nException );
176 INTREPID_TEST_COMMAND( hexBasis.getDofOrdinal(1,2,1), throwCounter, nException );
178 INTREPID_TEST_COMMAND( hexBasis.getDofOrdinal(0,4,1), throwCounter, nException );
180 INTREPID_TEST_COMMAND( hexBasis.getDofTag(28), throwCounter, nException );
182 INTREPID_TEST_COMMAND( hexBasis.getDofTag(-1), throwCounter, nException );
184 #ifdef HAVE_INTREPID_DEBUG
188 INTREPID_TEST_COMMAND( hexBasis.getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException );
192 INTREPID_TEST_COMMAND( hexBasis.getValues(vals, badPoints2, OPERATOR_VALUE), throwCounter, nException );
196 INTREPID_TEST_COMMAND( hexBasis.getValues(badVals1, hexNodes, OPERATOR_VALUE), throwCounter, nException );
200 INTREPID_TEST_COMMAND( hexBasis.getValues(badVals2, hexNodes, OPERATOR_GRAD), throwCounter, nException );
203 INTREPID_TEST_COMMAND( hexBasis.getValues(badVals2, hexNodes, OPERATOR_D1), throwCounter, nException );
206 INTREPID_TEST_COMMAND( hexBasis.getValues(badVals2, hexNodes, OPERATOR_D2), throwCounter, nException );
210 INTREPID_TEST_COMMAND( hexBasis.getValues(badVals3, hexNodes, OPERATOR_VALUE), throwCounter, nException );
214 INTREPID_TEST_COMMAND( hexBasis.getValues(badVals4, hexNodes, OPERATOR_VALUE), throwCounter, nException );
217 FieldContainer<double> badVals5(hexBasis.getCardinality(), hexNodes.dimension(0), hexBasis.getBaseCellTopology().getDimension() - 1);
218 INTREPID_TEST_COMMAND( hexBasis.getValues(badVals5, hexNodes, OPERATOR_GRAD), throwCounter, nException );
222 INTREPID_TEST_COMMAND( hexBasis.getValues(badVals6, hexNodes, OPERATOR_D2), throwCounter, nException );
226 INTREPID_TEST_COMMAND( hexBasis.getValues(badVals7, hexNodes, OPERATOR_D3), throwCounter, nException );
230 catch (
const std::logic_error & err) {
231 *outStream <<
"UNEXPECTED ERROR !!! ----------------------------------------------------------\n";
232 *outStream << err.what() <<
'\n';
233 *outStream <<
"-------------------------------------------------------------------------------" <<
"\n\n";
239 if (throwCounter != nException) {
241 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
246 <<
"===============================================================================\n"\
247 <<
"| TEST 2: correctness of tag to enum and enum to tag lookups |\n"\
248 <<
"===============================================================================\n";
251 std::vector<std::vector<int> > allTags = hexBasis.getAllDofTags();
254 for (
unsigned i = 0; i < allTags.size(); i++) {
255 int bfOrd = hexBasis.getDofOrdinal(allTags[i][0], allTags[i][1], allTags[i][2]);
259 std::vector<int> myTag = hexBasis.getDofTag(bfOrd);
260 if( !( (myTag[0] == allTags[i][0]) &&
261 (myTag[1] == allTags[i][1]) &&
262 (myTag[2] == allTags[i][2]) &&
263 (myTag[3] == allTags[i][3]) ) ) {
265 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
266 *outStream <<
" getDofOrdinal( {"
267 << allTags[i][0] <<
", "
268 << allTags[i][1] <<
", "
269 << allTags[i][2] <<
", "
270 << allTags[i][3] <<
"}) = " << bfOrd <<
" but \n";
271 *outStream <<
" getDofTag(" << bfOrd <<
") = { "
275 << myTag[3] <<
"}\n";
280 for(
int bfOrd = 0; bfOrd < hexBasis.getCardinality(); bfOrd++) {
281 std::vector<int> myTag = hexBasis.getDofTag(bfOrd);
282 int myBfOrd = hexBasis.getDofOrdinal(myTag[0], myTag[1], myTag[2]);
283 if( bfOrd != myBfOrd) {
285 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
286 *outStream <<
" getDofTag(" << bfOrd <<
") = { "
290 << myTag[3] <<
"} but getDofOrdinal({"
294 << myTag[3] <<
"} ) = " << myBfOrd <<
"\n";
298 catch (
const std::logic_error & err){
299 *outStream << err.what() <<
"\n\n";
306 <<
"===============================================================================\n"\
307 <<
"| TEST 3: correctness of basis function values |\n"\
308 <<
"===============================================================================\n";
310 outStream -> precision(20);
313 double basisValues[] = {
314 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
315 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
316 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
317 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
318 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
319 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
320 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
321 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
322 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
323 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
324 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
325 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
326 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
327 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
328 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
329 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
330 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
331 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
332 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, 0, \
333 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, 0, \
334 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, 0, \
335 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, 0, \
336 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, 0, \
337 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, 0, \
338 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, 0, \
339 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000, 0, \
340 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.000 };
344 std::string fileName;
345 std::ifstream dataFile;
348 std::vector<double> basisGrads;
350 fileName =
"./testdata/HEX_C2_GradVals.dat";
351 dataFile.open(fileName.c_str());
352 TEUCHOS_TEST_FOR_EXCEPTION( !dataFile.good(), std::logic_error,
353 ">>> ERROR (HGRAD_HEX_C2/test01): could not open GRAD values data file, test aborted.");
354 while (!dataFile.eof() ){
357 std::getline(dataFile, line);
358 stringstream data_line(line);
359 while(data_line >> temp){
360 basisGrads.push_back(temp);
371 std::vector<double> basisD2;
372 fileName =
"./testdata/HEX_C2_D2Vals.dat";
373 dataFile.open(fileName.c_str());
374 TEUCHOS_TEST_FOR_EXCEPTION( !dataFile.good(), std::logic_error,
375 ">>> ERROR (HGRAD_HEX_C2/test01): could not open D2 values data file, test aborted.");
376 while (!dataFile.eof() ){
379 std::getline(dataFile, line);
380 stringstream data_line(line);
381 while(data_line >> temp){
382 basisD2.push_back(temp);
390 std::vector<double> basisD3;
392 fileName =
"./testdata/HEX_C2_D3Vals.dat";
393 dataFile.open(fileName.c_str());
394 TEUCHOS_TEST_FOR_EXCEPTION( !dataFile.good(), std::logic_error,
395 ">>> ERROR (HGRAD_HEX_C2/test01): could not open D3 values data file, test aborted.");
397 while (!dataFile.eof() ){
400 std::getline(dataFile, line);
401 stringstream data_line(line);
402 while(data_line >> temp){
403 basisD3.push_back(temp);
411 std::vector<double> basisD4;
413 fileName =
"./testdata/HEX_C2_D4Vals.dat";
414 dataFile.open(fileName.c_str());
415 TEUCHOS_TEST_FOR_EXCEPTION( !dataFile.good(), std::logic_error,
416 ">>> ERROR (HGRAD_HEX_C2/test01): could not open D4 values data file, test aborted.");
418 while (!dataFile.eof() ){
421 std::getline(dataFile, line);
422 stringstream data_line(line);
423 while(data_line >> temp){
424 basisD4.push_back(temp);
434 int numFields = hexBasis.getCardinality();
435 int numPoints = hexNodes.dimension(0);
436 int spaceDim = hexBasis.getBaseCellTopology().getDimension();
442 vals.
resize(numFields, numPoints);
443 hexBasis.getValues(vals, hexNodes, OPERATOR_VALUE);
444 for (
int i = 0; i < numFields; i++) {
445 for (
int j = 0; j < numPoints; j++) {
446 int l = i + j * numFields;
447 if (std::abs(vals(i,j) - basisValues[l]) > INTREPID_TOL) {
449 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
452 *outStream <<
" At multi-index { ";
453 *outStream << i <<
" ";*outStream << j <<
" ";
454 *outStream <<
"} computed value: " << vals(i,j)
455 <<
" but reference value: " << basisValues[l] <<
"\n";
461 vals.
resize(numFields, numPoints, spaceDim);
462 hexBasis.getValues(vals, hexNodes, OPERATOR_GRAD);
463 for (
int i = 0; i < numFields; i++) {
464 for (
int j = 0; j < numPoints; j++) {
465 for (
int k = 0; k < spaceDim; k++) {
468 int l = k + j * spaceDim + i * spaceDim * numPoints;
469 if (std::abs(vals(i,j,k) - basisGrads[l]) > INTREPID_TOL) {
471 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
474 *outStream <<
" At multi-index { ";
475 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
476 *outStream <<
"} computed grad component: " << vals(i,j,k)
477 <<
" but reference grad component: " << basisGrads[l] <<
"\n";
484 hexBasis.getValues(vals, hexNodes, OPERATOR_D1);
485 for (
int i = 0; i < numFields; i++) {
486 for (
int j = 0; j < numPoints; j++) {
487 for (
int k = 0; k < spaceDim; k++) {
490 int l = k + j * spaceDim + i * spaceDim * numPoints;
491 if (std::abs(vals(i,j,k) - basisGrads[l]) > INTREPID_TOL) {
493 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
496 *outStream <<
" At multi-index { ";
497 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
498 *outStream <<
"} computed D1 component: " << vals(i,j,k)
499 <<
" but reference D1 component: " << basisGrads[l] <<
"\n";
507 int D2cardinality = Intrepid::getDkCardinality(OPERATOR_D2, spaceDim);
508 vals.
resize(numFields, numPoints, D2cardinality);
509 hexBasis.getValues(vals, hexNodes, OPERATOR_D2);
510 for (
int i = 0; i < numFields; i++) {
511 for (
int j = 0; j < numPoints; j++) {
512 for (
int k = 0; k < D2cardinality; k++) {
515 int l = k + j * D2cardinality + i * D2cardinality * numPoints;
516 if (std::abs(vals(i,j,k) - basisD2[l]) > INTREPID_TOL) {
518 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
521 *outStream <<
" At multi-index { ";
522 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
523 *outStream <<
"} computed D2 component: " << vals(i,j,k)
524 <<
" but reference D2 component: " << basisD2[l] <<
"\n";
532 int D3cardinality = Intrepid::getDkCardinality(OPERATOR_D3, spaceDim);
533 vals.
resize(numFields, numPoints, D3cardinality);
534 hexBasis.getValues(vals, hexNodes, OPERATOR_D3);
536 for (
int i = 0; i < numFields; i++) {
537 for (
int j = 0; j < numPoints; j++) {
538 for (
int k = 0; k < D3cardinality; k++) {
541 int l = k + j * D3cardinality + i * D3cardinality * numPoints;
542 if (std::abs(vals(i,j,k) - basisD3[l]) > INTREPID_TOL) {
544 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
547 *outStream <<
" At multi-index { ";
548 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
549 *outStream <<
"} computed D3 component: " << vals(i,j,k)
550 <<
" but reference D3 component: " << basisD3[l] <<
"\n";
558 int D4cardinality = Intrepid::getDkCardinality(OPERATOR_D4, spaceDim);
559 vals.
resize(numFields, numPoints, D4cardinality);
560 hexBasis.getValues(vals, hexNodes, OPERATOR_D4);
561 for (
int i = 0; i < numFields; i++) {
562 for (
int j = 0; j < numPoints; j++) {
563 for (
int k = 0; k < D4cardinality; k++) {
566 int l = k + j * D4cardinality + i * D4cardinality * numPoints;
567 if (std::abs(vals(i,j,k) - basisD4[l]) > INTREPID_TOL) {
569 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
572 *outStream <<
" At multi-index { ";
573 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
574 *outStream <<
"} computed D4 component: " << vals(i,j,k)
575 <<
" but reference D4 component: " << basisD2[l] <<
"\n";
584 for(EOperator op = OPERATOR_D7; op < OPERATOR_MAX; op++) {
587 int DkCardin = Intrepid::getDkCardinality(op, spaceDim);
588 vals.
resize(numFields, numPoints, DkCardin);
590 hexBasis.getValues(vals, hexNodes, op);
591 for (
int i = 0; i < vals.
size(); i++) {
592 if (std::abs(vals[i]) > INTREPID_TOL) {
594 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
597 std::vector<int> myIndex;
599 int ord = Intrepid::getOperatorOrder(op);
600 *outStream <<
" At multi-index { ";
601 for(
int j = 0; j < vals.
rank(); j++) {
602 *outStream << myIndex[j] <<
" ";
604 *outStream <<
"} computed D"<< ord <<
" component: " << vals[i]
605 <<
" but reference D" << ord <<
" component: 0 \n";
611 catch (
const std::logic_error & err) {
612 *outStream << err.what() <<
"\n\n";
617 std::cout <<
"End Result: TEST FAILED\n";
619 std::cout <<
"End Result: TEST PASSED\n";
622 std::cout.copyfmt(oldFormatState);
int size() const
Returns size of the FieldContainer defined as the product of its dimensions.
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...
Header file for the Intrepid::HGRAD_HEX_Cn_FEM class.
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 2 on Hexahedron cell...
void resize(const int dim0)
Resizes FieldContainer to a rank-1 container with the specified dimension, initialized by 0...