50 #include "Teuchos_oblackholestream.hpp"
51 #include "Teuchos_RCP.hpp"
52 #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_QUAD_Cn_FEM) |\n" \
95 <<
"| 1) Conversion of Dof tags into Dof ordinals and back |\n" \
96 <<
"| 2) Basis values for VALUE, GRAD, CURL, 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);
122 int throwCounter = 0;
126 quadNodes(0,0) = -1.0; quadNodes(0,1) = -1.0;
127 quadNodes(1,0) = 1.0; quadNodes(1,1) = -1.0;
128 quadNodes(2,0) = 1.0; quadNodes(2,1) = 1.0;
129 quadNodes(3,0) = -1.0; quadNodes(3,1) = 1.0;
131 quadNodes(4,0) = 0.0; quadNodes(4,1) = -1.0;
132 quadNodes(5,0) = 1.0; quadNodes(5,1) = 0.0;
133 quadNodes(6,0) = 0.0; quadNodes(6,1) = 1.0;
134 quadNodes(7,0) = -1.0; quadNodes(7,1) = 0.0;
136 quadNodes(8,0) = 0.0; quadNodes(8,1) = 0.0;
137 quadNodes(9,0) =1./3.; quadNodes(9,1) =-3./5.;
145 vals.
resize(quadBasis.getCardinality(), quadNodes.dimension(0));
146 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, quadNodes, OPERATOR_DIV), throwCounter, nException );
151 INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(3,0,0), throwCounter, nException );
153 INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(1,1,1), throwCounter, nException );
155 INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(0,4,0), throwCounter, nException );
157 INTREPID_TEST_COMMAND( quadBasis.getDofTag(10), throwCounter, nException );
159 INTREPID_TEST_COMMAND( quadBasis.getDofTag(-1), throwCounter, nException );
161 #ifdef HAVE_INTREPID_DEBUG
165 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException );
169 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, badPoints2, OPERATOR_VALUE), throwCounter, nException );
173 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals1, quadNodes, OPERATOR_VALUE), throwCounter, nException );
177 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_GRAD), throwCounter, nException );
180 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_CURL), throwCounter, nException );
183 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_D2), throwCounter, nException );
187 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals3, quadNodes, OPERATOR_VALUE), throwCounter, nException );
191 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals4, quadNodes, OPERATOR_VALUE), throwCounter, nException );
194 FieldContainer<double> badVals5(quadBasis.getCardinality(), quadNodes.dimension(0), quadBasis.getBaseCellTopology().getDimension() + 1);
195 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals5, quadNodes, OPERATOR_GRAD), throwCounter, nException );
199 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals6, quadNodes, OPERATOR_D2), throwCounter, nException );
202 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals6, quadNodes, OPERATOR_D3), throwCounter, nException );
206 catch (
const std::logic_error & err) {
207 *outStream <<
"UNEXPECTED ERROR !!! ----------------------------------------------------------\n";
208 *outStream << err.what() <<
'\n';
209 *outStream <<
"-------------------------------------------------------------------------------" <<
"\n\n";
214 if (throwCounter != nException) {
216 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
221 <<
"===============================================================================\n"\
222 <<
"| TEST 2: correctness of tag to enum and enum to tag lookups |\n"\
223 <<
"===============================================================================\n";
226 std::vector<std::vector<int> > allTags = quadBasis.getAllDofTags();
229 for (
unsigned i = 0; i < allTags.size(); i++) {
230 int bfOrd = quadBasis.getDofOrdinal(allTags[i][0], allTags[i][1], allTags[i][2]);
232 std::vector<int> myTag = quadBasis.getDofTag(bfOrd);
233 if( !( (myTag[0] == allTags[i][0]) &&
234 (myTag[1] == allTags[i][1]) &&
235 (myTag[2] == allTags[i][2]) &&
236 (myTag[3] == allTags[i][3]) ) ) {
238 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
239 *outStream <<
" getDofOrdinal( {"
240 << allTags[i][0] <<
", "
241 << allTags[i][1] <<
", "
242 << allTags[i][2] <<
", "
243 << allTags[i][3] <<
"}) = " << bfOrd <<
" but \n";
244 *outStream <<
" getDofTag(" << bfOrd <<
") = { "
248 << myTag[3] <<
"}\n";
253 for(
int bfOrd = 0; bfOrd < quadBasis.getCardinality(); bfOrd++) {
254 std::vector<int> myTag = quadBasis.getDofTag(bfOrd);
255 int myBfOrd = quadBasis.getDofOrdinal(myTag[0], myTag[1], myTag[2]);
256 if( bfOrd != myBfOrd) {
258 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
259 *outStream <<
" getDofTag(" << bfOrd <<
") = { "
263 << myTag[3] <<
"} but getDofOrdinal({"
267 << myTag[3] <<
"} ) = " << myBfOrd <<
"\n";
271 catch (
const std::logic_error & err){
272 *outStream << err.what() <<
"\n\n";
278 <<
"===============================================================================\n"\
279 <<
"| TEST 3: correctness of basis function values |\n"\
280 <<
"===============================================================================\n";
282 outStream -> precision(20);
285 double basisValues[] = {
286 1, 0, 0, 0, 0, 0, 0, 0, 0, -0.05333333333333334, \
287 0, 0, 0, 0, 1, 0, 0, 0, 0, 0.4266666666666667, \
288 0, 1, 0, 0, 0, 0, 0, 0, 0, 0.1066666666666667, \
289 0, 0, 0, 0, 0, 0, 0, 1, 0, -0.07111111111111112 , \
290 0, 0, 0, 0, 0, 0, 0, 0, 1, 0.5688888888888890, \
291 0, 0, 0, 0, 0, 1, 0, 0, 0, 0.1422222222222222 ,\
292 0, 0, 0, 1, 0, 0, 0, 0, 0, 0.01333333333333333, \
293 0, 0, 0, 0, 0, 0, 1, 0, 0, -0.1066666666666667, \
294 0, 0, 1, 0, 0, 0, 0, 0, 0, -0.02666666666666666 };
302 double basisGrads[] = {
304 -1.500000000000000, -1.500000000000000, 0.5000000000000000, 0, 0, 0, 0, 0.5000000000000000, -0.5000000000000000, 0, \
305 0, 0, 0, 0, 0, -0.5000000000000000, 0, 0, -0.08000000000000002, 0.1222222222222222, \
307 2.000000000000000, 0, -2.000000000000000, 0, 0, 0, 0, 0, 0, -1.500000000000000, \
308 0, 0, 0, 0.5000000000000000, 0, 0, 0, -0.5000000000000000, -0.3199999999999999, -0.9777777777777779, \
310 -0.5000000000000000, 0, 1.500000000000000, -1.500000000000000, 0, 0.5000000000000000, 0, 0, 0.5000000000000000, 0, \
311 0, -0.5000000000000000, 0, 0, 0, 0, 0, 0, 0.3999999999999999, -0.2444444444444444, \
313 0, 2.0, 0, 0, 0, 0, 0, -2.000000000000000, 0, 0, \
314 0.5000000000000000, 0, 0, 0, -1.50, 0, -0.50, 0, -0.1066666666666667, -0.1333333333333333, \
316 0, 0, 0, 0, 0, 0, 0, 0, 0, 2.0,\
317 -2.00, 0, 0, -2.0, 2.0, 0, 0, 0, -0.4266666666666667, 1.066666666666667 , \
319 0, 0, 0, 2.000000000000000, 0, -2.000000000000000, 0, 0, 0, 0, \
320 1.5, 0, 0, 0, -0.5, 0, 0.5000000000000000, 0, 0.5333333333333334, 0.2666666666666666 , \
322 0, -0.5000000000000000, 0, 0, 0.5000000000000000, 0, -1.500000000000000, 1.500000000000000, 0, 0, \
323 0, 0, -0.5000000000000000, 0, 0, 0.5000000000000000, 0, 0, 0.02000000000000000, 0.01111111111111112 , \
325 0, 0, 0, 0, -2.0, 0, 2.0, 0, 0, -0.50, \
326 0, 0, 0, 1.5, 0, 0, 0, 0.5000000000000000, 0.07999999999999997, -0.08888888888888888, \
328 0, 0, 0, -0.5000000000000000, 1.500000000000000, 1.500000000000000, -0.5000000000000000, 0, 0, 0, \
329 0, 0.5000000000000000, 0.5000000000000000, 0, 0, 0, 0, 0, -0.09999999999999998, -0.02222222222222221 \
337 1.0, 2.25, 1.0, 1.0, -0.75, 0, 0, 0.25, 0, 0, -0.75, 1.0, 1.0, 0.75, 0, 0, -0.25, 0, 0, -0.25, 0, 0, 0.75, 1.0, 0, 0.25, 0, 0.48, 0.1833333333333334, -0.1111111111111111,
339 -2.0, -3.0, 0, -2.0, 3.0, 0, 0, -1.0, 0, \
340 0, 1.0, 0, -2.0, 0, 1.0, 0, \
341 1.0, 0, 0, 0, 1.0, 0, -1.0, \
342 0, 0, 0, 1.0, -0.96, 0.7333333333333332, \
343 0.8888888888888890, \
345 1.0, 0.75, 0, 1.0, -2.25, 1.0, 0, 0.75, 1.0, 0, -0.25, 0, \
346 1.0, -0.75, 0, 0, -0.75, 1.0, 0, 0.25, 0, 0, 0.25, \
347 0, 0, -0.25, 0, 0.48, -0.9166666666666666, 0.2222222222222222,
349 0, -3.0, -2.0, 0, 1.0, 0, 0, -1.0, 0, 0, 3.0, \
350 -2.0, 0, -1.0, 0, 1.0, 0, 0, 0, 1.0, 0, 1.0, 0, -2.0, \
351 1.0, 0, 0, 0.6400000000000001, -0.2000000000000001, 0.2222222222222222, \
353 0, 4.0, 0, 0, -4.0, 0, 0, 4.0, 0, 0, -4.0, 0, 0, 0, \
354 -2.0, -2.0, 0, 0, 0, 0, -2.0, -2.0, 0, 0, -2.0, 0, \
355 -2.0, -1.280000000000000, -0.7999999999999998, -1.777777777777778 ,
357 0, -1.0, 0, 0, 3.0, -2.0, 0, -3.0, -2.0, 0, \
358 1.0, 0, 0, 1.0, 0, 1.0, 0, -2.0, 0, -1.0, 0, 1.0, 0, \
359 0, 1.0, 0, 0, 0.6400000000000001, 1.0, -0.4444444444444444, \
361 0, 0.75, 1.0, 0, -0.25, 0, 1.0, 0.75, 0, 1.0, -2.25, 1.0, 0, \
362 0.25, 0, 0, 0.25, 0, 1.0, -0.75, 0, 0, -0.75, 1.0, 0, \
363 -0.25, 0, -0.12, 0.01666666666666666, -0.1111111111111111, \
365 0, -1.0, 0, 0, 1.0, 0, -2.0, -3.0, 0, -2.0, 3.0, 0, 0, 0, 1.0, 0, -1.0, \
366 0, -2.0, 0, 1.0, 0, 1.0, 0, \
367 0, 0, 1.0, 0.24, 0.06666666666666665, 0.8888888888888890, \
369 0, 0.25, 0, 0, -0.75, 1.0, 1.0, 2.25, 1.0, 1.0, \
370 -0.75, 0, 0, -0.25, 0, 0, 0.75, 1.0, 1.0, \
371 0.75, 0, 0, -0.25, 0, 0, 0.25, 0, -0.12, -0.08333333333333331, 0.2222222222222222 \
376 0, -1.5, -1.5, 0, 0, -1.5, 0.5, 0, 0, 0.5,
377 0.5, 0, 0, 0.5, -1.5, 0, 0, -1.5, -0.5, 0,
378 0, -0.5, 0.5, 0, 0, 0.5, -0.5, 0, 0, -0.5,
379 -1.5, 0, 0, -0.5, -0.5, 0, 0, -1.1, -0.1666666666666667, 0,
381 0, 3.0, 2.0, 0, 0, 3.0, -2.0, 0, 0, -1.0,
382 -2.0, 0, 0, -1.0, 2.0, 0, 0, 3.0, 0, 0,
383 0, 1.0, -2.0, 0, 0, -1.0, 0, 0, 0, 1.0,
384 2.0, 0, 0, 1.0, 0, 0, 0, 2.2, -0.6666666666666665, 0,
386 0, -1.5, -0.5, 0, 0, -1.5, 1.5, 0, 0, 0.5,
387 1.5, 0, 0, 0.5, -0.5, 0, 0, -1.5, 0.5, 0,
388 0, -0.5, 1.5, 0, 0, 0.5, 0.5, 0, 0, -0.5,
389 -0.5, 0, 0, -0.5, 0.5, 0, 0, -1.1, 0.8333333333333333, 0,
391 0, 2.0, 3.0, 0, 0, 2.0, -1.0, 0, 0, -2.0,
392 -1.0, 0, 0, -2.0, 3.0, 0, 0, 2.0, 1.0, 0,
393 0, 0, -1.0, 0, 0, -2.0, 1.0, 0, 0, 0,
394 3.0, 0, 0, 0, 1.0, 0, 0, 1.2, 0.3333333333333334, 0,
396 0, -4.0, -4.0, 0, 0, -4.0, 4.0, 0, 0, 4.0,
397 4.0, 0, 0, 4.0, -4.0, 0, 0, -4.0, 0, 0,
398 0, 0, 4.0, 0, 0, 4.0, 0, 0, 0, 0,
399 -4.0, 0, 0, 0, 0, 0, 0, -2.40, 1.333333333333333, 0,
401 0, 2.0, 1.0, 0, 0, 2.0, -3.0, 0, 0, -2.0,
402 -3.0, 0, 0, -2.0, 1.0, 0, 0, 2.0, -1.0, 0,
403 0, 0, -3.0, 0, 0, -2.0, -1.0, 0, 0, 0,
404 1.0, 0, 0, 0, -1.0, 0, 0, 1.2, -1.666666666666667, 0 ,
406 0, -0.5, -1.5, 0, 0, -0.5, 0.5, 0, 0, 1.5,
407 0.5, 0, 0, 1.5, -1.5, 0, 0, -0.5, -0.5, 0,
408 0, 0.5, 0.5, 0, 0, 1.5, -0.5, 0, 0, 0.5,
409 -1.5, 0, 0, 0.5, -0.5, 0, 0, -0.09999999999999998, -0.1666666666666667, 0,
411 0, 1.0, 2.0, 0, 0, 1.0, -2.0, 0, 0, -3.0,
412 -2.0, 0, 0, -3.0, 2.0, 0, 0, 1.0, 0, 0,
413 0, -1.0, -2.0, 0, 0, -3.0, 0, 0, 0, -1.0,
414 2.0, 0, 0, -1.0, 0, 0, 0, 0.2, -0.6666666666666665, 0,
416 0, -0.5, -0.5, 0, 0, -0.5, 1.5, 0, 0, 1.5,
417 1.5, 0, 0, 1.5, -0.5, 0, 0, -0.5, 0.5, 0,
418 0, 0.5, 1.5, 0, 0, 1.5, 0.5, 0, 0, 0.5,
419 -0.5, 0, 0, 0.5, 0.5, 0, 0, -0.09999999999999998, 0.8333333333333333, 0
423 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
424 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
425 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
426 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
427 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
429 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
430 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
431 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
432 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
433 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
435 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
436 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
437 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
438 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
439 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
441 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
442 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
443 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
444 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
445 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
447 0, 0, 4.0, 0, 0, 0, 0, 4.0, 0, 0,
448 0, 0, 4.0, 0, 0, 0, 0, 4.0, 0, 0,
449 0, 0, 4.0, 0, 0, 0, 0, 4.0, 0, 0,
450 0, 0, 4.0, 0, 0, 0, 0, 4.0, 0, 0,
451 0, 0, 4.0, 0, 0, 0, 0, 4.0, 0, 0,
453 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
454 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
455 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
456 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
457 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
459 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
460 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
461 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
462 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
463 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
465 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
466 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
467 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
468 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
469 0, 0, -2.0, 0, 0, 0, 0, -2.0, 0, 0,
471 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
472 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
473 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
474 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0,
475 0, 0, 1.0, 0, 0, 0, 0, 1.0, 0, 0
481 int numFields = quadBasis.getCardinality();
482 int numPoints = quadNodes.dimension(0);
483 int spaceDim = quadBasis.getBaseCellTopology().getDimension();
489 vals.
resize(numFields, numPoints);
490 quadBasis.getValues(vals, quadNodes, OPERATOR_VALUE);
491 for (
int i = 0; i < numFields; i++) {
492 for (
int j = 0; j < numPoints; j++) {
495 int l = j + i * numPoints;
496 if (std::abs(vals(i,j) - basisValues[l]) > INTREPID_TOL) {
498 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
501 *outStream <<
" At multi-index { ";
502 *outStream << i <<
" ";*outStream << j <<
" ";
503 *outStream <<
"} computed value: " << vals(i,j)
504 <<
" but reference value: " << basisValues[l] <<
"\n";
510 vals.
resize(numFields, numPoints, spaceDim);
511 quadBasis.getValues(vals, quadNodes, OPERATOR_GRAD);
512 for (
int i = 0; i < numFields; i++) {
513 for (
int j = 0; j < numPoints; j++) {
514 for (
int k = 0; k < spaceDim; k++) {
517 int l = k + j * spaceDim + i * spaceDim * numPoints;
518 if (std::abs(vals(i,j,k) - basisGrads[l]) > INTREPID_TOL) {
520 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
523 *outStream <<
" At multi-index { ";
524 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
525 *outStream <<
"} computed grad component: " << vals(i,j,k)
526 <<
" but reference grad component: " << basisGrads[l] <<
"\n";
534 quadBasis.getValues(vals, quadNodes, OPERATOR_D1);
535 for (
int i = 0; i < numFields; i++) {
536 for (
int j = 0; j < numPoints; j++) {
537 for (
int k = 0; k < spaceDim; k++) {
540 int l = k + j * spaceDim + i * spaceDim * numPoints;
541 if (std::abs(vals(i,j,k) - basisGrads[l]) > INTREPID_TOL) {
543 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
546 *outStream <<
" At multi-index { ";
547 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
548 *outStream <<
"} computed D1 component: " << vals(i,j,k)
549 <<
" but reference D1 component: " << basisGrads[l] <<
"\n";
557 vals.
resize(numFields, numPoints, spaceDim);
558 quadBasis.getValues(vals, quadNodes, OPERATOR_CURL);
559 for (
int i = 0; i < numFields; i++) {
560 for (
int j = 0; j < numPoints; j++) {
562 int curl_0 = 1 + j * spaceDim + i * spaceDim * numPoints;
563 int curl_1 = 0 + j * spaceDim + i * spaceDim * numPoints;
565 double curl_value_0 = basisGrads[curl_0];
566 double curl_value_1 =-basisGrads[curl_1];
567 if (std::abs(vals(i,j,0) - curl_value_0) > INTREPID_TOL) {
569 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
571 *outStream <<
" At multi-index { ";
572 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << 0 <<
" ";
573 *outStream <<
"} computed curl component: " << vals(i,j,0)
574 <<
" but reference curl component: " << curl_value_0 <<
"\n";
576 if (std::abs(vals(i,j,1) - curl_value_1) > INTREPID_TOL) {
578 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
580 *outStream <<
" At multi-index { ";
581 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << 1 <<
" ";
582 *outStream <<
"} computed curl component: " << vals(i,j,1)
583 <<
" but reference curl component: " << curl_value_1 <<
"\n";
589 int D2cardinality = Intrepid::getDkCardinality(OPERATOR_D2, spaceDim);
590 vals.
resize(numFields, numPoints, D2cardinality);
591 quadBasis.getValues(vals, quadNodes, OPERATOR_D2);
592 for (
int i = 0; i < numFields; i++) {
593 for (
int j = 0; j < numPoints; j++) {
594 for (
int k = 0; k < D2cardinality; k++) {
597 int l = k + j * D2cardinality + i * D2cardinality * numPoints;
598 if (std::abs(vals(i,j,k) - basisD2[l]) > INTREPID_TOL) {
600 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
603 *outStream <<
" At multi-index { ";
604 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
605 *outStream <<
"} computed D2 component: " << vals(i,j,k)
606 <<
" but reference D2 component: " << basisD2[l] <<
"\n";
614 int D3cardinality = Intrepid::getDkCardinality(OPERATOR_D3, spaceDim);
615 vals.
resize(numFields, numPoints, D3cardinality);
616 quadBasis.getValues(vals, quadNodes, OPERATOR_D3);
617 for (
int i = 0; i < numFields; i++) {
618 for (
int j = 0; j < numPoints; j++) {
619 for (
int k = 0; k < D3cardinality; k++) {
622 int l = k + j * D3cardinality + i * D3cardinality * numPoints;
623 if (std::abs(vals(i,j,k) - basisD3[l]) > INTREPID_TOL) {
625 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
628 *outStream <<
" At multi-index { ";
629 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
630 *outStream <<
"} computed D3 component: " << vals(i,j,k)
631 <<
" but reference D3 component: " << basisD2[l] <<
"\n";
638 int D4cardinality = Intrepid::getDkCardinality(OPERATOR_D4, spaceDim);
639 vals.
resize(numFields, numPoints, D4cardinality);
640 quadBasis.getValues(vals, quadNodes, OPERATOR_D4);
641 for (
int i = 0; i < numFields; i++) {
642 for (
int j = 0; j < numPoints; j++) {
643 for (
int k = 0; k < D4cardinality; k++) {
646 int l = k + j * D4cardinality + i * D4cardinality * numPoints;
647 if (std::abs(vals(i,j,k) - basisD4[l]) > INTREPID_TOL) {
649 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
652 *outStream <<
" At multi-index { ";
653 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
654 *outStream <<
"} computed D4 component: " << vals(i,j,k)
655 <<
" but reference D4 component: " << basisD4[l] <<
"\n";
663 for(EOperator op = OPERATOR_D5; op <= OPERATOR_D6; op++) {
666 int DkCardin = Intrepid::getDkCardinality(op, spaceDim);
667 vals.
resize(numFields, numPoints, DkCardin);
669 quadBasis.getValues(vals, quadNodes, op);
670 for (
int i = 0; i < vals.
size(); i++) {
671 if (std::abs(vals[i]) > INTREPID_TOL) {
673 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
676 std::vector<int> myIndex;
678 int ord = Intrepid::getOperatorOrder(op);
679 *outStream <<
" At multi-index { ";
680 for(
int j = 0; j < vals.
rank(); j++) {
681 *outStream << myIndex[j] <<
" ";
683 *outStream <<
"} computed D"<< ord <<
" component: " << vals[i]
684 <<
" but reference D" << ord <<
" component: 0 \n";
691 catch (
const std::logic_error & err) {
692 *outStream << err.what() <<
"\n\n";
697 std::cout <<
"End Result: TEST FAILED\n";
699 std::cout <<
"End Result: TEST PASSED\n";
702 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...
int rank() const
Return rank of the FieldContainer = number of indices used to tag the multi-indexed value...
void resize(const int dim0)
Resizes FieldContainer to a rank-1 container with the specified dimension, initialized by 0...
Header file for the Intrepid::HGRAD_QUAD_Cn_FEM class.