|
Intrepid
|
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), 00025 // Denis Ridzal (dridzal@sandia.gov), 00026 // Kara Peterson (kjpeter@sandia.gov). 00027 // 00028 // ************************************************************************ 00029 // @HEADER 00030 00036 #include "Intrepid_FieldContainer.hpp" 00037 #include "Intrepid_HGRAD_LINE_Cn_FEM.hpp" 00038 #include "Intrepid_DefaultCubatureFactory.hpp" 00039 #include "Intrepid_ArrayTools.hpp" 00040 #include "Intrepid_PointTools.hpp" 00041 #include "Intrepid_FunctionSpaceTools.hpp" 00042 #include "Teuchos_oblackholestream.hpp" 00043 #include "Teuchos_RCP.hpp" 00044 #include "Teuchos_GlobalMPISession.hpp" 00045 00046 using namespace std; 00047 using namespace Intrepid; 00048 00049 #define INTREPID_TEST_COMMAND( S , throwCounter, nException ) \ 00050 { \ 00051 ++nException; \ 00052 try { \ 00053 S ; \ 00054 } \ 00055 catch (std::logic_error err) { \ 00056 ++throwCounter; \ 00057 *outStream << "Expected Error " << nException << " -------------------------------------------------------------\n"; \ 00058 *outStream << err.what() << '\n'; \ 00059 *outStream << "-------------------------------------------------------------------------------" << "\n\n"; \ 00060 }; \ 00061 } 00062 00063 00064 int main(int argc, char *argv[]) { 00065 00066 Teuchos::GlobalMPISession mpiSession(&argc, &argv); 00067 00068 // This little trick lets us print to std::cout only if 00069 // a (dummy) command-line argument is provided. 00070 int iprint = argc - 1; 00071 Teuchos::RCP<std::ostream> outStream; 00072 Teuchos::oblackholestream bhs; // outputs nothing 00073 if (iprint > 0) 00074 outStream = Teuchos::rcp(&std::cout, false); 00075 else 00076 outStream = Teuchos::rcp(&bhs, false); 00077 00078 // Save the format state of the original std::cout. 00079 Teuchos::oblackholestream oldFormatState; 00080 oldFormatState.copyfmt(std::cout); 00081 00082 *outStream \ 00083 << "===============================================================================\n" \ 00084 << "| |\n" \ 00085 << "| Unit Test (Basis_HGRAD_LINE_Cn_FEM) |\n" \ 00086 << "| |\n" \ 00087 << "| 1) Conversion of Dof tags into Dof ordinals and back |\n" \ 00088 << "| 2) Basis values for VALUE, GRAD, CURL, and Dk operators |\n" \ 00089 << "| |\n" \ 00090 << "| Questions? Contact Pavel Bochev (pbboche@sandia.gov), |\n" \ 00091 << "| Denis Ridzal (dridzal@sandia.gov), |\n" \ 00092 << "| Kara Peterson (kjpeter@sandia.gov). |\n" \ 00093 << "| |\n" \ 00094 << "| Intrepid's website: http://trilinos.sandia.gov/packages/intrepid |\n" \ 00095 << "| Trilinos website: http://trilinos.sandia.gov |\n" \ 00096 << "| |\n" \ 00097 << "===============================================================================\n"\ 00098 << "| TEST 1: Basis creation, exception testing |\n"\ 00099 << "===============================================================================\n"; 00100 00101 00102 // Define basis and error flag 00103 shards::CellTopology line(shards::getCellTopologyData< shards::Line<> >()); // create cell topology 00104 const int deg = 5; 00105 00106 // get the points for the basis 00107 FieldContainer<double> pts(PointTools::getLatticeSize(line,deg),1); 00108 PointTools::getLattice<double,FieldContainer<double> >(pts,line,deg); 00109 00110 Basis_HGRAD_LINE_Cn_FEM<double, FieldContainer<double> > lineBasis(deg, pts); 00111 int errorFlag = 0; 00112 00113 // Initialize throw counter for exception testing 00114 int nException = 0; 00115 int throwCounter = 0; 00116 00117 // Define array containing vertices of the reference Line and a few other points 00118 int numIntervals = 100; 00119 FieldContainer<double> lineNodes(numIntervals+1, 1); 00120 for (int i=0; i<numIntervals+1; i++) { 00121 lineNodes(i,0) = -1.0+(2.0*(double)i)/(double)numIntervals; 00122 } 00123 00124 // Generic array for the output values; needs to be properly resized depending on the operator type 00125 FieldContainer<double> vals; 00126 00127 try{ 00128 // exception #1: DIV cannot be applied to scalar functions 00129 // resize vals to rank-2 container with dimensions (num. points, num. basis functions) 00130 vals.resize(lineBasis.getCardinality(), lineNodes.dimension(0) ); 00131 //INTREPID_TEST_COMMAND( lineBasis.getValues(vals, lineNodes, OPERATOR_DIV), throwCounter, nException ); 00132 00133 // Exceptions 1-5: all bf tags/bf Ids below are wrong and should cause getDofOrdinal() and 00134 // getDofTag() to access invalid array elements thereby causing bounds check exception 00135 // exception #1 00136 INTREPID_TEST_COMMAND( lineBasis.getDofOrdinal(2,0,0), throwCounter, nException ); 00137 // exception #2 00138 INTREPID_TEST_COMMAND( lineBasis.getDofOrdinal(1,1,1), throwCounter, nException ); 00139 // exception #3 00140 INTREPID_TEST_COMMAND( lineBasis.getDofOrdinal(1,0,7), throwCounter, nException ); 00141 // exception #4 00142 INTREPID_TEST_COMMAND( lineBasis.getDofTag(6), throwCounter, nException ); 00143 // exception #5 00144 INTREPID_TEST_COMMAND( lineBasis.getDofTag(-1), throwCounter, nException ); 00145 // not an exception 00146 INTREPID_TEST_COMMAND( lineBasis.getDofTag(5), throwCounter, nException ); --nException; 00147 00148 #ifdef HAVE_INTREPID_DEBUG 00149 // Exceptions 6-14 test exception handling with incorrectly dimensioned input/output arrays 00150 // exception #6: input points array must be of rank-2 00151 FieldContainer<double> badPoints1(4, 5, 3); 00152 INTREPID_TEST_COMMAND( lineBasis.getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException ); 00153 00154 // exception #7: dimension 1 in the input point array must equal space dimension of the cell 00155 FieldContainer<double> badPoints2(4, 3); 00156 INTREPID_TEST_COMMAND( lineBasis.getValues(vals, badPoints2, OPERATOR_VALUE), throwCounter, nException ); 00157 00158 // exception #8: output values must be of rank-2 for OPERATOR_VALUE 00159 FieldContainer<double> badVals1(4, 3, 1); 00160 INTREPID_TEST_COMMAND( lineBasis.getValues(badVals1, lineNodes, OPERATOR_VALUE), throwCounter, nException ); 00161 00162 // exception #9: output values must be of rank-3 for OPERATOR_GRAD 00163 FieldContainer<double> badVals2(4, 3); 00164 INTREPID_TEST_COMMAND( lineBasis.getValues(badVals2, lineNodes, OPERATOR_GRAD), throwCounter, nException ); 00165 00166 // exception #10: output values must be of rank-2 for OPERATOR_D1 00167 INTREPID_TEST_COMMAND( lineBasis.getValues(badVals2, lineNodes, OPERATOR_D1), throwCounter, nException ); 00168 00169 // exception #11: incorrect 0th dimension of output array (must equal number of basis functions) 00170 FieldContainer<double> badVals3(lineBasis.getCardinality() + 1, lineNodes.dimension(0)); 00171 INTREPID_TEST_COMMAND( lineBasis.getValues(badVals3, lineNodes, OPERATOR_VALUE), throwCounter, nException ); 00172 00173 // exception #12: incorrect 1st dimension of output array (must equal number of points) 00174 FieldContainer<double> badVals4(lineBasis.getCardinality(), lineNodes.dimension(0) + 1); 00175 INTREPID_TEST_COMMAND( lineBasis.getValues(badVals4, lineNodes, OPERATOR_VALUE), throwCounter, nException ); 00176 00177 // exception #13: incorrect 2nd dimension of output array (must equal spatial dimension) 00178 FieldContainer<double> badVals5(lineBasis.getCardinality(), lineNodes.dimension(0), 2); 00179 INTREPID_TEST_COMMAND( lineBasis.getValues(badVals5, lineNodes, OPERATOR_GRAD), throwCounter, nException ); 00180 00181 00182 // not an exception 00183 FieldContainer<double> goodVals2(lineBasis.getCardinality(), lineNodes.dimension(0)); 00184 INTREPID_TEST_COMMAND( lineBasis.getValues(goodVals2, lineNodes, OPERATOR_VALUE), throwCounter, nException ); --nException; 00185 #endif 00186 00187 } 00188 catch (std::logic_error err) { 00189 *outStream << "UNEXPECTED ERROR !!! ----------------------------------------------------------\n"; 00190 *outStream << err.what() << '\n'; 00191 *outStream << "-------------------------------------------------------------------------------" << "\n\n"; 00192 errorFlag = -1000; 00193 }; 00194 00195 // Check if number of thrown exceptions matches the one we expect 00196 if (throwCounter != nException) { 00197 errorFlag++; 00198 *outStream << std::setw(70) << "FAILURE! Incorrect number of exceptions." << "\n"; 00199 } 00200 00201 *outStream \ 00202 << "\n" 00203 << "===============================================================================\n" \ 00204 << "| TEST 2: correctness of basis function values |\n" \ 00205 << "===============================================================================\n"; 00206 outStream -> precision(20); 00207 00208 00209 try { 00210 // Check Kronecker property for Lagrange polynomials. 00211 int maxorder = 4; 00212 00213 for (int ordi=1; ordi <= maxorder; ordi++) { 00214 FieldContainer<double> pts(PointTools::getLatticeSize(line,ordi),1); 00215 PointTools::getLattice<double,FieldContainer<double> >(pts,line,ordi); 00216 Basis_HGRAD_LINE_Cn_FEM<double, FieldContainer<double> > lineBasis(ordi, pts); 00217 FieldContainer<double> vals(lineBasis.getCardinality(),pts.dimension(0)); 00218 lineBasis.getValues(vals,pts,OPERATOR_VALUE); 00219 for (int i=0;i<lineBasis.getCardinality();i++) { 00220 for (int j=0;j<pts.dimension(0);j++) { 00221 if ( i==j && std::abs( vals(i,j) - 1.0 ) > INTREPID_TOL ) { 00222 errorFlag++; 00223 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n"; 00224 *outStream << " Basis function " << i << " does not have unit value at its node\n"; 00225 } 00226 if ( i!=j && std::abs( vals(i,j) ) > INTREPID_TOL ) { 00227 errorFlag++; 00228 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n"; 00229 *outStream << " Basis function " << i << " does not vanish at node " << j << "\n"; 00230 } 00231 } 00232 } 00233 } 00234 } 00235 // Catch unexpected errors 00236 catch (std::logic_error err) { 00237 *outStream << err.what() << "\n\n"; 00238 errorFlag = -1000; 00239 }; 00240 00241 if (errorFlag != 0) 00242 std::cout << "End Result: TEST FAILED\n"; 00243 else 00244 std::cout << "End Result: TEST PASSED\n"; 00245 00246 // reset format state of std::cout 00247 std::cout.copyfmt(oldFormatState); 00248 00249 return errorFlag; 00250 }
1.7.4