FEI Version of the Day
beam_main.cpp
00001 /*--------------------------------------------------------------------*/
00002 /*    Copyright 2005 Sandia Corporation.                              */
00003 /*    Under the terms of Contract DE-AC04-94AL85000, there is a       */
00004 /*    non-exclusive license for use of this work by or on behalf      */
00005 /*    of the U.S. Government.  Export of this program may require     */
00006 /*    a license from the United States Government.                    */
00007 /*--------------------------------------------------------------------*/
00008 
00009 //
00010 // This is a simple program to exercise FEI classes,
00011 // for the purposes of testing, code tuning and scaling studies.
00012 //
00013 
00014 #include <fei_macros.hpp> //includes std stuff like iostream, etc.
00015 
00016 //Including the header fei_base.hpp gets us the declaration for
00017 //the various fei:: and snl_fei:: classes.
00018 
00019 #include <fei_base.hpp>
00020 
00021 //Now make provision for using any one of several solver libraries. This is
00022 //handled by the code in LibraryFactory.{hpp,cpp}.
00023 
00024 #include <test_utils/LibraryFactory.hpp>
00025 
00026 //And finally, we need to include some 'utils' headers which are used for
00027 //setting up the data for the example problem.
00028 
00029 #include <test_utils/HexBeam.hpp>
00030 #include <test_utils/HexBeamCR.hpp>
00031 #include <snl_fei_Utils.hpp>
00032 #include <test_utils/fei_test_utils.hpp>
00033 //
00034 //Include definitions of macros to call functions and check the return code.
00035 //
00036 #undef fei_file
00037 #define fei_file "cube3.cpp"
00038 #include <fei_ErrMacros.hpp>
00039 
00040 //==============================================================================
00041 //Here's the main...
00042 //==============================================================================
00043 int beam_main(int argc, char** argv,
00044          MPI_Comm comm, int numProcs, int localProc){
00045 
00046   double start_time = fei::utils::cpu_time();
00047 
00048   std::vector<std::string> stdstrings;
00049   CHK_ERR( fei_test_utils::get_filename_and_read_input(argc, argv,
00050             comm, localProc,
00051             stdstrings) );
00052   fei::ParameterSet params;
00053   fei::utils::parse_strings(stdstrings, " ", params);
00054 
00055   std::string solverName;
00056   std::string datasource;
00057   std::string constraintform;
00058   int W = 0;
00059   int D = 0;
00060   int DofPerNode = 0;
00061 
00062   int errcode = 0;
00063   errcode += params.getStringParamValue("SOLVER_LIBRARY", solverName);
00064   errcode += params.getStringParamValue("DATA_SOURCE", datasource);
00065   errcode += params.getIntParamValue("W", W);
00066   errcode += params.getIntParamValue("D", D);
00067   errcode += params.getIntParamValue("DofPerNode", DofPerNode);
00068 
00069   if (errcode != 0) {
00070     fei::console_out() << "Failed to find one or more required parameters in input-file."
00071        << FEI_ENDL << "Required parameters:"<<FEI_ENDL
00072        << "SOLVER_LIBRARY" << FEI_ENDL
00073        << "DATA_SOURCE" << FEI_ENDL
00074        << "CONSTRAINT_FORM" << FEI_ENDL
00075        << "W" << FEI_ENDL << "D" << FEI_ENDL << "DofPerNode" << FEI_ENDL;
00076     return(-1);
00077   }
00078 
00079   params.getStringParamValue("CONSTRAINT_FORM", constraintform);
00080 
00081   bool slave_constraints = false;
00082   if ("slaves" == constraintform) {
00083     slave_constraints = true;
00084   }
00085 
00086   HexBeam* hexcubeptr = NULL;
00087   if (datasource == "HexBeam") {
00088     hexcubeptr = new HexBeam(W, D, DofPerNode,
00089            HexBeam::OneD, numProcs, localProc);
00090   }
00091   else {
00092     hexcubeptr = new HexBeamCR(W, D, DofPerNode,
00093              HexBeam::OneD, numProcs, localProc);
00094   }
00095 
00096   HexBeam& hexcube = *hexcubeptr;
00097 
00098   if (localProc == 0) {
00099     int numCRs = (W+1)*(W+1)* ((numProcs*2)-1);
00100     if (hexcube.getNumCRs() < 1) numCRs = 0;
00101     FEI_COUT << FEI_ENDL;
00102     FEI_COUT << "========================================================" 
00103    << FEI_ENDL;
00104     FEI_COUT << "FEI version: " << fei::utils::version() << FEI_ENDL;
00105     FEI_COUT << "--------------------------------------------------------"<<FEI_ENDL;
00106     FEI_COUT << "Size W: " << W << " (num-elements-along-side-of-cube)"<<FEI_ENDL;
00107     FEI_COUT << "Size D: " << D << " (num-elements-along-depth-of-cube)"<<FEI_ENDL;
00108     FEI_COUT << "DOF per node: " << DofPerNode <<FEI_ENDL;
00109     FEI_COUT << "Num local  elements: " << hexcube.localNumElems_ << FEI_ENDL;
00110     FEI_COUT << "Num global elements: " << hexcube.totalNumElems_ << FEI_ENDL;
00111     FEI_COUT << "Num local  DOF: " << hexcube.numLocalDOF_ << FEI_ENDL;
00112     FEI_COUT << "Num global DOF: " << hexcube.numGlobalDOF_ << FEI_ENDL;
00113     FEI_COUT << "Num global CRs: " << numCRs << FEI_ENDL;
00114     FEI_COUT << "========================================================" 
00115    << FEI_ENDL;
00116   }
00117 
00118   double start_init_time = fei::utils::cpu_time();
00119 
00120   fei::SharedPtr<fei::Factory> factory;
00121   try {
00122     factory = fei::create_fei_Factory(comm, solverName.c_str());
00123   }
00124   catch (...) {
00125     FEI_COUT << "library " << solverName << " not available."<<FEI_ENDL;
00126     return(-1);
00127   }
00128   if (factory.get() == NULL) {
00129     FEI_COUT << "snl_fei::Factory allocation failed." << FEI_ENDL;
00130     return(-1);
00131   }
00132 
00133   factory->parameters(params);
00134 
00135   fei::SharedPtr<fei::VectorSpace> nodeSpace =
00136     factory->createVectorSpace(comm, "cube3");
00137 
00138   fei::SharedPtr<fei::VectorSpace> dummy;
00139   fei::SharedPtr<fei::MatrixGraph> matrixGraph =
00140     factory->createMatrixGraph(nodeSpace, dummy, "cube3");
00141 
00142   //load some solver-control parameters.
00143   nodeSpace->setParameters(params);
00144   matrixGraph->setParameters(params);
00145 
00146   int fieldID = 0;
00147   int fieldSize = hexcube.numDofPerNode();
00148   int nodeIDType = 0;
00149   int crIDType = 1;
00150 
00151   nodeSpace->defineFields( 1, &fieldID, &fieldSize );
00152   nodeSpace->defineIDTypes(1, &nodeIDType );
00153   nodeSpace->defineIDTypes(1, &crIDType );
00154 
00155   CHK_ERR( HexBeam_Functions::
00156      init_elem_connectivities( matrixGraph.get(), hexcube ) );
00157 
00158   CHK_ERR( HexBeam_Functions::
00159      init_shared_nodes( matrixGraph.get(), hexcube ) );
00160 
00161   int firstLocalCRID = 0;
00162   if (slave_constraints) {
00163     CHK_ERR( HexBeam_Functions::
00164        init_slave_constraints( matrixGraph.get(), hexcube) );
00165   }
00166   else {
00167     CHK_ERR( HexBeam_Functions::
00168        init_constraints( matrixGraph.get(), hexcube, localProc, firstLocalCRID) );
00169   }
00170 
00171   CHK_ERR( matrixGraph->initComplete() );
00172 
00173   double fei_init_time = fei::utils::cpu_time() - start_init_time;
00174 
00175   if (localProc == 0) {
00176     FEI_COUT.setf(IOS_FIXED, IOS_FLOATFIELD);
00177     FEI_COUT << "Initialization cpu time:   " << fei_init_time << FEI_ENDL;
00178   }
00179 
00180   //Now the initialization phase is complete. Next we'll do the load phase,
00181   //which for this problem just consists of loading the element data
00182   //(element-wise stiffness arrays and load vectors) and the boundary
00183   //condition data.
00184 
00185   double fei_creatematrix_start_time = fei::utils::cpu_time();
00186 
00187   fei::SharedPtr<fei::Matrix> mat = factory->createMatrix(matrixGraph);
00188 
00189   double fei_creatematrix_time = fei::utils::cpu_time() - fei_creatematrix_start_time;
00190   if (localProc == 0) {
00191     FEI_COUT.setf(IOS_FIXED, IOS_FLOATFIELD);
00192     FEI_COUT << "Create-Matrix cpu time:   " << fei_creatematrix_time << FEI_ENDL;
00193   }
00194 
00195   double start_load_time = fei::utils::cpu_time();
00196 
00197   fei::SharedPtr<fei::Vector> solnVec = factory->createVector(matrixGraph, true);
00198   fei::SharedPtr<fei::Vector> rhsVec  = factory->createVector(matrixGraph);
00199   fei::SharedPtr<fei::LinearSystem> linSys =
00200     factory->createLinearSystem(matrixGraph);
00201 
00202   CHK_ERR( linSys->parameters(params));
00203 
00204   linSys->setMatrix(mat);
00205   linSys->setSolutionVector(solnVec);
00206   linSys->setRHS(rhsVec);
00207 
00208   CHK_ERR( HexBeam_Functions::
00209      load_elem_data(matrixGraph.get(), mat.get(), rhsVec.get(), hexcube) );
00210 
00211   //temporarily disable boundary-conditions
00212   //  CHK_ERR( load_BC_data(linSys, hexcube) );
00213 
00214   if (!slave_constraints) {
00215     CHK_ERR( HexBeam_Functions::
00216        load_constraints(linSys.get(), hexcube, firstLocalCRID) );
00217   }
00218 
00219   CHK_ERR( linSys->loadComplete() );
00220 
00221   double fei_load_time = fei::utils::cpu_time() - start_load_time;
00222 
00223   if (localProc == 0) {
00224     //IOS macros are defined in fei_macros.h
00225     FEI_COUT.setf(IOS_FIXED, IOS_FLOATFIELD);
00226     FEI_COUT << "Coef. loading cpu time:    " << fei_load_time << FEI_ENDL;
00227     FEI_COUT << "Total assembly wall time:   "
00228    << fei_init_time + fei_creatematrix_time + fei_load_time << FEI_ENDL;
00229   }
00230 
00231   //
00232   //now the load phase is complete, so we're ready to launch the underlying
00233   //solver and solve Ax=b
00234   //
00235 
00236   //First, check whether the params (set from the input .i file) specify
00237   //a "Trilinos_Solver" string (which may be Amesos...). If not, it won't
00238   //matter but if so, the factory may switch based on this value, when
00239   //creating the solver object instance.
00240 
00241   std::string solver_name_value;
00242   params.getStringParamValue("Trilinos_Solver", solver_name_value);
00243 
00244   const char* charptr_solvername =
00245     solver_name_value.empty() ? 0 : solver_name_value.c_str();
00246 
00247   fei::SharedPtr<fei::Solver> solver = factory->createSolver(charptr_solvername);
00248 
00249   int status;
00250   int itersTaken = 0;
00251 
00252   if (localProc==0) FEI_COUT << "solve..." << FEI_ENDL;
00253   double start_solve_time = fei::utils::cpu_time();
00254 
00255   int err = solver->solve(linSys.get(),
00256         NULL, //preconditioningMatrix
00257         params,
00258         itersTaken,
00259         status);
00260 
00261   double solve_time = fei::utils::cpu_time()-start_solve_time;
00262 
00263   if (err!=0) {
00264     if (localProc==0) FEI_COUT << "solve returned err: " << err <<", status: "
00265          << status << FEI_ENDL;
00266   }
00267 
00268   if (localProc==0) {
00269     FEI_COUT << " cpu-time in solve: " << solve_time << FEI_ENDL;
00270   }
00271 
00272   CHK_ERR( solnVec->scatterToOverlap() );
00273 
00274   delete hexcubeptr;
00275 
00276   double elapsed_cpu_time = fei::utils::cpu_time() - start_time;
00277   int returnValue = 0;
00278 
00279   //The following IOS_... macros are defined in base/fei_macros.h
00280   FEI_COUT.setf(IOS_FIXED, IOS_FLOATFIELD);
00281   if (localProc==0) {
00282     FEI_COUT << "Proc0 cpu times (seconds):" << FEI_ENDL
00283    << "   FEI initialize:    " << fei_init_time << FEI_ENDL
00284          << "   FEI create-matrix: " << fei_creatematrix_time << FEI_ENDL
00285          << "   FEI load:          " << fei_load_time << FEI_ENDL
00286          << "      solve:          " << solve_time << FEI_ENDL
00287          << "Total program time:   " << elapsed_cpu_time << FEI_ENDL;
00288 
00289 #if defined(FEI_PLATFORM) && defined(FEI_OPT_LEVEL)
00290     double benchmark = fei_init_time;
00291 
00292     std::string slavestr;
00293     if (!constraintform.empty()) {
00294       slavestr = constraintform;
00295     }
00296     if (slavestr.size() > 0) slavestr += "_";
00297 
00298     FEI_OSTRINGSTREAM testname_init;
00299     testname_init << "cube3_init_"<<slavestr<<W<<"_"<<D<<"_"<<DofPerNode<<"_"
00300       <<solverName<<"_np"<<numProcs<<"_"
00301       <<FEI_PLATFORM<<"_"<<FEI_OPT_LEVEL;
00302     FEI_OSTRINGSTREAM testname_create;
00303     testname_create << "cube3_creatematrix_"<<slavestr<<W<<"_"<<D<<"_"<<DofPerNode
00304         <<"_"<<solverName<<"_np"<<numProcs<<"_"
00305         <<FEI_PLATFORM<<"_"<<FEI_OPT_LEVEL;
00306     FEI_OSTRINGSTREAM testname_load;
00307     testname_load << "cube3_load_"<<slavestr<<W<<"_"<<D<<"_"<<DofPerNode<<"_"
00308       <<solverName<<"_np"<<numProcs<<"_"
00309       <<FEI_PLATFORM<<"_"<<FEI_OPT_LEVEL;
00310 
00311     double file_init, file_create, file_load;
00312     bool file_benchmarks_available = true;
00313     try {
00314       file_init = fei_test_utils::get_file_benchmark("./cube3_timings.txt",
00315                 testname_init.str().c_str());
00316       file_create = fei_test_utils::get_file_benchmark("./cube3_timings.txt",
00317             testname_create.str().c_str());
00318       file_load = fei_test_utils::get_file_benchmark("./cube3_timings.txt",
00319                 testname_load.str().c_str());
00320     }
00321     catch (std::runtime_error& exc) {
00322       file_benchmarks_available = false;
00323     }
00324 
00325     if (file_benchmarks_available) {
00326 
00327       bool init_test_passed =
00328   fei_test_utils::check_and_cout_test_result(testname_init.str(), benchmark,
00329               file_init, 10);
00330 
00331       benchmark = fei_creatematrix_time;
00332       bool create_test_passed =
00333   fei_test_utils::check_and_cout_test_result(testname_create.str(), benchmark,
00334               file_create, 10);
00335 
00336       benchmark = fei_load_time;
00337       bool load_test_passed =
00338   fei_test_utils::check_and_cout_test_result(testname_load.str(), benchmark,
00339               file_load, 10);
00340 
00341       returnValue = init_test_passed&&create_test_passed&&load_test_passed ? 0 : 1;
00342     }
00343 #endif
00344 
00345   }
00346 
00347   bool testPassed = returnValue==0;
00348   if (testPassed && localProc == 0) {
00349     //A string for the SIERRA runtest tool to search for in test output...
00350     FEI_COUT << "FEI test successful" << FEI_ENDL;
00351   }
00352 
00353   return(returnValue);
00354 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends