50 #include "Teuchos_oblackholestream.hpp"
51 #include "Teuchos_RCP.hpp"
52 #include "Teuchos_GlobalMPISession.hpp"
55 using 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"; \
71 int 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_HCURL_TRI_I1_FEM) |\n" \
94 <<
"| 1) Conversion of Dof tags into Dof ordinals and back |\n" \
95 <<
"| 2) Basis values for VALUE and CURL 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";
114 int throwCounter = 0;
118 triNodes(0,0) = 0.0; triNodes(0,1) = 0.0;
119 triNodes(1,0) = 1.0; triNodes(1,1) = 0.0;
120 triNodes(2,0) = 0.0; triNodes(2,1) = 1.0;
122 triNodes(3,0) = 0.5; triNodes(3,1) = 0.0;
123 triNodes(4,0) = 0.5; triNodes(4,1) = 0.5;
124 triNodes(5,0) = 0.0; triNodes(5,1) = 0.5;
126 triNodes(6,0) = 0.25; triNodes(6,1) = 0.25;
136 INTREPID_TEST_COMMAND( triBasis.
getValues(vals, triNodes, OPERATOR_GRAD), throwCounter, nException );
141 INTREPID_TEST_COMMAND( triBasis.
getValues(vals, triNodes, OPERATOR_DIV), throwCounter, nException );
146 INTREPID_TEST_COMMAND( triBasis.
getDofOrdinal(3,0,0), throwCounter, nException );
148 INTREPID_TEST_COMMAND( triBasis.
getDofOrdinal(1,1,1), throwCounter, nException );
150 INTREPID_TEST_COMMAND( triBasis.
getDofOrdinal(0,4,1), throwCounter, nException );
152 INTREPID_TEST_COMMAND( triBasis.
getDofTag(12), throwCounter, nException );
154 INTREPID_TEST_COMMAND( triBasis.
getDofTag(-1), throwCounter, nException );
156 #ifdef HAVE_INTREPID_DEBUG
160 INTREPID_TEST_COMMAND( triBasis.
getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException );
164 INTREPID_TEST_COMMAND( triBasis.
getValues(vals, badPoints2, OPERATOR_VALUE), throwCounter, nException );
168 INTREPID_TEST_COMMAND( triBasis.
getValues(badVals1, triNodes, OPERATOR_VALUE), throwCounter, nException );
172 INTREPID_TEST_COMMAND( triBasis.
getValues(badCurls1, triNodes, OPERATOR_CURL), throwCounter, nException );
176 INTREPID_TEST_COMMAND( triBasis.
getValues(badVals2, triNodes, OPERATOR_VALUE), throwCounter, nException ) ;
180 INTREPID_TEST_COMMAND( triBasis.
getValues(badVals3, triNodes, OPERATOR_VALUE), throwCounter, nException ) ;
184 INTREPID_TEST_COMMAND( triBasis.
getValues(badVals4, triNodes, OPERATOR_VALUE), throwCounter, nException ) ;
189 triNodes.dimension(0),
191 INTREPID_TEST_COMMAND( triBasis.
getValues(vals, triNodes, OPERATOR_D2), throwCounter, nException );
195 catch (
const std::logic_error & err) {
196 *outStream <<
"UNEXPECTED ERROR !!! ----------------------------------------------------------\n";
197 *outStream << err.what() <<
'\n';
198 *outStream <<
"-------------------------------------------------------------------------------" <<
"\n\n";
204 if (throwCounter != nException) {
206 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
211 <<
"===============================================================================\n"\
212 <<
"| TEST 2: correctness of tag to enum and enum to tag lookups |\n"\
213 <<
"===============================================================================\n";
216 std::vector<std::vector<int> > allTags = triBasis.
getAllDofTags();
219 for (
unsigned i = 0; i < allTags.size(); i++) {
220 int bfOrd = triBasis.
getDofOrdinal(allTags[i][0], allTags[i][1], allTags[i][2]);
222 std::vector<int> myTag = triBasis.
getDofTag(bfOrd);
223 if( !( (myTag[0] == allTags[i][0]) &&
224 (myTag[1] == allTags[i][1]) &&
225 (myTag[2] == allTags[i][2]) &&
226 (myTag[3] == allTags[i][3]) ) ) {
228 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
229 *outStream <<
" getDofOrdinal( {"
230 << allTags[i][0] <<
", "
231 << allTags[i][1] <<
", "
232 << allTags[i][2] <<
", "
233 << allTags[i][3] <<
"}) = " << bfOrd <<
" but \n";
234 *outStream <<
" getDofTag(" << bfOrd <<
") = { "
238 << myTag[3] <<
"}\n";
244 std::vector<int> myTag = triBasis.
getDofTag(bfOrd);
245 int myBfOrd = triBasis.
getDofOrdinal(myTag[0], myTag[1], myTag[2]);
246 if( bfOrd != myBfOrd) {
248 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
249 *outStream <<
" getDofTag(" << bfOrd <<
") = { "
253 << myTag[3] <<
"} but getDofOrdinal({"
257 << myTag[3] <<
"} ) = " << myBfOrd <<
"\n";
261 catch (
const std::logic_error & err){
262 *outStream << err.what() <<
"\n\n";
268 <<
"===============================================================================\n"\
269 <<
"| TEST 3: correctness of basis function values |\n"\
270 <<
"===============================================================================\n";
272 outStream -> precision(20);
275 double basisValues[] = {
276 1.000, 0, 0, 0, 0, -1.000, 1.000, 1.000, 0, 1.000, 0, 0, 0, 0, \
277 -1.000, 0, -1.000, -1.000, 1.000, 0.5000, 0, 0.5000, 0, -0.5000, \
278 0.5000, 0.5000, -0.5000, 0.5000, -0.5000, -0.5000, 0.5000, 0, \
279 -0.5000, 0, -0.5000, -1.000, 0.7500, 0.2500, -0.2500, 0.2500, \
283 double basisCurls[] = {
297 int numPoints = triNodes.dimension(0);
304 vals.
resize(numFields, numPoints, spaceDim);
305 triBasis.
getValues(vals, triNodes, OPERATOR_VALUE);
306 for (
int i = 0; i < numFields; i++) {
307 for (
int j = 0; j < numPoints; j++) {
308 for (
int k = 0; k < spaceDim; k++) {
311 int l = k + i * spaceDim + j * spaceDim * numFields;
312 if (std::abs(vals(i,j,k) - basisValues[l]) > INTREPID_TOL) {
314 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
317 *outStream <<
" At multi-index { ";
318 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
319 *outStream <<
"} computed value: " << vals(i,j,k)
320 <<
" but reference value: " << basisValues[l] <<
"\n";
327 vals.
resize(numFields, numPoints);
328 triBasis.
getValues(vals, triNodes, OPERATOR_CURL);
329 for (
int i = 0; i < numFields; i++) {
330 for (
int j = 0; j < numPoints; j++) {
331 int l = i + j * numFields;
332 if (std::abs(vals(i,j) - basisCurls[l]) > INTREPID_TOL) {
334 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
337 *outStream <<
" At multi-index { ";
338 *outStream << i <<
" ";*outStream << j <<
" ";
339 *outStream <<
"} computed curl component: " << vals(i,j)
340 <<
" but reference curl component: " << basisCurls[l] <<
"\n";
347 catch (
const std::logic_error & err) {
348 *outStream << err.what() <<
"\n\n";
354 <<
"===============================================================================\n"\
355 <<
"| TEST 4: correctness of DoF locations |\n"\
356 <<
"===============================================================================\n";
359 Teuchos::RCP<Basis<double, FieldContainer<double> > > basis =
361 Teuchos::RCP<DofCoordsInterface<FieldContainer<double> > > coord_iface =
369 #ifdef HAVE_INTREPID_DEBUG
371 INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
373 INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
375 INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
378 INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException ); nException--;
380 if (throwCounter != nException) {
382 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
388 tangents(0,0) = 1.0; tangents(0,1) = 0.0;
389 tangents(1,0) = -1.0; tangents(1,1) = 1.0;
390 tangents(2,0) = 0.0; tangents(2,1) = -1.0;
392 basis->getValues(bvals, cvals, OPERATOR_VALUE);
394 for (
int i=0; i<bvals.dimension(0); i++) {
395 for (
int j=0; j<bvals.dimension(1); j++) {
397 double tangent = 0.0;
398 for(
int d=0;d<spaceDim;d++)
399 tangent += bvals(i,j,d)*tangents(j,d);
401 if ((i != j) && (std::abs(tangent - 0.0) > INTREPID_TOL)) {
403 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), tangent, 0.0);
404 *outStream << buffer;
406 else if ((i == j) && (std::abs(tangent - 1.0) > INTREPID_TOL)) {
408 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), tangent, 1.0);
409 *outStream << buffer;
415 catch (
const std::logic_error & err){
416 *outStream << err.what() <<
"\n\n";
422 std::cout <<
"End Result: TEST FAILED\n";
424 std::cout <<
"End Result: TEST PASSED\n";
427 std::cout.copyfmt(oldFormatState);
Implementation of the default H(curl)-compatible FEM basis of degree 1 on Triangle cell...
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 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
Evaluation of a FEM basis on a reference Triangle 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...
Header file for the Intrepid::HCURL_TRI_I1_FEM class.