|
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) or 00025 // Denis Ridzal (dridzal@sandia.gov). 00026 // 00027 // ************************************************************************ 00028 // @HEADER 00029 00035 #include "Intrepid_ArrayTools.hpp" 00036 #include "Intrepid_FieldContainer.hpp" 00037 #include "Intrepid_RealSpaceTools.hpp" 00038 #include "Teuchos_oblackholestream.hpp" 00039 #include "Teuchos_RCP.hpp" 00040 #include "Teuchos_ScalarTraits.hpp" 00041 #include "Teuchos_GlobalMPISession.hpp" 00042 00043 using namespace std; 00044 using namespace Intrepid; 00045 00046 #define INTREPID_TEST_COMMAND( S ) \ 00047 { \ 00048 try { \ 00049 S ; \ 00050 } \ 00051 catch (std::logic_error err) { \ 00052 *outStream << "Expected Error ----------------------------------------------------------------\n"; \ 00053 *outStream << err.what() << '\n'; \ 00054 *outStream << "-------------------------------------------------------------------------------" << "\n\n"; \ 00055 }; \ 00056 } 00057 00058 00059 int main(int argc, char *argv[]) { 00060 00061 Teuchos::GlobalMPISession mpiSession(&argc, &argv); 00062 00063 // This little trick lets us print to std::cout only if 00064 // a (dummy) command-line argument is provided. 00065 int iprint = argc - 1; 00066 Teuchos::RCP<std::ostream> outStream; 00067 Teuchos::oblackholestream bhs; // outputs nothing 00068 if (iprint > 0) 00069 outStream = Teuchos::rcp(&std::cout, false); 00070 else 00071 outStream = Teuchos::rcp(&bhs, false); 00072 00073 // Save the format state of the original std::cout. 00074 Teuchos::oblackholestream oldFormatState; 00075 oldFormatState.copyfmt(std::cout); 00076 00077 *outStream \ 00078 << "===============================================================================\n" \ 00079 << "| |\n" \ 00080 << "| Unit Test (ArrayTools) |\n" \ 00081 << "| |\n" \ 00082 << "| 1) Array operations: dot multiply |\n" \ 00083 << "| |\n" \ 00084 << "| Questions? Contact Pavel Bochev (pbboche@sandia.gov) or |\n" \ 00085 << "| Denis Ridzal (dridzal@sandia.gov). |\n" \ 00086 << "| |\n" \ 00087 << "| Intrepid's website: http://trilinos.sandia.gov/packages/intrepid |\n" \ 00088 << "| Trilinos website: http://trilinos.sandia.gov |\n" \ 00089 << "| |\n" \ 00090 << "===============================================================================\n"; 00091 00092 00093 int errorFlag = 0; 00094 #ifdef HAVE_INTREPID_DEBUG 00095 int beginThrowNumber = TestForException_getThrowNumber(); 00096 int endThrowNumber = beginThrowNumber + 36; 00097 #endif 00098 00099 typedef ArrayTools art; 00100 typedef RealSpaceTools<double> rst; 00101 #ifdef HAVE_INTREPID_DEBUG 00102 ArrayTools atools; 00103 #endif 00104 00105 *outStream \ 00106 << "\n" 00107 << "===============================================================================\n"\ 00108 << "| TEST 1: exceptions |\n"\ 00109 << "===============================================================================\n"; 00110 00111 try{ 00112 00113 #ifdef HAVE_INTREPID_DEBUG 00114 FieldContainer<double> a_2(2); 00115 FieldContainer<double> a_2_2(2, 2); 00116 FieldContainer<double> a_3_2(3, 2); 00117 FieldContainer<double> a_2_2_2(2, 2, 2); 00118 FieldContainer<double> a_2_2_2_2(2, 2, 2, 2); 00119 FieldContainer<double> a_10_1(10, 1); 00120 FieldContainer<double> a_10_2(10, 2); 00121 FieldContainer<double> a_10_3(10, 3); 00122 FieldContainer<double> a_10_1_2(10, 1, 2); 00123 FieldContainer<double> a_10_2_2(10, 2, 2); 00124 FieldContainer<double> a_10_2_2_2(10, 2, 2, 2); 00125 FieldContainer<double> a_10_2_2_2_2(10, 2, 2, 2, 2); 00126 FieldContainer<double> a_9_2_2(9, 2, 2); 00127 FieldContainer<double> a_10_3_2(10, 3, 2); 00128 FieldContainer<double> a_10_2_3(10, 2, 3); 00129 FieldContainer<double> a_10_2_2_3(10, 2, 2, 3); 00130 FieldContainer<double> a_10_1_2_2(10, 1, 2, 2); 00131 00132 *outStream << "-> dotMultiplyDataField:\n"; 00133 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_2_2, a_2, a_2_2) ); 00134 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_2_2, a_10_2_2, a_10_2_2_2_2) ); 00135 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2, a_10_2_2, a_10_2_2_2) ); 00136 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_9_2_2, a_10_2_2_2) ); 00137 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_3_2, a_10_2_2_2) ); 00138 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_2_3, a_10_2_2_2) ); 00139 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_2_2_3, a_10_2_2_2_2) ); 00140 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_2_2_2, a_10_2_2_2, a_10_2_2_2_2) ); 00141 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_3_2, a_10_2_2_2, a_10_2_2_2_2) ); 00142 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_3, a_10_2_2_2, a_10_2_2_2_2) ); 00143 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_2_2_2, a_10_2_2_2_2) ); 00144 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_1_2_2, a_10_2_2_2_2) ); 00145 // 00146 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_2_2, a_2, a_2) ); 00147 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_2_2, a_10_2_2, a_10_2) ); 00148 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2, a_10_2_2, a_10_2_2) ); 00149 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_2_2, a_10_3_2) ); 00150 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_2_2, a_10_2_2) ); 00151 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_3, a_10_2_2, a_2_2_2) ); 00152 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_9_2_2, a_2_2_2) ); 00153 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_2_3, a_2_2_2) ); 00154 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_2_2_3, a_2_2_2_2) ); 00155 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_2_2_2, a_2_2_2_2) ); 00156 INTREPID_TEST_COMMAND( atools.dotMultiplyDataField<double>(a_10_2_2, a_10_1, a_2_2) ); 00157 00158 *outStream << "-> dotMultiplyDataData:\n"; 00159 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_2_2, a_2, a_2_2) ); 00160 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_2_2, a_10_2_2, a_10_2_2_2) ); 00161 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2_2, a_10_2_2, a_10_2_2) ); 00162 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_9_2_2, a_10_2_2) ); 00163 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_3_2, a_10_2_2) ); 00164 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2_3, a_10_2_2) ); 00165 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2_2_2, a_10_2_2_3) ); 00166 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_9_2_2, a_9_2_2) ); 00167 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_3_2, a_10_3_2) ); 00168 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2, a_10_2) ); 00169 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2_2, a_10_2_2) ); 00170 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2_2_2, a_10_2_2_2) ); 00171 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_1, a_10_2) ); 00172 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_1_2, a_10_2_2) ); 00173 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_1_2_2, a_10_2_2_2) ); 00174 // 00175 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_2_2, a_10_2_2_2_2, a_10_2) ); 00176 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_2_2, a_10_2_2, a_2) ); 00177 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2_2, a_10_2_2, a_10_2) ); 00178 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2_2, a_10_2) ); 00179 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_3, a_10_2_2, a_2_2) ); 00180 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_9_2_2, a_2_2) ); 00181 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2_3, a_2_2) ); 00182 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2_2_3, a_2_2_2) ); 00183 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2, a_2) ); 00184 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2_2, a_2_2) ); 00185 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_2_2_2, a_2_2_2) ); 00186 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_1, a_2) ); 00187 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_1_2, a_2_2) ); 00188 INTREPID_TEST_COMMAND( atools.dotMultiplyDataData<double>(a_10_2, a_10_1_2_2, a_2_2_2) ); 00189 #endif 00190 00191 } 00192 catch (std::logic_error err) { 00193 *outStream << "UNEXPECTED ERROR !!! ----------------------------------------------------------\n"; 00194 *outStream << err.what() << '\n'; 00195 *outStream << "-------------------------------------------------------------------------------" << "\n\n"; 00196 errorFlag = -1000; 00197 }; 00198 00199 #ifdef HAVE_INTREPID_DEBUG 00200 if (TestForException_getThrowNumber() != endThrowNumber) 00201 errorFlag++; 00202 #endif 00203 00204 *outStream \ 00205 << "\n" 00206 << "===============================================================================\n"\ 00207 << "| TEST 2: correctness of math operations |\n"\ 00208 << "===============================================================================\n"; 00209 00210 outStream->precision(20); 00211 00212 try { 00213 { // start scope 00214 *outStream << "\n************ Checking dotMultiplyDataField, (branch 1) ************\n"; 00215 00216 int c=5, p=9, f=7, d1=6, d2=14; 00217 00218 FieldContainer<double> in_c_f_p(c, f, p); 00219 FieldContainer<double> in_c_f_p_d(c, f, p, d1); 00220 FieldContainer<double> in_c_f_p_d_d(c, f, p, d1, d2); 00221 FieldContainer<double> data_c_p(c, p); 00222 FieldContainer<double> data_c_p_d(c, p, d1); 00223 FieldContainer<double> data_c_p_d_d(c, p, d1, d2); 00224 FieldContainer<double> data_c_1(c, 1); 00225 FieldContainer<double> data_c_1_d(c, 1, d1); 00226 FieldContainer<double> data_c_1_d_d(c, 1, d1, d2); 00227 FieldContainer<double> out_c_f_p(c, f, p); 00228 FieldContainer<double> outSM_c_f_p(c, f, p); 00229 FieldContainer<double> outDM_c_f_p(c, f, p); 00230 double zero = INTREPID_TOL*10000.0; 00231 00232 // fill with random numbers 00233 for (int i=0; i<in_c_f_p.size(); i++) { 00234 in_c_f_p[i] = Teuchos::ScalarTraits<double>::random(); 00235 } 00236 // fill with alternating 1's and -1's 00237 for (int i=0; i<in_c_f_p_d.size(); i++) { 00238 in_c_f_p_d[i] = std::pow((double)-1.0, (int)i); 00239 } 00240 // fill with 1's 00241 for (int i=0; i<in_c_f_p_d_d.size(); i++) { 00242 in_c_f_p_d_d[i] = 1.0; 00243 } 00244 // fill with random numbers 00245 for (int i=0; i<data_c_p.size(); i++) { 00246 data_c_p[i] = Teuchos::ScalarTraits<double>::random(); 00247 } 00248 // fill with 1's 00249 for (int i=0; i<data_c_1.size(); i++) { 00250 data_c_1[i] = Teuchos::ScalarTraits<double>::random(); 00251 } 00252 // fill with 1's 00253 for (int i=0; i<data_c_p_d.size(); i++) { 00254 data_c_p_d[i] = 1.0; 00255 } 00256 // fill with 1's 00257 for (int i=0; i<data_c_1_d.size(); i++) { 00258 data_c_1_d[i] = 1.0; 00259 } 00260 // fill with 1's 00261 for (int i=0; i<data_c_p_d_d.size(); i++) { 00262 data_c_p_d_d[i] = 1.0; 00263 } 00264 // fill with 1's 00265 for (int i=0; i<data_c_1_d_d.size(); i++) { 00266 data_c_1_d_d[i] = 1.0; 00267 } 00268 00269 art::scalarMultiplyDataField<double>(outSM_c_f_p, data_c_p, in_c_f_p); 00270 art::dotMultiplyDataField<double>(outDM_c_f_p, data_c_p, in_c_f_p); 00271 rst::subtract(out_c_f_p, outSM_c_f_p, outDM_c_f_p); 00272 if (rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_ONE) > zero) { 00273 *outStream << "\n\nINCORRECT dotMultiplyDataField (1): check dot multiply for scalars vs. scalar multiply\n\n"; 00274 errorFlag = -1000; 00275 } 00276 art::dotMultiplyDataField<double>(out_c_f_p, data_c_p_d, in_c_f_p_d); 00277 if (rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_ONE) > zero) { 00278 *outStream << "\n\nINCORRECT dotMultiplyDataField (2): check dot multiply of orthogonal vectors\n\n"; 00279 errorFlag = -1000; 00280 } 00281 art::dotMultiplyDataField<double>(out_c_f_p, data_c_p_d_d, in_c_f_p_d_d); 00282 if ((rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_INF) - d1*d2) > zero) { 00283 *outStream << "\n\nINCORRECT dotMultiplyDataField (3): check dot multiply for tensors of 1s\n\n"; 00284 errorFlag = -1000; 00285 } 00286 art::scalarMultiplyDataField<double>(outSM_c_f_p, data_c_1, in_c_f_p); 00287 art::dotMultiplyDataField<double>(outDM_c_f_p, data_c_1, in_c_f_p); 00288 rst::subtract(out_c_f_p, outSM_c_f_p, outDM_c_f_p); 00289 if (rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_ONE) > zero) { 00290 *outStream << "\n\nINCORRECT dotMultiplyDataField (4): check dot multiply for scalars vs. scalar multiply\n\n"; 00291 errorFlag = -1000; 00292 } 00293 art::dotMultiplyDataField<double>(out_c_f_p, data_c_1_d, in_c_f_p_d); 00294 if (rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_ONE) > zero) { 00295 *outStream << "\n\nINCORRECT dotMultiplyDataField (5): check dot multiply of orthogonal vectors\n\n"; 00296 errorFlag = -1000; 00297 } 00298 art::dotMultiplyDataField<double>(out_c_f_p, data_c_1_d_d, in_c_f_p_d_d); 00299 if ((rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_INF) - d1*d2) > zero) { 00300 *outStream << "\n\nINCORRECT dotMultiplyDataField (6): check dot multiply for tensors of 1s\n\n"; 00301 errorFlag = -1000; 00302 } 00303 } // end scope 00304 00305 00306 { // start scope 00307 *outStream << "\n************ Checking dotMultiplyDataField, (branch 2) ************\n"; 00308 00309 int c=5, p=9, f=7, d1=6, d2=14; 00310 00311 FieldContainer<double> in_f_p(f, p); 00312 FieldContainer<double> in_f_p_d(f, p, d1); 00313 FieldContainer<double> in_f_p_d_d(f, p, d1, d2); 00314 FieldContainer<double> data_c_p(c, p); 00315 FieldContainer<double> data_c_p_d(c, p, d1); 00316 FieldContainer<double> data_c_p_d_d(c, p, d1, d2); 00317 FieldContainer<double> data_c_1(c, 1); 00318 FieldContainer<double> data_c_1_d(c, 1, d1); 00319 FieldContainer<double> data_c_1_d_d(c, 1, d1, d2); 00320 FieldContainer<double> out_c_f_p(c, f, p); 00321 FieldContainer<double> outSM_c_f_p(c, f, p); 00322 FieldContainer<double> outDM_c_f_p(c, f, p); 00323 double zero = INTREPID_TOL*10000.0; 00324 00325 // fill with random numbers 00326 for (int i=0; i<in_f_p.size(); i++) { 00327 in_f_p[i] = Teuchos::ScalarTraits<double>::random(); 00328 } 00329 // fill with alternating 1's and -1's 00330 for (int i=0; i<in_f_p_d.size(); i++) { 00331 in_f_p_d[i] = std::pow((double)-1.0, (int)i); 00332 } 00333 // fill with 1's 00334 for (int i=0; i<in_f_p_d_d.size(); i++) { 00335 in_f_p_d_d[i] = 1.0; 00336 } 00337 // fill with random numbers 00338 for (int i=0; i<data_c_p.size(); i++) { 00339 data_c_p[i] = Teuchos::ScalarTraits<double>::random(); 00340 } 00341 // fill with 1's 00342 for (int i=0; i<data_c_1.size(); i++) { 00343 data_c_1[i] = Teuchos::ScalarTraits<double>::random(); 00344 } 00345 // fill with 1's 00346 for (int i=0; i<data_c_p_d.size(); i++) { 00347 data_c_p_d[i] = 1.0; 00348 } 00349 // fill with 1's 00350 for (int i=0; i<data_c_1_d.size(); i++) { 00351 data_c_1_d[i] = 1.0; 00352 } 00353 // fill with 1's 00354 for (int i=0; i<data_c_p_d_d.size(); i++) { 00355 data_c_p_d_d[i] = 1.0; 00356 } 00357 // fill with 1's 00358 for (int i=0; i<data_c_1_d_d.size(); i++) { 00359 data_c_1_d_d[i] = 1.0; 00360 } 00361 00362 art::scalarMultiplyDataField<double>(outSM_c_f_p, data_c_p, in_f_p); 00363 art::dotMultiplyDataField<double>(outDM_c_f_p, data_c_p, in_f_p); 00364 rst::subtract(out_c_f_p, outSM_c_f_p, outDM_c_f_p); 00365 if (rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_ONE) > zero) { 00366 *outStream << "\n\nINCORRECT dotMultiplyDataField (7): check dot multiply for scalars vs. scalar multiply\n\n"; 00367 errorFlag = -1000; 00368 } 00369 art::dotMultiplyDataField<double>(out_c_f_p, data_c_p_d, in_f_p_d); 00370 if (rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_ONE) > zero) { 00371 *outStream << "\n\nINCORRECT dotMultiplyDataField (8): check dot multiply of orthogonal vectors\n\n"; 00372 errorFlag = -1000; 00373 } 00374 art::dotMultiplyDataField<double>(out_c_f_p, data_c_p_d_d, in_f_p_d_d); 00375 if ((rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_INF) - d1*d2) > zero) { 00376 *outStream << "\n\nINCORRECT dotMultiplyDataField (9): check dot multiply for tensors of 1s\n\n"; 00377 errorFlag = -1000; 00378 } 00379 art::scalarMultiplyDataField<double>(outSM_c_f_p, data_c_1, in_f_p); 00380 art::dotMultiplyDataField<double>(outDM_c_f_p, data_c_1, in_f_p); 00381 rst::subtract(out_c_f_p, outSM_c_f_p, outDM_c_f_p); 00382 if (rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_ONE) > zero) { 00383 *outStream << "\n\nINCORRECT dotMultiplyDataField (10): check dot multiply for scalars vs. scalar multiply\n\n"; 00384 errorFlag = -1000; 00385 } 00386 art::dotMultiplyDataField<double>(out_c_f_p, data_c_1_d, in_f_p_d); 00387 if (rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_ONE) > zero) { 00388 *outStream << "\n\nINCORRECT dotMultiplyDataField (11): check dot multiply of orthogonal vectors\n\n"; 00389 errorFlag = -1000; 00390 } 00391 art::dotMultiplyDataField<double>(out_c_f_p, data_c_1_d_d, in_f_p_d_d); 00392 if ((rst::vectorNorm(&out_c_f_p[0], out_c_f_p.size(), NORM_INF) - d1*d2) > zero) { 00393 *outStream << "\n\nINCORRECT dotMultiplyDataField (12): check dot multiply for tensors of 1s\n\n"; 00394 errorFlag = -1000; 00395 } 00396 } // end scope 00397 00398 00399 00400 00401 { // start scope 00402 *outStream << "\n************ Checking dotMultiplyDataData, (branch 1) ************\n"; 00403 00404 int c=5, p=9, d1=6, d2=14; 00405 00406 FieldContainer<double> in_c_p(c, p); 00407 FieldContainer<double> in_c_p_d(c, p, d1); 00408 FieldContainer<double> in_c_p_d_d(c, p, d1, d2); 00409 FieldContainer<double> data_c_p(c, p); 00410 FieldContainer<double> data_c_p_d(c, p, d1); 00411 FieldContainer<double> data_c_p_d_d(c, p, d1, d2); 00412 FieldContainer<double> data_c_1(c, 1); 00413 FieldContainer<double> data_c_1_d(c, 1, d1); 00414 FieldContainer<double> data_c_1_d_d(c, 1, d1, d2); 00415 FieldContainer<double> out_c_p(c, p); 00416 FieldContainer<double> outSM_c_p(c, p); 00417 FieldContainer<double> outDM_c_p(c, p); 00418 double zero = INTREPID_TOL*10000.0; 00419 00420 // fill with random numbers 00421 for (int i=0; i<in_c_p.size(); i++) { 00422 in_c_p[i] = Teuchos::ScalarTraits<double>::random(); 00423 } 00424 // fill with alternating 1's and -1's 00425 for (int i=0; i<in_c_p_d.size(); i++) { 00426 in_c_p_d[i] = std::pow((double)-1.0, (int)i); 00427 } 00428 // fill with 1's 00429 for (int i=0; i<in_c_p_d_d.size(); i++) { 00430 in_c_p_d_d[i] = 1.0; 00431 } 00432 // fill with random numbers 00433 for (int i=0; i<data_c_p.size(); i++) { 00434 data_c_p[i] = Teuchos::ScalarTraits<double>::random(); 00435 } 00436 // fill with 1's 00437 for (int i=0; i<data_c_1.size(); i++) { 00438 data_c_1[i] = Teuchos::ScalarTraits<double>::random(); 00439 } 00440 // fill with 1's 00441 for (int i=0; i<data_c_p_d.size(); i++) { 00442 data_c_p_d[i] = 1.0; 00443 } 00444 // fill with 1's 00445 for (int i=0; i<data_c_1_d.size(); i++) { 00446 data_c_1_d[i] = 1.0; 00447 } 00448 // fill with 1's 00449 for (int i=0; i<data_c_p_d_d.size(); i++) { 00450 data_c_p_d_d[i] = 1.0; 00451 } 00452 // fill with 1's 00453 for (int i=0; i<data_c_1_d_d.size(); i++) { 00454 data_c_1_d_d[i] = 1.0; 00455 } 00456 00457 art::scalarMultiplyDataData<double>(outSM_c_p, data_c_p, in_c_p); 00458 art::dotMultiplyDataData<double>(outDM_c_p, data_c_p, in_c_p); 00459 rst::subtract(out_c_p, outSM_c_p, outDM_c_p); 00460 if (rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_ONE) > zero) { 00461 *outStream << "\n\nINCORRECT dotMultiplyDataData (1): check dot multiply for scalars vs. scalar multiply\n\n"; 00462 errorFlag = -1000; 00463 } 00464 art::dotMultiplyDataData<double>(out_c_p, data_c_p_d, in_c_p_d); 00465 if (rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_ONE) > zero) { 00466 *outStream << "\n\nINCORRECT dotMultiplyDataData (2): check dot multiply of orthogonal vectors\n\n"; 00467 errorFlag = -1000; 00468 } 00469 art::dotMultiplyDataData<double>(out_c_p, data_c_p_d_d, in_c_p_d_d); 00470 if ((rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_INF) - d1*d2) > zero) { 00471 *outStream << "\n\nINCORRECT dotMultiplyDataData (3): check dot multiply for tensors of 1s\n\n"; 00472 errorFlag = -1000; 00473 } 00474 art::scalarMultiplyDataData<double>(outSM_c_p, data_c_1, in_c_p); 00475 art::dotMultiplyDataData<double>(outDM_c_p, data_c_1, in_c_p); 00476 rst::subtract(out_c_p, outSM_c_p, outDM_c_p); 00477 if (rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_ONE) > zero) { 00478 *outStream << "\n\nINCORRECT dotMultiplyDataData (4): check dot multiply for scalars vs. scalar multiply\n\n"; 00479 errorFlag = -1000; 00480 } 00481 art::dotMultiplyDataData<double>(out_c_p, data_c_1_d, in_c_p_d); 00482 if (rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_ONE) > zero) { 00483 *outStream << "\n\nINCORRECT dotMultiplyDataData (5): check dot multiply of orthogonal vectors\n\n"; 00484 errorFlag = -1000; 00485 } 00486 art::dotMultiplyDataData<double>(out_c_p, data_c_1_d_d, in_c_p_d_d); 00487 if ((rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_INF) - d1*d2) > zero) { 00488 *outStream << "\n\nINCORRECT dotMultiplyDataData (6): check dot multiply for tensors of 1s\n\n"; 00489 errorFlag = -1000; 00490 } 00491 } // end scope 00492 00493 00494 { // start scope 00495 *outStream << "\n************ Checking dotMultiplyDataData, (branch 2) ************\n"; 00496 00497 int c=5, p=9, d1=6, d2=14; 00498 00499 FieldContainer<double> in_p(p); 00500 FieldContainer<double> in_p_d(p, d1); 00501 FieldContainer<double> in_p_d_d(p, d1, d2); 00502 FieldContainer<double> data_c_p(c, p); 00503 FieldContainer<double> data_c_p_d(c, p, d1); 00504 FieldContainer<double> data_c_p_d_d(c, p, d1, d2); 00505 FieldContainer<double> data_c_1(c, 1); 00506 FieldContainer<double> data_c_1_d(c, 1, d1); 00507 FieldContainer<double> data_c_1_d_d(c, 1, d1, d2); 00508 FieldContainer<double> out_c_p(c, p); 00509 FieldContainer<double> outSM_c_p(c, p); 00510 FieldContainer<double> outDM_c_p(c, p); 00511 double zero = INTREPID_TOL*10000.0; 00512 00513 // fill with random numbers 00514 for (int i=0; i<in_p.size(); i++) { 00515 in_p[i] = Teuchos::ScalarTraits<double>::random(); 00516 } 00517 // fill with alternating 1's and -1's 00518 for (int i=0; i<in_p_d.size(); i++) { 00519 in_p_d[i] = std::pow((double)-1.0, (int)i); 00520 } 00521 // fill with 1's 00522 for (int i=0; i<in_p_d_d.size(); i++) { 00523 in_p_d_d[i] = 1.0; 00524 } 00525 // fill with random numbers 00526 for (int i=0; i<data_c_p.size(); i++) { 00527 data_c_p[i] = Teuchos::ScalarTraits<double>::random(); 00528 } 00529 // fill with 1's 00530 for (int i=0; i<data_c_1.size(); i++) { 00531 data_c_1[i] = Teuchos::ScalarTraits<double>::random(); 00532 } 00533 // fill with 1's 00534 for (int i=0; i<data_c_p_d.size(); i++) { 00535 data_c_p_d[i] = 1.0; 00536 } 00537 // fill with 1's 00538 for (int i=0; i<data_c_1_d.size(); i++) { 00539 data_c_1_d[i] = 1.0; 00540 } 00541 // fill with 1's 00542 for (int i=0; i<data_c_p_d_d.size(); i++) { 00543 data_c_p_d_d[i] = 1.0; 00544 } 00545 // fill with 1's 00546 for (int i=0; i<data_c_1_d_d.size(); i++) { 00547 data_c_1_d_d[i] = 1.0; 00548 } 00549 00550 art::scalarMultiplyDataData<double>(outSM_c_p, data_c_p, in_p); 00551 art::dotMultiplyDataData<double>(outDM_c_p, data_c_p, in_p); 00552 rst::subtract(out_c_p, outSM_c_p, outDM_c_p); 00553 if (rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_ONE) > zero) { 00554 *outStream << "\n\nINCORRECT dotMultiplyDataData (7): check dot multiply for scalars vs. scalar multiply\n\n"; 00555 errorFlag = -1000; 00556 } 00557 art::dotMultiplyDataData<double>(out_c_p, data_c_p_d, in_p_d); 00558 if (rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_ONE) > zero) { 00559 *outStream << "\n\nINCORRECT dotMultiplyDataData (8): check dot multiply of orthogonal vectors\n\n"; 00560 errorFlag = -1000; 00561 } 00562 art::dotMultiplyDataData<double>(out_c_p, data_c_p_d_d, in_p_d_d); 00563 if ((rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_INF) - d1*d2) > zero) { 00564 *outStream << "\n\nINCORRECT dotMultiplyDataData (9): check dot multiply for tensors of 1s\n\n"; 00565 errorFlag = -1000; 00566 } 00567 art::scalarMultiplyDataData<double>(outSM_c_p, data_c_1, in_p); 00568 art::dotMultiplyDataData<double>(outDM_c_p, data_c_1, in_p); 00569 rst::subtract(out_c_p, outSM_c_p, outDM_c_p); 00570 if (rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_ONE) > zero) { 00571 *outStream << "\n\nINCORRECT dotMultiplyDataData (10): check dot multiply for scalars vs. scalar multiply\n\n"; 00572 errorFlag = -1000; 00573 } 00574 art::dotMultiplyDataData<double>(out_c_p, data_c_1_d, in_p_d); 00575 if (rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_ONE) > zero) { 00576 *outStream << "\n\nINCORRECT dotMultiplyDataData (11): check dot multiply of orthogonal vectors\n\n"; 00577 errorFlag = -1000; 00578 } 00579 art::dotMultiplyDataData<double>(out_c_p, data_c_1_d_d, in_p_d_d); 00580 if ((rst::vectorNorm(&out_c_p[0], out_c_p.size(), NORM_INF) - d1*d2) > zero) { 00581 *outStream << "\n\nINCORRECT dotMultiplyDataData (12): check dot multiply for tensors of 1s\n\n"; 00582 errorFlag = -1000; 00583 } 00584 } // end scope 00585 00586 /******************************************/ 00587 *outStream << "\n"; 00588 } 00589 catch (std::logic_error err) { 00590 *outStream << "UNEXPECTED ERROR !!! ----------------------------------------------------------\n"; 00591 *outStream << err.what() << '\n'; 00592 *outStream << "-------------------------------------------------------------------------------" << "\n\n"; 00593 errorFlag = -1000; 00594 }; 00595 00596 00597 if (errorFlag != 0) 00598 std::cout << "End Result: TEST FAILED\n"; 00599 else 00600 std::cout << "End Result: TEST PASSED\n"; 00601 00602 // reset format state of std::cout 00603 std::cout.copyfmt(oldFormatState); 00604 00605 return errorFlag; 00606 }
1.7.4