FEI Version of the Day
Loading...
Searching...
No Matches
beam_main.cpp
1/*--------------------------------------------------------------------*/
2/* Copyright 2005 Sandia Corporation. */
3/* Under the terms of Contract DE-AC04-94AL85000, there is a */
4/* non-exclusive license for use of this work by or on behalf */
5/* of the U.S. Government. Export of this program may require */
6/* a license from the United States Government. */
7/*--------------------------------------------------------------------*/
8
9//
10// This is a simple program to exercise FEI classes,
11// for the purposes of testing, code tuning and scaling studies.
12//
13
14#include <fei_macros.hpp> //includes std stuff like iostream, etc.
15
16//Including the header fei_base.hpp gets us the declaration for
17//the various fei:: and snl_fei:: classes.
18
19#include <fei_base.hpp>
20
21//Now make provision for using any one of several solver libraries. This is
22//handled by the code in LibraryFactory.{hpp,cpp}.
23
24#include <test_utils/LibraryFactory.hpp>
25
26//And finally, we need to include some 'utils' headers which are used for
27//setting up the data for the example problem.
28
29#include <test_utils/HexBeam.hpp>
30#include <test_utils/HexBeamCR.hpp>
31#include <snl_fei_Utils.hpp>
32#include <test_utils/fei_test_utils.hpp>
33//
34//Include definitions of macros to call functions and check the return code.
35//
36#undef fei_file
37#define fei_file "cube3.cpp"
38#include <fei_ErrMacros.hpp>
39
40//==============================================================================
41//Here's the main...
42//==============================================================================
43int beam_main(int argc, char** argv,
44 MPI_Comm comm, int numProcs, int localProc){
45
46 double start_time = fei::utils::cpu_time();
47
48 std::vector<std::string> stdstrings;
50 comm, localProc,
51 stdstrings) );
52 fei::ParameterSet params;
53 fei::utils::parse_strings(stdstrings, " ", params);
54
55 std::string solverName;
56 std::string datasource;
57 std::string constraintform;
58 int W = 0;
59 int D = 0;
60 int DofPerNode = 0;
61
62 int errcode = 0;
63 errcode += params.getStringParamValue("SOLVER_LIBRARY", solverName);
64 errcode += params.getStringParamValue("DATA_SOURCE", datasource);
65 errcode += params.getIntParamValue("W", W);
66 errcode += params.getIntParamValue("D", D);
67 errcode += params.getIntParamValue("DofPerNode", DofPerNode);
68
69 if (errcode != 0) {
70 fei::console_out() << "Failed to find one or more required parameters in input-file."
71 << FEI_ENDL << "Required parameters:"<<FEI_ENDL
72 << "SOLVER_LIBRARY" << FEI_ENDL
73 << "DATA_SOURCE" << FEI_ENDL
74 << "CONSTRAINT_FORM" << FEI_ENDL
75 << "W" << FEI_ENDL << "D" << FEI_ENDL << "DofPerNode" << FEI_ENDL;
76 return(-1);
77 }
78
79 params.getStringParamValue("CONSTRAINT_FORM", constraintform);
80
81 bool slave_constraints = false;
82 if ("slaves" == constraintform) {
83 slave_constraints = true;
84 }
85
86 HexBeam* hexcubeptr = NULL;
87 if (datasource == "HexBeam") {
88 hexcubeptr = new HexBeam(W, D, DofPerNode,
89 HexBeam::OneD, numProcs, localProc);
90 }
91 else {
92 hexcubeptr = new HexBeamCR(W, D, DofPerNode,
93 HexBeam::OneD, numProcs, localProc);
94 }
95
96 HexBeam& hexcube = *hexcubeptr;
97
98 if (localProc == 0) {
99 int numCRs = (W+1)*(W+1)* ((numProcs*2)-1);
100 if (hexcube.getNumCRs() < 1) numCRs = 0;
101 FEI_COUT << FEI_ENDL;
102 FEI_COUT << "========================================================"
103 << FEI_ENDL;
104 FEI_COUT << "FEI version: " << fei::utils::version() << FEI_ENDL;
105 FEI_COUT << "--------------------------------------------------------"<<FEI_ENDL;
106 FEI_COUT << "Size W: " << W << " (num-elements-along-side-of-cube)"<<FEI_ENDL;
107 FEI_COUT << "Size D: " << D << " (num-elements-along-depth-of-cube)"<<FEI_ENDL;
108 FEI_COUT << "DOF per node: " << DofPerNode <<FEI_ENDL;
109 FEI_COUT << "Num local elements: " << hexcube.localNumElems_ << FEI_ENDL;
110 FEI_COUT << "Num global elements: " << hexcube.totalNumElems_ << FEI_ENDL;
111 FEI_COUT << "Num local DOF: " << hexcube.numLocalDOF_ << FEI_ENDL;
112 FEI_COUT << "Num global DOF: " << hexcube.numGlobalDOF_ << FEI_ENDL;
113 FEI_COUT << "Num global CRs: " << numCRs << FEI_ENDL;
114 FEI_COUT << "========================================================"
115 << FEI_ENDL;
116 }
117
118 double start_init_time = fei::utils::cpu_time();
119
121 try {
122 factory = fei::create_fei_Factory(comm, solverName.c_str());
123 }
124 catch (...) {
125 FEI_COUT << "library " << solverName << " not available."<<FEI_ENDL;
126 return(-1);
127 }
128 if (factory.get() == NULL) {
129 FEI_COUT << "snl_fei::Factory allocation failed." << FEI_ENDL;
130 return(-1);
131 }
132
133 factory->parameters(params);
134
136 factory->createVectorSpace(comm, "cube3");
137
140 factory->createMatrixGraph(nodeSpace, dummy, "cube3");
141
142 //load some solver-control parameters.
143 nodeSpace->setParameters(params);
144 matrixGraph->setParameters(params);
145
146 int fieldID = 0;
147 int fieldSize = hexcube.numDofPerNode();
148 int nodeIDType = 0;
149 int crIDType = 1;
150
151 nodeSpace->defineFields( 1, &fieldID, &fieldSize );
152 nodeSpace->defineIDTypes(1, &nodeIDType );
153 nodeSpace->defineIDTypes(1, &crIDType );
154
155 CHK_ERR( HexBeam_Functions::
156 init_elem_connectivities( matrixGraph.get(), hexcube ) );
157
158 CHK_ERR( HexBeam_Functions::
159 init_shared_nodes( matrixGraph.get(), hexcube ) );
160
161 int firstLocalCRID = 0;
162 if (slave_constraints) {
163 CHK_ERR( HexBeam_Functions::
164 init_slave_constraints( matrixGraph.get(), hexcube) );
165 }
166 else {
167 CHK_ERR( HexBeam_Functions::
168 init_constraints( matrixGraph.get(), hexcube, localProc, firstLocalCRID) );
169 }
170
171 CHK_ERR( matrixGraph->initComplete() );
172
173 double fei_init_time = fei::utils::cpu_time() - start_init_time;
174
175 if (localProc == 0) {
176 FEI_COUT.setf(IOS_FIXED, IOS_FLOATFIELD);
177 FEI_COUT << "Initialization cpu time: " << fei_init_time << FEI_ENDL;
178 }
179
180 //Now the initialization phase is complete. Next we'll do the load phase,
181 //which for this problem just consists of loading the element data
182 //(element-wise stiffness arrays and load vectors) and the boundary
183 //condition data.
184
185 double fei_creatematrix_start_time = fei::utils::cpu_time();
186
187 fei::SharedPtr<fei::Matrix> mat = factory->createMatrix(matrixGraph);
188
189 double fei_creatematrix_time = fei::utils::cpu_time() - fei_creatematrix_start_time;
190 if (localProc == 0) {
191 FEI_COUT.setf(IOS_FIXED, IOS_FLOATFIELD);
192 FEI_COUT << "Create-Matrix cpu time: " << fei_creatematrix_time << FEI_ENDL;
193 }
194
195 double start_load_time = fei::utils::cpu_time();
196
197 fei::SharedPtr<fei::Vector> solnVec = factory->createVector(matrixGraph, true);
198 fei::SharedPtr<fei::Vector> rhsVec = factory->createVector(matrixGraph);
200 factory->createLinearSystem(matrixGraph);
201
202 CHK_ERR( linSys->parameters(params));
203
204 linSys->setMatrix(mat);
205 linSys->setSolutionVector(solnVec);
206 linSys->setRHS(rhsVec);
207
208 CHK_ERR( HexBeam_Functions::
209 load_elem_data(matrixGraph.get(), mat.get(), rhsVec.get(), hexcube) );
210
211 //temporarily disable boundary-conditions
212 // CHK_ERR( load_BC_data(linSys, hexcube) );
213
214 if (!slave_constraints) {
215 CHK_ERR( HexBeam_Functions::
216 load_constraints(linSys.get(), hexcube, firstLocalCRID) );
217 }
218
219 CHK_ERR( linSys->loadComplete() );
220
221 double fei_load_time = fei::utils::cpu_time() - start_load_time;
222
223 if (localProc == 0) {
224 //IOS macros are defined in fei_macros.h
225 FEI_COUT.setf(IOS_FIXED, IOS_FLOATFIELD);
226 FEI_COUT << "Coef. loading cpu time: " << fei_load_time << FEI_ENDL;
227 FEI_COUT << "Total assembly wall time: "
228 << fei_init_time + fei_creatematrix_time + fei_load_time << FEI_ENDL;
229 }
230
231 //
232 //now the load phase is complete, so we're ready to launch the underlying
233 //solver and solve Ax=b
234 //
235
236 //First, check whether the params (set from the input .i file) specify
237 //a "Trilinos_Solver" string (which may be Amesos...). If not, it won't
238 //matter but if so, the factory may switch based on this value, when
239 //creating the solver object instance.
240
241 std::string solver_name_value;
242 params.getStringParamValue("Trilinos_Solver", solver_name_value);
243
244 const char* charptr_solvername =
245 solver_name_value.empty() ? 0 : solver_name_value.c_str();
246
247 fei::SharedPtr<fei::Solver> solver = factory->createSolver(charptr_solvername);
248
249 int status;
250 int itersTaken = 0;
251
252 if (localProc==0) FEI_COUT << "solve..." << FEI_ENDL;
253 double start_solve_time = fei::utils::cpu_time();
254
255 int err = solver->solve(linSys.get(),
256 NULL, //preconditioningMatrix
257 params,
258 itersTaken,
259 status);
260
261 double solve_time = fei::utils::cpu_time()-start_solve_time;
262
263 if (err!=0) {
264 if (localProc==0) FEI_COUT << "solve returned err: " << err <<", status: "
265 << status << FEI_ENDL;
266 }
267
268 if (localProc==0) {
269 FEI_COUT << " cpu-time in solve: " << solve_time << FEI_ENDL;
270 }
271
272 CHK_ERR( solnVec->scatterToOverlap() );
273
274 delete hexcubeptr;
275
276 double elapsed_cpu_time = fei::utils::cpu_time() - start_time;
277 int returnValue = 0;
278
279 //The following IOS_... macros are defined in base/fei_macros.h
280 FEI_COUT.setf(IOS_FIXED, IOS_FLOATFIELD);
281 if (localProc==0) {
282 FEI_COUT << "Proc0 cpu times (seconds):" << FEI_ENDL
283 << " FEI initialize: " << fei_init_time << FEI_ENDL
284 << " FEI create-matrix: " << fei_creatematrix_time << FEI_ENDL
285 << " FEI load: " << fei_load_time << FEI_ENDL
286 << " solve: " << solve_time << FEI_ENDL
287 << "Total program time: " << elapsed_cpu_time << FEI_ENDL;
288
289#if defined(FEI_PLATFORM) && defined(FEI_OPT_LEVEL)
290 double benchmark = fei_init_time;
291
292 std::string slavestr;
293 if (!constraintform.empty()) {
294 slavestr = constraintform;
295 }
296 if (slavestr.size() > 0) slavestr += "_";
297
298 FEI_OSTRINGSTREAM testname_init;
299 testname_init << "cube3_init_"<<slavestr<<W<<"_"<<D<<"_"<<DofPerNode<<"_"
300 <<solverName<<"_np"<<numProcs<<"_"
301 <<FEI_PLATFORM<<"_"<<FEI_OPT_LEVEL;
302 FEI_OSTRINGSTREAM testname_create;
303 testname_create << "cube3_creatematrix_"<<slavestr<<W<<"_"<<D<<"_"<<DofPerNode
304 <<"_"<<solverName<<"_np"<<numProcs<<"_"
305 <<FEI_PLATFORM<<"_"<<FEI_OPT_LEVEL;
306 FEI_OSTRINGSTREAM testname_load;
307 testname_load << "cube3_load_"<<slavestr<<W<<"_"<<D<<"_"<<DofPerNode<<"_"
308 <<solverName<<"_np"<<numProcs<<"_"
309 <<FEI_PLATFORM<<"_"<<FEI_OPT_LEVEL;
310
311 double file_init, file_create, file_load;
312 bool file_benchmarks_available = true;
313 try {
314 file_init = fei_test_utils::get_file_benchmark("./cube3_timings.txt",
315 testname_init.str().c_str());
316 file_create = fei_test_utils::get_file_benchmark("./cube3_timings.txt",
317 testname_create.str().c_str());
318 file_load = fei_test_utils::get_file_benchmark("./cube3_timings.txt",
319 testname_load.str().c_str());
320 }
321 catch (std::runtime_error& exc) {
322 file_benchmarks_available = false;
323 }
324
325 if (file_benchmarks_available) {
326
327 bool init_test_passed =
328 fei_test_utils::check_and_cout_test_result(testname_init.str(), benchmark,
329 file_init, 10);
330
331 benchmark = fei_creatematrix_time;
332 bool create_test_passed =
333 fei_test_utils::check_and_cout_test_result(testname_create.str(), benchmark,
334 file_create, 10);
335
336 benchmark = fei_load_time;
337 bool load_test_passed =
338 fei_test_utils::check_and_cout_test_result(testname_load.str(), benchmark,
339 file_load, 10);
340
341 returnValue = init_test_passed&&create_test_passed&&load_test_passed ? 0 : 1;
342 }
343#endif
344
345 }
346
347 bool testPassed = returnValue==0;
348 if (testPassed && localProc == 0) {
349 //A string for the SIERRA runtest tool to search for in test output...
350 FEI_COUT << "FEI test successful" << FEI_ENDL;
351 }
352
353 return(returnValue);
354}
int getIntParamValue(const char *name, int &paramValue) const
int getStringParamValue(const char *name, std::string &paramValue) const
T * get() const
void parse_strings(std::vector< std::string > &stdstrings, const char *separator_string, fei::ParameterSet &paramset)
Definition: fei_utils.cpp:191
const char * version()
Definition: fei_utils.hpp:53
double cpu_time()
Definition: fei_utils.cpp:46
int get_filename_and_read_input(int argc, char **argv, MPI_Comm comm, int localProc, std::vector< std::string > &stdstrings)
double get_file_benchmark(const char *filename, const char *testname)
fei::SharedPtr< fei::Factory > create_fei_Factory(MPI_Comm comm, const char *libraryName)
std::ostream & console_out()
int numProcs(MPI_Comm comm)