50#include "Teuchos_oblackholestream.hpp"
51#include "Teuchos_RCP.hpp"
52#include "Teuchos_GlobalMPISession.hpp"
55using namespace Intrepid;
57#define INTREPID_TEST_COMMAND( S , throwCounter, nException ) \
63 catch (const std::logic_error & err) { \
65 *outStream << "Expected Error " << nException << " -------------------------------------------------------------\n"; \
66 *outStream << err.what() << '\n'; \
67 *outStream << "-------------------------------------------------------------------------------" << "\n\n"; \
71int main(
int argc,
char *argv[]) {
73 Teuchos::GlobalMPISession mpiSession(&argc, &argv);
77 int iprint = argc - 1;
78 Teuchos::RCP<std::ostream> outStream;
79 Teuchos::oblackholestream bhs;
81 outStream = Teuchos::rcp(&std::cout,
false);
83 outStream = Teuchos::rcp(&bhs,
false);
86 Teuchos::oblackholestream oldFormatState;
87 oldFormatState.copyfmt(std::cout);
90 <<
"===============================================================================\n" \
92 <<
"| Unit Test (Basis_HDIV_QUAD_I1_FEM) |\n" \
94 <<
"| 1) Conversion of Dof tags into Dof ordinals and back |\n" \
95 <<
"| 2) Basis values for VALUE and DIV operators |\n" \
97 <<
"| Questions? Contact Pavel Bochev (pbboche@sandia.gov), |\n" \
98 <<
"| Denis Ridzal (dridzal@sandia.gov), |\n" \
99 <<
"| Kara Peterson (kjpeter@sandia.gov). |\n" \
101 <<
"| Intrepid's website: http://trilinos.sandia.gov/packages/intrepid |\n" \
102 <<
"| Trilinos website: http://trilinos.sandia.gov |\n" \
104 <<
"===============================================================================\n"\
105 <<
"| TEST 1: Basis creation, exception testing |\n"\
106 <<
"===============================================================================\n";
109 Basis_HDIV_QUAD_I1_FEM<double, FieldContainer<double> > quadBasis;
114 int throwCounter = 0;
117 FieldContainer<double> quadNodes(9, 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;
123 quadNodes(4,0) = 0.0; quadNodes(4,1) = 0.0;
124 quadNodes(5,0) = 0.0; quadNodes(5,1) = -0.5;
125 quadNodes(6,0) = 0.0; quadNodes(6,1) = 0.5;
126 quadNodes(7,0) = -0.5; quadNodes(7,1) = 0.0;
127 quadNodes(8,0) = 0.5; quadNodes(8,1) = 0.0;
130 FieldContainer<double> vals;
134 int spaceDim = quadBasis.getBaseCellTopology().getDimension();
138 vals.resize(quadBasis.getCardinality(), quadNodes.dimension(0), spaceDim );
139 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, quadNodes, OPERATOR_GRAD), throwCounter, nException );
142 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, quadNodes, OPERATOR_CURL), throwCounter, nException );
147 INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(3,0,0), throwCounter, nException );
149 INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(1,1,1), throwCounter, nException );
151 INTREPID_TEST_COMMAND( quadBasis.getDofOrdinal(0,4,1), throwCounter, nException );
153 INTREPID_TEST_COMMAND( quadBasis.getDofTag(12), throwCounter, nException );
155 INTREPID_TEST_COMMAND( quadBasis.getDofTag(-1), throwCounter, nException );
157#ifdef HAVE_INTREPID_DEBUG
160 FieldContainer<double> badPoints1(4, 5, 3);
161 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException );
164 FieldContainer<double> badPoints2(4, spaceDim + 1);
165 INTREPID_TEST_COMMAND( quadBasis.getValues(vals, badPoints2, OPERATOR_VALUE), throwCounter, nException );
168 FieldContainer<double> badVals1(4, 5);
169 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals1, quadNodes, OPERATOR_VALUE), throwCounter, nException );
172 FieldContainer<double> badVals2(4, 5, spaceDim);
173 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals2, quadNodes, OPERATOR_DIV), throwCounter, nException );
176 FieldContainer<double> badVals3(quadBasis.getCardinality() + 1, quadNodes.dimension(0), spaceDim);
177 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals3, quadNodes, OPERATOR_VALUE), throwCounter, nException );
180 FieldContainer<double> badVals4(quadBasis.getCardinality() + 1, quadNodes.dimension(0));
181 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals4, quadNodes, OPERATOR_DIV), throwCounter, nException );
184 FieldContainer<double> badVals5(quadBasis.getCardinality(), quadNodes.dimension(0) + 1, spaceDim);
185 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals5, quadNodes, OPERATOR_VALUE), throwCounter, nException );
188 FieldContainer<double> badVals6(quadBasis.getCardinality(), quadNodes.dimension(0) + 1);
189 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals6, quadNodes, OPERATOR_DIV), throwCounter, nException );
192 FieldContainer<double> badVals7(quadBasis.getCardinality(), quadNodes.dimension(0), spaceDim + 1);
193 INTREPID_TEST_COMMAND( quadBasis.getValues(badVals7, quadNodes, OPERATOR_VALUE), throwCounter, nException );
197 catch (
const std::logic_error & err) {
198 *outStream <<
"UNEXPECTED ERROR !!! ----------------------------------------------------------\n";
199 *outStream << err.what() <<
'\n';
200 *outStream <<
"-------------------------------------------------------------------------------" <<
"\n\n";
205 if (throwCounter != nException) {
207 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
212 <<
"===============================================================================\n"\
213 <<
"| TEST 2: correctness of tag to enum and enum to tag lookups |\n"\
214 <<
"===============================================================================\n";
217 std::vector<std::vector<int> > allTags = quadBasis.getAllDofTags();
220 for (
unsigned i = 0; i < allTags.size(); i++) {
221 int bfOrd = quadBasis.getDofOrdinal(allTags[i][0], allTags[i][1], allTags[i][2]);
223 std::vector<int> myTag = quadBasis.getDofTag(bfOrd);
224 if( !( (myTag[0] == allTags[i][0]) &&
225 (myTag[1] == allTags[i][1]) &&
226 (myTag[2] == allTags[i][2]) &&
227 (myTag[3] == allTags[i][3]) ) ) {
229 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
230 *outStream <<
" getDofOrdinal( {"
231 << allTags[i][0] <<
", "
232 << allTags[i][1] <<
", "
233 << allTags[i][2] <<
", "
234 << allTags[i][3] <<
"}) = " << bfOrd <<
" but \n";
235 *outStream <<
" getDofTag(" << bfOrd <<
") = { "
239 << myTag[3] <<
"}\n";
244 for(
int bfOrd = 0; bfOrd < quadBasis.getCardinality(); bfOrd++) {
245 std::vector<int> myTag = quadBasis.getDofTag(bfOrd);
246 int myBfOrd = quadBasis.getDofOrdinal(myTag[0], myTag[1], myTag[2]);
247 if( bfOrd != myBfOrd) {
249 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
250 *outStream <<
" getDofTag(" << bfOrd <<
") = { "
254 << myTag[3] <<
"} but getDofOrdinal({"
258 << myTag[3] <<
"} ) = " << myBfOrd <<
"\n";
262 catch (
const std::logic_error & err){
263 *outStream << err.what() <<
"\n\n";
269 <<
"===============================================================================\n"\
270 <<
"| TEST 3: correctness of basis function values |\n"\
271 <<
"===============================================================================\n";
273 outStream -> precision(20);
276 double basisValues[] = {
277 0, -0.500000, 0, 0, 0, 0, -0.500000, 0, 0, -0.500000, 0.500000, 0, 0, \
278 0, 0, 0, 0, 0, 0.500000, 0, 0, 0.500000, 0, 0, 0, 0, 0, 0, 0, \
279 0.500000, -0.500000, 0, 0, -0.250000, 0.250000, 0, 0, 0.250000, \
280 -0.250000, 0, 0, -0.375000, 0.250000, 0, 0, 0.125000, -0.250000, 0, \
281 0, -0.125000, 0.250000, 0, 0, 0.375000, -0.250000, 0, 0, -0.250000, \
282 0.125000, 0, 0, 0.250000, -0.375000, 0, 0, -0.250000, 0.375000, 0, 0, \
283 0.250000, -0.125000, 0 };
286 double basisDivs[] = {
287 0.25, 0.25, 0.25, 0.25,
288 0.25, 0.25, 0.25, 0.25,
289 0.25, 0.25, 0.25, 0.25,
290 0.25, 0.25, 0.25, 0.25,
291 0.25, 0.25, 0.25, 0.25,
292 0.25, 0.25, 0.25, 0.25,
293 0.25, 0.25, 0.25, 0.25,
294 0.25, 0.25, 0.25, 0.25,
295 0.25, 0.25, 0.25, 0.25,
301 int numPoints = quadNodes.dimension(0);
302 int numFields = quadBasis.getCardinality();
303 int spaceDim = quadBasis.getBaseCellTopology().getDimension();
306 FieldContainer<double> vals;
309 vals.resize(numFields, numPoints, spaceDim);
310 quadBasis.getValues(vals, quadNodes, OPERATOR_VALUE);
311 for (
int i = 0; i < numFields; i++) {
312 for (
int j = 0; j < numPoints; j++) {
313 for (
int k = 0; k < spaceDim; k++) {
316 int l = k + i * spaceDim + j * spaceDim * numFields;
317 if (std::abs(vals(i,j,k) - basisValues[l]) > INTREPID_TOL) {
319 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
322 *outStream <<
" At multi-index { ";
323 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
324 *outStream <<
"} computed value: " << vals(i,j,k)
325 <<
" but reference value: " << basisValues[l] <<
"\n";
332 vals.resize(numFields, numPoints);
333 quadBasis.getValues(vals, quadNodes, OPERATOR_DIV);
334 for (
int i = 0; i < numFields; i++) {
335 for (
int j = 0; j < numPoints; j++) {
336 int l = i + j * numFields;
337 if (std::abs(vals(i,j) - basisDivs[l]) > INTREPID_TOL) {
339 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
342 *outStream <<
" At multi-index { ";
343 *outStream << i <<
" ";*outStream << j <<
" ";
344 *outStream <<
"} computed divergence component: " << vals(i,j)
345 <<
" but reference divergence component: " << basisDivs[l] <<
"\n";
353 catch (
const std::logic_error & err) {
354 *outStream << err.what() <<
"\n\n";
359 std::cout <<
"End Result: TEST FAILED\n";
361 std::cout <<
"End Result: TEST PASSED\n";
364 std::cout.copyfmt(oldFormatState);
Header file for utility class to provide multidimensional containers.
Header file for the Intrepid::HDIV_QUAD_I1_FEM class.