FEI Version of the Day
Loading...
Searching...
No Matches
beam_oldfei_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 the new FEI (version 3?) classes,
11// for the purposes of testing, code tuning and scaling studies.
12//
13
14#include <fei_macros.hpp> //simply includes std stuff like iostream, etc.
15 //users could include that stuff themselves rather
16 //than using this header.
17
18//Including the header fei_base.hpp gets us the declaration for
19//various fei classes and interfaces.
20
21#include <fei_base.hpp>
22#include <FEI_Implementation.hpp>
23
24//Now make provision for using any one of several solver libraries. This is
25//handled by the code in LibraryFactory.[hC].
26
27#include <test_utils/LibraryFactory.hpp>
28
29//And finally, we need to include some 'utils' headers which are used for
30//setting up the data for the example problem.
31
32#include <test_utils/HexBeam.hpp>
33#include <test_utils/HexBeamCR.hpp>
34
35#include <test_utils/fei_test_utils.hpp>
36#include <snl_fei_Utils.hpp>
37
38#undef fei_file
39#define fei_file "cube_main.cpp"
40#include <fei_ErrMacros.hpp>
41
42//==============================================================================
43//Here's the main...
44//==============================================================================
45int beam_oldfei_main(int argc, char** argv,
46 MPI_Comm comm, int numProcs, int localProc){
47
48 std::vector<std::string> stdstrings;
50 comm, localProc,
51 stdstrings) );
52
53 const char** params = NULL;
54 int numParams = 0;
55 fei::utils::strings_to_char_ptrs(stdstrings, numParams, params);
56
57 fei::ParameterSet paramset;
58 fei::utils::parse_strings(stdstrings, " ", paramset);
59
60 std::string whichFEI;
61 std::string solverName;
62 std::string datasource;
63 int W = 0;
64 int D = 0;
65 int DofPerNode = 0;
66
67 int errcode = 0;
68 errcode += paramset.getStringParamValue("SOLVER_LIBRARY", solverName);
69 errcode += paramset.getStringParamValue("DATA_SOURCE", datasource);
70 errcode += paramset.getStringParamValue("WHICH_FEI", whichFEI);
71 errcode += paramset.getIntParamValue("W", W);
72 errcode += paramset.getIntParamValue("D", D);
73 errcode += paramset.getIntParamValue("DofPerNode", DofPerNode);
74
75 if (errcode != 0) {
76 fei::console_out() << "Failed to find one or more required parameters in input-file."
77 << FEI_ENDL << "Required parameters:"<<FEI_ENDL
78 << "SOLVER_LIBRARY" << FEI_ENDL
79 << "DATA_SOURCE" << FEI_ENDL
80 << "CONSTRAINT_FORM" << FEI_ENDL
81 << "W" << FEI_ENDL << "D" << FEI_ENDL << "DofPerNode" << FEI_ENDL;
82 return(-1);
83 }
84
85 HexBeam* hexcubeptr = NULL;
86 if (datasource == "HexBeam") {
87 hexcubeptr = new HexBeam(W, D, DofPerNode,
88 HexBeam::OneD, numProcs, localProc);
89 }
90 else {
91 hexcubeptr = new HexBeamCR(W, D, DofPerNode,
92 HexBeam::OneD, numProcs, localProc);
93 }
94
95 HexBeam& hexcube = *hexcubeptr;
96
97 if (localProc == 0) {
98 int numCRs = (W+1)*(W+1)*(numProcs*2)-1;
99 if (hexcube.getNumCRs() < 1) numCRs = 0;
100 FEI_COUT << FEI_ENDL;
101 FEI_COUT << "========================================================"
102 << FEI_ENDL;
103 FEI_COUT << "Size W: " << W << " (num-elements-along-side-of-cube)"<<FEI_ENDL;
104 FEI_COUT << "Size D: " << D << " (num-elements-along-depth-of-cube)"<<FEI_ENDL;
105 FEI_COUT << "DOF per node: " << DofPerNode <<FEI_ENDL;
106 FEI_COUT << "Num local elements: " << hexcube.localNumElems_ << FEI_ENDL;
107 FEI_COUT << "Num global elements: " << hexcube.totalNumElems_ << FEI_ENDL;
108 FEI_COUT << "Num local DOF: " << hexcube.numLocalDOF_ << FEI_ENDL;
109 FEI_COUT << "Num global DOF: " << hexcube.numGlobalDOF_ << FEI_ENDL;
110 FEI_COUT << "Num global CRs: " << numCRs << FEI_ENDL;
111 FEI_COUT << "========================================================"
112 << FEI_ENDL;
113 }
114
115 double start_init_time = fei::utils::cpu_time();
116
117 //CHK_ERR( print_cube_data(hexcube, numProcs, localProc) );
118
122
123 if (whichFEI == "OLDFEI") {
124 try {
125 wrapper = fei::create_LibraryWrapper(comm, solverName.c_str());
126 }
127 catch (std::runtime_error& exc) {
128 fei::console_out() << exc.what() << FEI_ENDL;
129 ERReturn(-1);
130 }
131 fei.reset(new FEI_Implementation(wrapper, comm));
132 }
133 else if (whichFEI == "fei::FEI_Impl") {
134 try {
135 factory = fei::create_fei_Factory(comm, solverName.c_str());
136 }
137 catch (std::runtime_error& exc) {
138 fei::console_out() << exc.what() << FEI_ENDL;
139 ERReturn(-1);
140 }
141 fei = factory->createFEI(comm);
142 }
143 else {
144 fei::console_out() << "cube ERROR, value of 'WHICH_FEI' must be 'OLDFEI' or 'fei::FEI_Impl'"<< FEI_ENDL;
145 ERReturn(-1);
146 }
147
148 //load some control parameters.
149 CHK_ERR( fei->parameters(numParams, params) );
150
151 delete [] params;
152
153 CHK_ERR( fei->setSolveType(FEI_SINGLE_SYSTEM) );
154
155 int fieldID = 0;
156 int fieldSize = hexcube.numDofPerNode();
157
158 CHK_ERR( fei->initFields( 1, &fieldSize, &fieldID ) );
159
160 CHK_ERR( HexBeam_Functions::init_elem_connectivities( fei.get(), hexcube ) );
161
162 CHK_ERR( HexBeam_Functions::init_shared_nodes( fei.get(), hexcube ) );
163
164 int firstLocalCRID;
165 CHK_ERR( HexBeam_Functions::
166 init_constraints( fei.get(), hexcube, firstLocalCRID) );
167
168 CHK_ERR( fei->initComplete() );
169
170 double fei_init_time = fei::utils::cpu_time() - start_init_time;
171
172 if (localProc == 0) {
173 FEI_COUT.setf(IOS_FIXED, IOS_FLOATFIELD);
174 FEI_COUT << "Initialization time: " << fei_init_time << FEI_ENDL;
175 }
176
177 //Now the initialization phase is complete. Next we'll do the load phase,
178 //which for this problem just consists of loading the element data
179 //(element-wise stiffness arrays and load vectors) and the boundary
180 //condition data.
181
182 double start_load_time = fei::utils::cpu_time();
183
184 CHK_ERR( HexBeam_Functions::load_elem_data(fei.get(), hexcube) );
185
186 CHK_ERR( HexBeam_Functions::load_constraints(fei.get(), hexcube, firstLocalCRID) );
187
188 //FEI_COUT << "calling load_BC_data"<<FEI_ENDL;
189 //CHK_ERR( load_BC_data(fei, hexcube) );
190
191 fei->loadComplete();
192
193 double fei_load_time = fei::utils::cpu_time() - start_load_time;
194
195 delete hexcubeptr;
196
197 if (localProc == 0) {
198 //IOS macros are defined in fei_macros.h
199 FEI_COUT.setf(IOS_FIXED, IOS_FLOATFIELD);
200 FEI_COUT << "Coef. loading time: " << fei_load_time << FEI_ENDL;
201 FEI_COUT << "Total assembly time: " << fei_init_time + fei_load_time << FEI_ENDL;
202 }
203
204 //
205 //now the load phase is complete, so we're ready to launch the underlying
206 //solver and solve Ax=b
207 //
208
209 int err;
210 int status;
211 double start_solve_time = fei::utils::cpu_time();
212
213 err = fei->solve(status);
214
215 if (err || status) {
216 if (localProc==0) FEI_COUT << "solve returned status: " << status << FEI_ENDL;
217 }
218
219 double solve_time = fei::utils::cpu_time() - start_solve_time;
220
221 if (localProc == 0) {
222 FEI_COUT << "Solver time: " << solve_time << FEI_ENDL;
223 }
224
225 int returnValue = 0;
226 if (localProc == 0) {
227#if defined(FEI_PLATFORM) && defined(FEI_OPT_LEVEL)
228 double benchmark = fei_init_time;
229
230 FEI_OSTRINGSTREAM testname_init;
231 testname_init << "cube_"<<whichFEI<<"_init_"<<W<<"_"<<D<<"_"<<DofPerNode<<"_"
232 <<solverName<<"_np"<<numProcs<<"_"
233 <<FEI_PLATFORM<<"_"<<FEI_OPT_LEVEL;
234
235 FEI_OSTRINGSTREAM testname_load;
236 testname_load << "cube_"<<whichFEI<<"_load_"<<W<<"_"<<D<<"_"<<DofPerNode<<"_"
237 <<solverName<<"_np"<<numProcs<<"_"
238 <<FEI_PLATFORM<<"_"<<FEI_OPT_LEVEL;
239
240 double file_init, file_load;
241 bool file_benchmarks_available = true;
242 try {
243 file_init = fei_test_utils::get_file_benchmark("./cube_timings.txt",
244 testname_init.str().c_str());
245 file_load = fei_test_utils::get_file_benchmark("./cube_timings.txt",
246 testname_load.str().c_str());
247 }
248 catch (std::runtime_error& exc) {
249 file_benchmarks_available = false;
250 }
251
252 if (file_benchmarks_available) {
253
254 bool init_test_passed =
255 fei_test_utils::check_and_cout_test_result(testname_init.str(), benchmark,
256 file_init, 10);
257
258 benchmark = fei_load_time;
259 bool load_test_passed =
260 fei_test_utils::check_and_cout_test_result(testname_load.str(), benchmark,
261 file_load, 10);
262
263 returnValue = init_test_passed&&load_test_passed ? 0 : 1;
264 }
265#endif
266 }
267
268 bool testPassed = returnValue==0;
269 if (testPassed && localProc == 0) {
270 //A string for the SIERRA runtest tool to search for in test output...
271 FEI_COUT << "cube execution successful" << FEI_ENDL;
272 }
273
274 fei.reset();
275
276 return(0);
277}
int getIntParamValue(const char *name, int &paramValue) const
int getStringParamValue(const char *name, std::string &paramValue) const
void parse_strings(std::vector< std::string > &stdstrings, const char *separator_string, fei::ParameterSet &paramset)
Definition: fei_utils.cpp:191
double cpu_time()
Definition: fei_utils.cpp:46
void strings_to_char_ptrs(std::vector< std::string > &stdstrings, int &numStrings, const char **&charPtrs)
Definition: fei_utils.cpp:178
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()
fei::SharedPtr< LibraryWrapper > create_LibraryWrapper(MPI_Comm comm, const char *libraryName)
int numProcs(MPI_Comm comm)