Intrepid: http://trilinos.sandia.gov/packages/docs/r10.0/packages/intrepid/test/Discretization/Integration/test

00001 // @HEADER
00002 // ************************************************************************
00003 //
00004 //                           Intrepid Package
00005 //                 Copyright (2007) Sandia Corporation
00006 //
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 //
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Pavel Bochev (pbboche@sandia.gov) or
00025 //                    Denis Ridzal (dridzal@sandia.gov).
00026 //
00027 // ************************************************************************
00028 // @HEADER
00029 
00030 
00037 #include "Intrepid_CubatureDirectLineGauss.hpp"
00038 #include "Intrepid_CubatureDirectTriDefault.hpp"
00039 #include "Intrepid_CubatureDirectTetDefault.hpp"
00040 #include "Intrepid_CubatureTensor.hpp"
00041 #include "Shards_CellTopology.hpp"
00042 #include "Teuchos_oblackholestream.hpp"
00043 #include "Teuchos_RCP.hpp"
00044 #include "Teuchos_GlobalMPISession.hpp"
00045 
00046 using namespace Intrepid;
00047 
00048 #define INTREPID_TEST_COMMAND( S )                                                                                  \
00049 {                                                                                                                   \
00050   try {                                                                                                             \
00051     S ;                                                                                                             \
00052   }                                                                                                                 \
00053   catch (std::logic_error err) {                                                                                    \
00054       *outStream << "Expected Error ----------------------------------------------------------------\n";            \
00055       *outStream << err.what() << '\n';                                                                             \
00056       *outStream << "-------------------------------------------------------------------------------" << "\n\n";    \
00057   };                                                                                                                \
00058 }
00059 
00060 /*
00061   Computes volume of a given reference cell.
00062 */
00063 double computeRefVolume(shards::CellTopology & cellTopology, int cubDegree) {
00064   Teuchos::RCP< Cubature<double> > myCub;
00065   double vol = 0.0;
00066 
00067   switch (cellTopology.getBaseTopology()->key) {
00068 
00069     case shards::Line<>::key:
00070         myCub = Teuchos::rcp(new CubatureDirectLineGauss<double>(cubDegree));
00071       break;
00072     case shards::Triangle<>::key:
00073         myCub = Teuchos::rcp(new CubatureDirectTriDefault<double>(cubDegree));
00074       break;
00075     case shards::Tetrahedron<>::key:
00076         myCub = Teuchos::rcp(new CubatureDirectTetDefault<double>(cubDegree));
00077       break;
00078     case shards::Quadrilateral<>::key: { 
00079         std::vector< Teuchos::RCP< Cubature<double> > > lineCubs(2);
00080         lineCubs[0] = Teuchos::rcp(new CubatureDirectLineGauss<double>(cubDegree));
00081         lineCubs[1] = Teuchos::rcp(new CubatureDirectLineGauss<double>((cubDegree+7) % INTREPID_CUBATURE_LINE_GAUSS_MAX));
00082         myCub = Teuchos::rcp(new CubatureTensor<double>(lineCubs));
00083         }
00084       break;
00085     case shards::Hexahedron<>::key: { 
00086         std::vector< Teuchos::RCP< Cubature<double> > > lineCubs(2);
00087         std::vector< Teuchos::RCP< Cubature<double> > > miscCubs(2);
00088         lineCubs[0] = Teuchos::rcp(new CubatureDirectLineGauss<double>(cubDegree));
00089         lineCubs[1] = Teuchos::rcp(new CubatureDirectLineGauss<double>((cubDegree+7) % INTREPID_CUBATURE_LINE_GAUSS_MAX));
00090         miscCubs[0] = Teuchos::rcp(new CubatureTensor<double>(lineCubs));
00091         miscCubs[1] = Teuchos::rcp(new CubatureDirectLineGauss<double>((cubDegree+25) % INTREPID_CUBATURE_LINE_GAUSS_MAX));
00092         myCub = Teuchos::rcp(new CubatureTensor<double>(miscCubs));
00093         }
00094       break;
00095     case shards::Wedge<>::key: {
00096         Teuchos::RCP<CubatureDirect<double> > triCub  = Teuchos::rcp(new CubatureDirectTriDefault<double>(cubDegree));
00097         Teuchos::RCP<CubatureDirect<double> > lineCub = Teuchos::rcp(new CubatureDirectLineGauss<double>(cubDegree));
00098         myCub = Teuchos::rcp(new CubatureTensor<double>(triCub,lineCub));
00099         }
00100       break;
00101 
00102     default:
00103       TEST_FOR_EXCEPTION( ( (cellTopology.getBaseTopology()->key != shards::Line<>::key),
00104                             (cellTopology.getBaseTopology()->key != shards::Triangle<>::key),
00105                             (cellTopology.getBaseTopology()->key != shards::Tetrahedron<>::key),
00106                             (cellTopology.getBaseTopology()->key != shards::Quadrilateral<>::key),
00107                             (cellTopology.getBaseTopology()->key != shards::Hexahedron<>::key),
00108                             (cellTopology.getBaseTopology()->key != shards::Wedge<>::key) ),
00109                           std::invalid_argument,
00110                           ">>> ERROR (Unit Test -- Cubature -- Volume): Invalid cell type.");
00111   } // end switch
00112 
00113   int numCubPoints = myCub->getNumPoints();
00114 
00115   FieldContainer<double> cubPoints( numCubPoints, myCub->getDimension() );
00116   FieldContainer<double> cubWeights( numCubPoints );
00117   myCub->getCubature(cubPoints, cubWeights);
00118 
00119   for (int i=0; i<numCubPoints; i++)
00120     vol += cubWeights[i];
00121 
00122   return vol;
00123 }
00124 
00125 
00126 int main(int argc, char *argv[]) {
00127 
00128   Teuchos::GlobalMPISession mpiSession(&argc, &argv);
00129 
00130   // This little trick lets us print to std::cout only if
00131   // a (dummy) command-line argument is provided.
00132   int iprint     = argc - 1;
00133   Teuchos::RCP<std::ostream> outStream;
00134   Teuchos::oblackholestream bhs; // outputs nothing
00135   if (iprint > 0)
00136     outStream = Teuchos::rcp(&std::cout, false);
00137   else
00138     outStream = Teuchos::rcp(&bhs, false);
00139 
00140   // Save the format state of the original std::cout.
00141   Teuchos::oblackholestream oldFormatState;
00142   oldFormatState.copyfmt(std::cout);
00143  
00144   *outStream \
00145   << "===============================================================================\n" \
00146   << "|                                                                             |\n" \
00147   << "|                  Unit Test (CubatureDirect,CubatureTensor)                  |\n" \
00148   << "|                                                                             |\n" \
00149   << "|     1) Computing volumes of reference cells                                 |\n" \
00150   << "|                                                                             |\n" \
00151   << "|  Questions? Contact  Pavel Bochev (pbboche@sandia.gov) or                   |\n" \
00152   << "|                      Denis Ridzal (dridzal@sandia.gov).                     |\n" \
00153   << "|                                                                             |\n" \
00154   << "|  Intrepid's website: http://trilinos.sandia.gov/packages/intrepid           |\n" \
00155   << "|  Trilinos website:   http://trilinos.sandia.gov                             |\n" \
00156   << "|                                                                             |\n" \
00157   << "===============================================================================\n"\
00158   << "| TEST 1: construction and basic functionality                                |\n"\
00159   << "===============================================================================\n";
00160 
00161   int errorFlag  = 0;
00162 
00163   int beginThrowNumber = TestForException_getThrowNumber();
00164   int endThrowNumber = beginThrowNumber + 7;  
00165 
00166   try {
00167     /* Line cubature. */
00168     INTREPID_TEST_COMMAND( CubatureDirectLineGauss<double> lineCub(-1) );
00169     INTREPID_TEST_COMMAND( CubatureDirectLineGauss<double> lineCub(INTREPID_CUBATURE_LINE_GAUSS_MAX+1) );
00170     INTREPID_TEST_COMMAND( CubatureDirectLineGauss<double> lineCub;
00171                            std::string testName    = "INTREPID_CUBATURE_LINE_GAUSS";
00172                            std::string lineCubName = lineCub.getName();
00173                            *outStream << "\nComparing strings: " << testName << " and " << lineCubName << "\n\n";
00174                            TEST_FOR_EXCEPTION( (testName != lineCubName), std::logic_error, "Name mismatch!" ) );
00175     INTREPID_TEST_COMMAND( CubatureDirectLineGauss<double> lineCub;
00176                            std::vector<int> accuracy;
00177                            lineCub.getAccuracy(accuracy);
00178                            TEST_FOR_EXCEPTION( (accuracy[0] != 0), std::logic_error, "Check member getAccuracy!" ) );
00179     INTREPID_TEST_COMMAND( CubatureDirectLineGauss<double> lineCub(55);
00180                            TEST_FOR_EXCEPTION( (lineCub.getNumPoints() != 28), std::logic_error, "Check member getNumPoints!" ) );
00181     INTREPID_TEST_COMMAND( CubatureDirectLineGauss<double> lineCub;
00182                            TEST_FOR_EXCEPTION( (lineCub.getDimension() != 1),
00183                                                std::logic_error,
00184                                                "Check member dimension!" ) );
00185     /* Triangle cubature. */
00186     INTREPID_TEST_COMMAND( CubatureDirectTriDefault<double> triCub(-1) );
00187     INTREPID_TEST_COMMAND( CubatureDirectTriDefault<double> triCub(INTREPID_CUBATURE_TRI_DEFAULT_MAX+1) );
00188     INTREPID_TEST_COMMAND( CubatureDirectTriDefault<double> triCub;
00189                            std::string testName    = "INTREPID_CUBATURE_TRI_DEFAULT";
00190                            std::string triCubName = triCub.getName();
00191                            *outStream << "\nComparing strings: " << testName << " and " << triCubName << "\n\n";
00192                            TEST_FOR_EXCEPTION( (testName != triCubName), std::logic_error, "Name mismatch!" ) );
00193     INTREPID_TEST_COMMAND( CubatureDirectTriDefault<double> triCub;
00194                            std::vector<int> accuracy;
00195                            triCub.getAccuracy(accuracy);
00196                            TEST_FOR_EXCEPTION( (accuracy[0] != 0), std::logic_error, "Check member getAccuracy!" ) );
00197     INTREPID_TEST_COMMAND( CubatureDirectTriDefault<double> triCub(17);
00198                            TEST_FOR_EXCEPTION( (triCub.getNumPoints() != 61), std::logic_error, "Check member getNumPoints!" ) );
00199     INTREPID_TEST_COMMAND( CubatureDirectTriDefault<double> triCub;
00200                            TEST_FOR_EXCEPTION( (triCub.getDimension() != 2),
00201                                                std::logic_error,
00202                                                "Check member dimension!" ) );
00203     /* Tetrahedron cubature. */
00204     INTREPID_TEST_COMMAND( CubatureDirectTetDefault<double> tetCub(-1) );
00205     INTREPID_TEST_COMMAND( CubatureDirectTetDefault<double> tetCub(INTREPID_CUBATURE_TET_DEFAULT_MAX+1) );
00206     INTREPID_TEST_COMMAND( CubatureDirectTetDefault<double> tetCub;
00207                            std::string testName    = "INTREPID_CUBATURE_TET_DEFAULT";
00208                            std::string tetCubName = tetCub.getName();
00209                            *outStream << "\nComparing strings: " << testName << " and " << tetCubName << "\n\n";
00210         std::vector< Teuchos::RCP< Cubature<double> > > lineCubs(2);
00211                            TEST_FOR_EXCEPTION( (testName != tetCubName), std::logic_error, "Name mismatch!" ) );
00212     INTREPID_TEST_COMMAND( CubatureDirectTetDefault<double> tetCub;
00213                            std::vector<int> accuracy;
00214                            tetCub.getAccuracy(accuracy);
00215                            TEST_FOR_EXCEPTION( (accuracy[0] != 0), std::logic_error, "Check member getAccuracy!" ) );
00216     INTREPID_TEST_COMMAND( CubatureDirectTetDefault<double> tetCub(17);
00217                            TEST_FOR_EXCEPTION( (tetCub.getNumPoints() != 495), std::logic_error, "Check member getNumPoints!" ) );
00218     INTREPID_TEST_COMMAND( CubatureDirectTetDefault<double> tetCub;
00219                            TEST_FOR_EXCEPTION( (tetCub.getDimension() != 3),
00220                                                std::logic_error,
00221                                                "Check member getCellTopology!" ) );
00222     /* Tensor cubature. */
00223     INTREPID_TEST_COMMAND( std::vector< Teuchos::RCP< Cubature<double> > > lineCubs(0);
00224                            CubatureTensor<double> quadCub(lineCubs) );
00225     INTREPID_TEST_COMMAND( std::vector< Teuchos::RCP< Cubature<double> > > miscCubs(3);
00226                            std::vector< Teuchos::RCP< Cubature<double> > > lineCubs(2);
00227                            lineCubs[0] = Teuchos::rcp(new CubatureDirectLineGauss<double>(3));
00228                            lineCubs[1] = Teuchos::rcp(new CubatureDirectLineGauss<double>(16));
00229                            miscCubs[0] = Teuchos::rcp(new CubatureTensor<double>(lineCubs));
00230                            miscCubs[1] = Teuchos::rcp(new CubatureDirectTriDefault<double>);
00231                            miscCubs[2] = Teuchos::rcp(new CubatureDirectTetDefault<double>(19));
00232                            CubatureTensor<double> tensorCub(miscCubs);
00233                            std::vector<int> a(4); a[0]=3; a[1]=16; a[2]=0; a[3]=19;
00234                            std::vector<int> atest(4);
00235                            tensorCub.getAccuracy(atest);
00236                            TEST_FOR_EXCEPTION( (a != atest), std::logic_error, "Check member getAccuracy!" ) );
00237     INTREPID_TEST_COMMAND( std::vector< Teuchos::RCP< Cubature<double> > > lineCubs(2);
00238                            lineCubs[0] = Teuchos::rcp(new CubatureDirectLineGauss<double>(15));
00239                            lineCubs[1] = Teuchos::rcp(new CubatureDirectLineGauss<double>(11));
00240                            CubatureTensor<double> tensorCub(lineCubs);
00241                            TEST_FOR_EXCEPTION( (tensorCub.getNumPoints() != 48), std::logic_error, "Check member getNumPoints!" ) );
00242     INTREPID_TEST_COMMAND( std::vector< Teuchos::RCP< Cubature<double> > > miscCubs(3);
00243                            miscCubs[0] = Teuchos::rcp(new CubatureDirectLineGauss<double>);
00244                            miscCubs[1] = Teuchos::rcp(new CubatureDirectTriDefault<double>);
00245                            miscCubs[2] = Teuchos::rcp(new CubatureDirectTetDefault<double>);
00246                            CubatureTensor<double> tensorCub(miscCubs);
00247                            TEST_FOR_EXCEPTION( (tensorCub.getDimension() != 6), std::logic_error, "Check member dimension!" ) );
00248     INTREPID_TEST_COMMAND( std::vector< Teuchos::RCP< Cubature<double> > > miscCubs(3);
00249                            miscCubs[0] = Teuchos::rcp(new CubatureDirectLineGauss<double>(3));
00250                            miscCubs[1] = Teuchos::rcp(new CubatureDirectLineGauss<double>(7));
00251                            miscCubs[2] = Teuchos::rcp(new CubatureDirectLineGauss<double>(5));
00252                            CubatureTensor<double> tensorCub(miscCubs);
00253                            FieldContainer<double> points(tensorCub.getNumPoints(), tensorCub.getDimension());
00254                            FieldContainer<double> weights(tensorCub.getNumPoints());
00255                            tensorCub.getCubature(points, weights)
00256                          )
00257     INTREPID_TEST_COMMAND( Teuchos::RCP< CubatureDirect<double> > lineCub = Teuchos::rcp(new CubatureDirectLineGauss<double>(15));
00258                            Teuchos::RCP< CubatureDirect<double> > triCub = Teuchos::rcp(new CubatureDirectTriDefault<double>(12));
00259                            CubatureTensor<double> tensorCub(lineCub, triCub);
00260                            std::vector<int> a(2); a[0] = 15; a[1] = 12;
00261                            std::vector<int> atest(2);
00262                            tensorCub.getAccuracy(atest);
00263                            TEST_FOR_EXCEPTION( (tensorCub.getDimension() != 3) || (a != atest),
00264                                                std::logic_error,
00265                                                "Check constructormembers dimension and getAccuracy!" ) );
00266     INTREPID_TEST_COMMAND( Teuchos::RCP< CubatureDirect<double> > lineCub = Teuchos::rcp(new CubatureDirectLineGauss<double>(15));
00267                            Teuchos::RCP< CubatureDirect<double> > triCub = Teuchos::rcp(new CubatureDirectTriDefault<double>(12));
00268                            CubatureTensor<double> tensorCub(triCub, lineCub, triCub);
00269                            std::vector<int> a(3); a[0] = 12; a[1] = 15; a[2] = 12;
00270                            std::vector<int> atest(3);
00271                            tensorCub.getAccuracy(atest);
00272                            TEST_FOR_EXCEPTION( (tensorCub.getDimension() != 5) || (a != atest),
00273                                                std::logic_error,
00274                                                "Check constructor and members dimension and getAccuracy!" ) );
00275     INTREPID_TEST_COMMAND( Teuchos::RCP< CubatureDirect<double> > triCub = Teuchos::rcp(new CubatureDirectTriDefault<double>(12));
00276                            CubatureTensor<double> tensorCub(triCub, 5);
00277                            std::vector<int> a(5); a[0] = 12; a[1] = 12; a[2] = 12; a[3] = 12; a[4] = 12;
00278                            std::vector<int> atest(5);
00279                            tensorCub.getAccuracy(atest);
00280                            TEST_FOR_EXCEPTION( (tensorCub.getDimension() != 10) || (a != atest),
00281                                                std::logic_error,
00282                                                "Check constructor and members dimension and getAccuracy!" ) );
00283     if (TestForException_getThrowNumber() != endThrowNumber) {
00284       errorFlag = -1000;
00285     }
00286   }
00287   catch (std::logic_error err) {
00288     *outStream << err.what() << "\n";
00289     errorFlag = -1000;
00290   };
00291  
00292   *outStream \
00293   << "===============================================================================\n"\
00294   << "| TEST 2: volume computations                                                 |\n"\
00295   << "===============================================================================\n";
00296 
00297   double testVol = 0.0;
00298   double tol     = 100.0 * INTREPID_TOL;
00299 
00300   // list of analytic volume values, listed in the enumerated reference cell order up to CELL_HEXAPRISM
00301   double volumeList[] = {0.0, 2.0, 1.0/2.0, 4.0, 1.0/6.0, 8.0, 1.0, 32.0};
00302 
00303   *outStream << "\nReference cell volumes:\n\n";
00304 
00305   try {
00306     shards::CellTopology line(shards::getCellTopologyData< shards::Line<> >());
00307     for (int deg=0; deg<=INTREPID_CUBATURE_LINE_GAUSS_MAX; deg++) {
00308       testVol = computeRefVolume(line, deg);
00309       *outStream << std::setw(30) << "Line volume --> " << std::setw(10) << std::scientific << testVol <<
00310                     std::setw(10) << "diff = " << std::setw(10) << std::scientific << std::abs(testVol - volumeList[1]) << "\n";
00311       if (std::abs(testVol - volumeList[1]) > tol) {
00312         errorFlag = 1;
00313         *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
00314       }
00315     }
00316     *outStream << "\n\n";
00317     shards::CellTopology tri(shards::getCellTopologyData< shards::Triangle<> >());
00318     for (int deg=0; deg<=INTREPID_CUBATURE_TRI_DEFAULT_MAX; deg++) {
00319       testVol = computeRefVolume(tri, deg);
00320       *outStream << std::setw(30) << "Triangle volume --> " << std::setw(10) << std::scientific << testVol <<
00321                     std::setw(10) << "diff = " << std::setw(10) << std::scientific << std::abs(testVol - volumeList[2]) << "\n";
00322       if (std::abs(testVol - volumeList[2]) > tol) {
00323         errorFlag = 1;
00324         *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
00325       }
00326     }
00327     *outStream << "\n\n";
00328     shards::CellTopology quad(shards::getCellTopologyData< shards::Quadrilateral<> >());
00329     for (int deg=0; deg<=INTREPID_CUBATURE_LINE_GAUSS_MAX; deg++) {
00330       testVol = computeRefVolume(quad, deg);
00331       *outStream << std::setw(30) << "Quadrilateral volume --> " << std::setw(10) << std::scientific << testVol <<
00332                     std::setw(10) << "diff = " << std::setw(10) << std::scientific << std::abs(testVol - volumeList[3]) << "\n";
00333       if (std::abs(testVol - volumeList[3]) > tol) {
00334         errorFlag = 1;
00335         *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
00336       }
00337     }
00338     *outStream << "\n\n";
00339     shards::CellTopology tet(shards::getCellTopologyData< shards::Tetrahedron<> >());
00340     for (int deg=0; deg<=INTREPID_CUBATURE_TET_DEFAULT_MAX; deg++) {
00341       testVol = computeRefVolume(tet, deg);
00342       *outStream << std::setw(30) << "Tetrahedron volume --> " << std::setw(10) << std::scientific << testVol <<
00343                     std::setw(10) << "diff = " << std::setw(10) << std::scientific << std::abs(testVol - volumeList[4]) << "\n";
00344       if (std::abs(testVol - volumeList[4]) > tol) {
00345         errorFlag = 1;
00346         *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
00347       }
00348     }
00349     *outStream << "\n\n";
00350     shards::CellTopology hex(shards::getCellTopologyData< shards::Hexahedron<> >());
00351     for (int deg=0; deg<=INTREPID_CUBATURE_LINE_GAUSS_MAX; deg++) {
00352       testVol = computeRefVolume(hex, deg);
00353       *outStream << std::setw(30) << "Hexahedron volume --> " << std::setw(10) << std::scientific << testVol <<
00354                     std::setw(10) << "diff = " << std::setw(10) << std::scientific << std::abs(testVol - volumeList[5]) << "\n";
00355       if (std::abs(testVol - volumeList[5]) > tol) {
00356         errorFlag = 1;
00357         *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
00358       }
00359     }
00360     *outStream << "\n\n";
00361     shards::CellTopology wedge(shards::getCellTopologyData< shards::Wedge<> >());
00362     for (int deg=0; deg<=std::min(INTREPID_CUBATURE_LINE_GAUSS_MAX,INTREPID_CUBATURE_TRI_DEFAULT_MAX); deg++) {
00363       testVol = computeRefVolume(wedge, deg);
00364       *outStream << std::setw(30) << "Wedge volume --> " << std::setw(10) << std::scientific << testVol <<
00365                     std::setw(10) << "diff = " << std::setw(10) << std::scientific << std::abs(testVol - volumeList[6]) << "\n";
00366       if (std::abs(testVol - volumeList[6]) > tol) {
00367         errorFlag = 1;
00368         *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
00369       }
00370     }
00371     *outStream << "\n\n";
00372     for (int deg=0; deg<=20; deg++) {
00373       Teuchos::RCP<CubatureDirectLineGauss<double> > lineCub = Teuchos::rcp(new CubatureDirectLineGauss<double>(deg));
00374       CubatureTensor<double> hypercubeCub(lineCub, 5);
00375       int numCubPoints = hypercubeCub.getNumPoints();
00376       FieldContainer<double> cubPoints( numCubPoints, hypercubeCub.getDimension() );
00377       FieldContainer<double> cubWeights( numCubPoints );
00378       hypercubeCub.getCubature(cubPoints, cubWeights);
00379       testVol = 0;
00380       for (int i=0; i<numCubPoints; i++)
00381         testVol += cubWeights[i];
00382       *outStream << std::setw(30) << "5-D Hypercube volume --> " << std::setw(10) << std::scientific << testVol <<
00383                     std::setw(10) << "diff = " << std::setw(10) << std::scientific << std::abs(testVol - volumeList[7]) << "\n";
00384       if (std::abs(testVol - volumeList[7])/std::abs(testVol) > tol) {
00385         errorFlag = 1;
00386         *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
00387       }
00388     }
00389   }
00390   catch (std::logic_error err) {
00391     *outStream << err.what() << "\n";
00392     errorFlag = -1;
00393   };
00394 
00395 
00396   if (errorFlag != 0)
00397     std::cout << "End Result: TEST FAILED\n";
00398   else
00399     std::cout << "End Result: TEST PASSED\n";
00400 
00401   // reset format state of std::cout
00402   std::cout.copyfmt(oldFormatState);
00403 
00404   return errorFlag;
00405 }
http://trilinos.sandia.gov/packages/docs/r10.0/packages/intrepid/test/Discretization/Integration/test_01.cpp
