FEI Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
beam_oldfei.cpp
Go to the documentation of this file.
1/*
2// @HEADER
3// ************************************************************************
4// FEI: Finite Element Interface to Linear Solvers
5// Copyright (2005) Sandia Corporation.
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, the
8// U.S. Government retains certain rights in this software.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Alan Williams (william@sandia.gov)
38//
39// ************************************************************************
40// @HEADER
41*/
42
43
44//
45// This is a simple program to exercise the FEI interface
46// for the purposes of testing, code tuning and scaling studies.
47//
48
49#include <fei_iostream.hpp>
50
51//Including the header fei_base.hpp pulls in the declaration for
52//various fei classes and interfaces.
53
54#include <fei_base.hpp>
56
57//Now make provision for using any one of several solver libraries. This is
58//handled by the code in LibraryFactory.{hpp,cpp}.
59
61
62
63//we need to include some 'utils' headers which are used for
64//setting up the data for the example problem.
65
68
70#include <snl_fei_Utils.hpp>
71
72#undef fei_file
73#define fei_file "beam_oldfei.cpp"
74#include <fei_ErrMacros.hpp>
75
76//==============================================================================
77//Here's the main...
78//==============================================================================
79int main(int argc, char** argv)
80{
81 MPI_Comm comm;
82 int numProcs, localProc;
83 CHK_ERR( fei_test_utils::initialize_mpi(argc, argv, localProc, numProcs) );
84 comm = MPI_COMM_WORLD;
85
86 std::vector<std::string> stdstrings;
88 comm, localProc,
89 stdstrings) );
90
91 const char** params = NULL;
92 int numParams = 0;
93 fei::utils::strings_to_char_ptrs(stdstrings, numParams, params);
94
95 fei::ParameterSet paramset;
96 fei::utils::parse_strings(stdstrings, " ", paramset);
97
98 std::string whichFEI;
99 std::string solverName;
100 std::string datasource;
101 int W = 0;
102 int D = 0;
103 int DofPerNode = 0;
104
105 int errcode = 0;
106 errcode += paramset.getStringParamValue("SOLVER_LIBRARY", solverName);
107 errcode += paramset.getStringParamValue("DATA_SOURCE", datasource);
108 errcode += paramset.getStringParamValue("WHICH_FEI", whichFEI);
109 errcode += paramset.getIntParamValue("W", W);
110 errcode += paramset.getIntParamValue("D", D);
111 errcode += paramset.getIntParamValue("DofPerNode", DofPerNode);
112
113 if (errcode != 0) {
114 fei::console_out() << "Failed to find one or more required parameters in input-file."
115 << std::endl << "Required parameters:"<<std::endl
116 << "SOLVER_LIBRARY" << std::endl
117 << "DATA_SOURCE" << std::endl
118 << "WHICH_FEI" << std::endl
119 << "W" << std::endl << "D" << std::endl << "DofPerNode" << std::endl;
120 return(-1);
121 }
122
123 HexBeam* hexcubeptr = NULL;
124 if (datasource == "HexBeam") {
125 hexcubeptr = new HexBeam(W, D, DofPerNode,
126 HexBeam::OneD, numProcs, localProc);
127 }
128 else {
129 hexcubeptr = new HexBeamCR(W, D, DofPerNode,
130 HexBeam::OneD, numProcs, localProc);
131 }
132
133 HexBeam& hexcube = *hexcubeptr;
134
135 if (localProc == 0) {
136 int numCRs = (W+1)*(W+1)*(numProcs*2)-1;
137 if (hexcube.getNumCRs() < 1) numCRs = 0;
138 std::cout << std::endl;
139 std::cout << "========================================================"
140 << std::endl;
141 std::cout << "Size W: " << W << " (num-elements-along-side-of-cube)"<<std::endl;
142 std::cout << "Size D: " << D << " (num-elements-along-depth-of-cube)"<<std::endl;
143 std::cout << "DOF per node: " << DofPerNode <<std::endl;
144 std::cout << "Num local elements: " << hexcube.localNumElems_ << std::endl;
145 std::cout << "Num global elements: " << hexcube.totalNumElems_ << std::endl;
146 std::cout << "Num local DOF: " << hexcube.numLocalDOF_ << std::endl;
147 std::cout << "Num global DOF: " << hexcube.numGlobalDOF_ << std::endl;
148 std::cout << "Num global CRs: " << numCRs << std::endl;
149 std::cout << "========================================================"
150 << std::endl;
151 }
152
153 double start_init_time = fei::utils::cpu_time();
154
155 //CHK_ERR( print_cube_data(hexcube, numProcs, localProc) );
156
160
161 if (whichFEI == "OLDFEI") {
162 try {
163 wrapper = fei::create_LibraryWrapper(comm, solverName.c_str());
164 }
165 catch (std::runtime_error& exc) {
166 fei::console_out() << exc.what() << std::endl;
167 ERReturn(-1);
168 }
169 fei.reset(new FEI_Implementation(wrapper, comm));
170 }
171 else if (whichFEI == "fei::FEI_Impl") {
172 try {
173 factory = fei::create_fei_Factory(comm, solverName.c_str());
174 }
175 catch (std::runtime_error& exc) {
176 fei::console_out() << exc.what() << std::endl;
177 ERReturn(-1);
178 }
179 fei = factory->createFEI(comm);
180 }
181 else {
182 fei::console_out() << "beam ERROR, value of 'WHICH_FEI' must be 'OLDFEI' or 'fei::FEI_Impl'"<< std::endl;
183 ERReturn(-1);
184 }
185
186 //load some control parameters.
187 CHK_ERR( fei->parameters(numParams, params) );
188
189 delete [] params;
190
191 CHK_ERR( fei->setSolveType(FEI_SINGLE_SYSTEM) );
192
193 int fieldID = 0;
194 int fieldSize = hexcube.numDofPerNode();
195
196 CHK_ERR( fei->initFields( 1, &fieldSize, &fieldID ) );
197
199
201
202 int firstLocalCRID;
204 init_constraints( fei.get(), hexcube, firstLocalCRID) );
205
206 CHK_ERR( fei->initComplete() );
207
208 double fei_init_time = fei::utils::cpu_time() - start_init_time;
209
210 if (localProc == 0) {
211 std::cout.setf(IOS_FIXED, IOS_FLOATFIELD);
212 std::cout << "Initialization time: " << fei_init_time << std::endl;
213 }
214
215 //Now the initialization phase is complete. Next we'll do the load phase,
216 //which for this problem just consists of loading the element data
217 //(element-wise stiffness arrays and load vectors) and the boundary
218 //condition data.
219
220 double start_load_time = fei::utils::cpu_time();
221
223
224 CHK_ERR( HexBeam_Functions::load_constraints(fei.get(), hexcube, firstLocalCRID) );
225
226 //std::cout << "calling load_BC_data"<<std::endl;
227 //CHK_ERR( load_BC_data(fei, hexcube) );
228
229 fei->loadComplete();
230
231 double fei_load_time = fei::utils::cpu_time() - start_load_time;
232
233 delete hexcubeptr;
234
235 if (localProc == 0) {
236 //IOS macros are defined in fei_macros.h
237 std::cout.setf(IOS_FIXED, IOS_FLOATFIELD);
238 std::cout << "Coef. loading time: " << fei_load_time << std::endl;
239 std::cout << "Total assembly time: " << fei_init_time + fei_load_time << std::endl;
240 }
241
242 //
243 //now the load phase is complete, so we're ready to launch the underlying
244 //solver and solve Ax=b
245 //
246
247 int err;
248 int status;
249 double start_solve_time = fei::utils::cpu_time();
250
251 err = fei->solve(status);
252
253 if (err || status) {
254 if (localProc==0) std::cout << "solve returned status: " << status << std::endl;
255 }
256
257 double solve_time = fei::utils::cpu_time() - start_solve_time;
258
259 if (localProc == 0) {
260 std::cout << "Solver time: " << solve_time << std::endl;
261 }
262
263 int returnValue = 0;
264 if (localProc == 0) {
265#if defined(FEI_PLATFORM) && defined(FEI_OPT_LEVEL)
266 double benchmark = fei_init_time;
267
268 FEI_OSTRINGSTREAM testname_init;
269 testname_init << "cube_"<<whichFEI<<"_init_"<<W<<"_"<<D<<"_"<<DofPerNode<<"_"
270 <<solverName<<"_np"<<numProcs<<"_"
271 <<FEI_PLATFORM<<"_"<<FEI_OPT_LEVEL;
272
273 FEI_OSTRINGSTREAM testname_load;
274 testname_load << "cube_"<<whichFEI<<"_load_"<<W<<"_"<<D<<"_"<<DofPerNode<<"_"
275 <<solverName<<"_np"<<numProcs<<"_"
276 <<FEI_PLATFORM<<"_"<<FEI_OPT_LEVEL;
277
278 double file_init, file_load;
279 bool file_benchmarks_available = true;
280 try {
281 file_init = fei_test_utils::get_file_benchmark("./cube_timings.txt",
282 testname_init.str().c_str());
283 file_load = fei_test_utils::get_file_benchmark("./cube_timings.txt",
284 testname_load.str().c_str());
285 }
286 catch (std::runtime_error& exc) {
287 file_benchmarks_available = false;
288 }
289
290 if (file_benchmarks_available) {
291
292 bool init_test_passed =
293 fei_test_utils::check_and_cout_test_result(testname_init.str(), benchmark,
294 file_init, 10);
295
296 benchmark = fei_load_time;
297 bool load_test_passed =
298 fei_test_utils::check_and_cout_test_result(testname_load.str(), benchmark,
299 file_load, 10);
300
301 returnValue = init_test_passed&&load_test_passed ? 0 : 1;
302 }
303#endif
304 }
305
306 bool testPassed = returnValue==0;
307 if (testPassed && localProc == 0) {
308 //A string for the SIERRA runtest tool to search for in test output...
309 std::cout << "beam execution successful" << std::endl;
310 }
311
312 fei.reset();
313
314#ifndef FEI_SER
315 MPI_Finalize();
316#endif
317
318 return(0);
319}
int main(int argc, char **argv)
Definition: beam_oldfei.cpp:79
@ OneD
Definition: HexBeam.hpp:26
virtual int numDofPerNode()
Definition: HexBeam.hpp:36
int numLocalDOF_
Definition: HexBeam.hpp:91
int totalNumElems_
Definition: HexBeam.hpp:75
virtual int getNumCRs()
Definition: HexBeam.hpp:63
int localNumElems_
Definition: HexBeam.hpp:77
int numGlobalDOF_
Definition: HexBeam.hpp:92
int getIntParamValue(const char *name, int &paramValue) const
int getStringParamValue(const char *name, std::string &paramValue) const
#define ERReturn(a)
#define CHK_ERR(a)
#define FEI_SINGLE_SYSTEM
Definition: fei_defs.h:65
#define IOS_FIXED
#define IOS_FLOATFIELD
#define MPI_COMM_WORLD
Definition: fei_mpi.h:58
#define MPI_Comm
Definition: fei_mpi.h:56
#define FEI_OSTRINGSTREAM
Definition: fei_sstream.hpp:32
#define FEI_Implementation
Definition: fei_version.h:66
int init_elem_connectivities(FEI *fei, HexBeam &hexcube)
Definition: HexBeam.cpp:280
int load_elem_data(FEI *fei, HexBeam &hexcube)
Definition: HexBeam.cpp:440
int load_constraints(FEI *fei, HexBeam &hexcube, int firstLocalCRID)
Definition: HexBeam.cpp:397
int init_shared_nodes(FEI *fei, HexBeam &hexcube)
Definition: HexBeam.cpp:328
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 initialize_mpi(int argc, char **argv, int &localProc, int &numProcs)
int get_filename_and_read_input(int argc, char **argv, MPI_Comm comm, int localProc, std::vector< std::string > &stdstrings)
bool check_and_cout_test_result(std::string testname, double value, double file_value, unsigned margin)
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)