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_HDIV_QUAD_In_FEM) |\n" \
95 <<
"| 1) Conversion of Dof tags into Dof ordinals and back |\n" \
96 <<
"| 2) Basis values for VALUE and DIV operators |\n" \
98 <<
"| Questions? Contact Pavel Bochev (pbboche@sandia.gov), |\n" \
99 <<
"| Denis Ridzal (dridzal@sandia.gov), |\n" \
100 <<
"| Kara Peterson (kjpeter@sandia.gov). |\n" \
102 <<
"| Intrepid's website: http://trilinos.sandia.gov/packages/intrepid |\n" \
103 <<
"| Trilinos website: http://trilinos.sandia.gov |\n" \
105 <<
"===============================================================================\n"\
106 <<
"| TEST 1: Basis creation, exception testing |\n"\
107 <<
"===============================================================================\n";
112 shards::CellTopology line(shards::getCellTopologyData< shards::Line<> >());
115 PointTools::getLattice<double,FieldContainer<double> >(closedPts,line,deg);
116 PointTools::getLattice<double,FieldContainer<double> >(openPts,line,deg+1,1);
123 int throwCounter = 0;
127 quadNodes(0,0) = -1.0; quadNodes(0,1) = -1.0;
128 quadNodes(1,0) = 0.0; quadNodes(1,1) = -1.0;
129 quadNodes(2,0) = 1.0; quadNodes(2,1) = -1.0;
130 quadNodes(3,0) = -1.0; quadNodes(3,1) = 0.0;
131 quadNodes(4,0) = 0.0; quadNodes(4,1) = 0.0;
132 quadNodes(5,0) = 1.0; quadNodes(5,1) = 0.0;
133 quadNodes(6,0) = -1.0; quadNodes(6,1) = 1.0;
134 quadNodes(7,0) = 0.0; quadNodes(7,1) = 1.0;
135 quadNodes(8,0) = 1.0; quadNodes(8,1) = 1.0;
142 int spaceDim = quadBasis.getBaseCellTopology().getDimension();
146 vals.
resize(quadBasis.getCardinality(), quadNodes.dimension(0), spaceDim );
147 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, quadNodes, OPERATOR_GRAD), throwCounter, nException );
150 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, quadNodes, OPERATOR_CURL), throwCounter, nException );
155 INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(3,0,0), throwCounter, nException );
157 INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(1,1,1), throwCounter, nException );
159 INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(0,4,1), throwCounter, nException );
161 INTREPID_TEST_COMMAND( quadBasis.getDofTag(12), throwCounter, nException );
163 INTREPID_TEST_COMMAND( quadBasis.getDofTag(-1), throwCounter, nException );
165 #ifdef HAVE_INTREPID_DEBUG
169 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException );
173 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, badPoints2, OPERATOR_VALUE), throwCounter, nException );
177 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals1, quadNodes, OPERATOR_VALUE), throwCounter, nException );
181 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_DIV), throwCounter, nException );
185 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals3, quadNodes, OPERATOR_VALUE), throwCounter, nException );
189 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals4, quadNodes, OPERATOR_DIV), throwCounter, nException );
193 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals5, quadNodes, OPERATOR_VALUE), throwCounter, nException );
197 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals6, quadNodes, OPERATOR_DIV), throwCounter, nException );
201 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals7, quadNodes, OPERATOR_VALUE), throwCounter, nException );
205 catch (
const std::logic_error & err) {
206 *outStream <<
"UNEXPECTED ERROR !!! ----------------------------------------------------------\n";
207 *outStream << err.what() <<
'\n';
208 *outStream <<
"-------------------------------------------------------------------------------" <<
"\n\n";
213 if (throwCounter != nException) {
215 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
220 <<
"===============================================================================\n"\
221 <<
"| TEST 2: correctness of tag to enum and enum to tag lookups |\n"\
222 <<
"===============================================================================\n";
225 std::vector<std::vector<int> > allTags = quadBasis.getAllDofTags();
228 for (
unsigned i = 0; i < allTags.size(); i++) {
229 int bfOrd = quadBasis.getDofOrdinal(allTags[i][0], allTags[i][1], allTags[i][2]);
231 std::vector<int> myTag = quadBasis.getDofTag(bfOrd);
232 if( !( (myTag[0] == allTags[i][0]) &&
233 (myTag[1] == allTags[i][1]) &&
234 (myTag[2] == allTags[i][2]) &&
235 (myTag[3] == allTags[i][3]) ) ) {
237 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
238 *outStream <<
" getDofOrdinal( {"
239 << allTags[i][0] <<
", "
240 << allTags[i][1] <<
", "
241 << allTags[i][2] <<
", "
242 << allTags[i][3] <<
"}) = " << bfOrd <<
" but \n";
243 *outStream <<
" getDofTag(" << bfOrd <<
") = { "
247 << myTag[3] <<
"}\n";
252 for(
int bfOrd = 0; bfOrd < quadBasis.getCardinality(); bfOrd++) {
253 std::vector<int> myTag = quadBasis.getDofTag(bfOrd);
254 int myBfOrd = quadBasis.getDofOrdinal(myTag[0], myTag[1], myTag[2]);
255 if( bfOrd != myBfOrd) {
257 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
258 *outStream <<
" getDofTag(" << bfOrd <<
") = { "
262 << myTag[3] <<
"} but getDofOrdinal({"
266 << myTag[3] <<
"} ) = " << myBfOrd <<
"\n";
270 catch (
const std::logic_error & err){
271 *outStream << err.what() <<
"\n\n";
277 <<
"===============================================================================\n"\
278 <<
"| TEST 3: correctness of basis function values |\n"\
279 <<
"===============================================================================\n";
281 outStream -> precision(20);
284 double basisValues[] = {
286 1.0, 0.0, 0.5, 0.0, 0.0, 0.0,
288 1.0, 0.0, 0.5, 0.0, 0.0, 0.0,
290 1.0, 0.0, 0.5, 0.0, 0.0, 0.0,
292 0.0, 0.0, 0.5, 0.0, 1.0, 0.0,
294 0.0, 0.0, 0.5, 0.0, 1.0, 0.0,
296 0.0, 0.0, 0.5, 0.0, 1.0, 0.0,
298 0.0, 1.0, 0.0, 1.0, 0.0, 1.0,
300 0.0, 0.5, 0.0, 0.5, 0.0, 0.5,
302 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
304 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
306 0.0, 0.5, 0.0, 0.5, 0.0, 0.5,
308 0.0, 1.0, 0.0, 1.0, 0.0, 1.0
312 double basisDivs[] = {
314 -0.5, -0.5, -0.5, -0.5, -0.5, -0.5, -0.5, -0.5, -0.5,
316 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5,
318 -0.5, -0.5, -0.5, -0.5, -0.5, -0.5, -0.5, -0.5, -0.5,
320 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5,
326 int numPoints = quadNodes.dimension(0);
327 int numFields = quadBasis.getCardinality();
328 int spaceDim = quadBasis.getBaseCellTopology().getDimension();
334 vals.
resize(numFields, numPoints, spaceDim);
335 quadBasis.getValues(vals, quadNodes, OPERATOR_VALUE);
336 for (
int i = 0; i < numFields; i++) {
337 for (
int j = 0; j < numPoints; j++) {
338 for (
int k = 0; k < spaceDim; k++) {
340 int l = i * numPoints * spaceDim + j * spaceDim + k;
341 if (std::abs(vals(i,j,k) - basisValues[l]) > INTREPID_TOL) {
343 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
346 *outStream <<
" At multi-index { ";
347 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
348 *outStream <<
"} computed value: " << vals(i,j,k)
349 <<
" but reference value: " << basisValues[l] <<
"\n";
356 vals.
resize(numFields, numPoints);
357 quadBasis.getValues(vals, quadNodes, OPERATOR_DIV);
358 for (
int i = 0; i < numFields; i++) {
359 for (
int j = 0; j < numPoints; j++) {
360 int l = j + i * numPoints;
361 if (std::abs(vals(i,j) - basisDivs[l]) > INTREPID_TOL) {
363 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
366 *outStream <<
" At multi-index { ";
367 *outStream << i <<
" ";*outStream << j <<
" ";
368 *outStream <<
"} computed divergence component: " << vals(i,j)
369 <<
" but reference divergence component: " << basisDivs[l] <<
"\n";
377 catch (
const std::logic_error & err) {
378 *outStream << err.what() <<
"\n\n";
383 std::cout <<
"End Result: TEST FAILED\n";
385 std::cout <<
"End Result: TEST PASSED\n";
388 std::cout.copyfmt(oldFormatState);
Header file for utility class to provide multidimensional containers.
Header file for the Intrepid::HDIV_QUAD_In_FEM class.
void resize(const int dim0)
Resizes FieldContainer to a rank-1 container with the specified dimension, initialized by 0...
Implementation of the default H(div)-compatible FEM basis of degree 1 on Quadrilateral cell...