FEI Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
fei_Factory_Trilinos.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 #include "fei_trilinos_macros.hpp"
45 
46 #include <fei_Factory_Trilinos.hpp>
47 
48 #include <fei_VectorReducer.hpp>
49 #include <fei_MatrixReducer.hpp>
50 #include <fei_Matrix_Local.hpp>
51 #include <fei_Vector_Local.hpp>
52 
53 #ifdef HAVE_FEI_AZTECOO
54 #include <fei_Solver_AztecOO.hpp>
55 #endif
56 #ifdef HAVE_FEI_BELOS
57 #include <fei_Solver_Belos.hpp>
58 #endif
59 #ifdef HAVE_FEI_AMESOS
60 #include <fei_Solver_Amesos.hpp>
61 #endif
62 
63 #ifdef HAVE_FEI_EPETRA
64 
66  : fei::Factory(comm),
67  comm_(comm),
68  reducer_(),
69  lpm_epetrabasic_(),
70  use_lpm_epetrabasic_(false),
71  useAmesos_(false),
72  useBelos_(false),
73  use_feiMatrixLocal_(false),
74  blockEntryMatrix_(false),
75  orderRowsWithLocalColsFirst_(false),
76  outputLevel_(0)
77 {
78 }
79 
81 {
82 }
83 
84 int Factory_Trilinos::parameters(int numParams,
85  const char* const* paramStrings)
86 {
87  std::vector<std::string> stdstrings;
88  fei::utils::char_ptrs_to_strings(numParams, paramStrings, stdstrings);
89 
90  fei::ParameterSet paramset;
91  fei::utils::parse_strings(stdstrings, " ", paramset);
92 
93  parameters(paramset);
94  return(0);
95 }
96 
97 void Factory_Trilinos::parameters(const fei::ParameterSet& parameterset)
98 {
99  fei::Factory::parameters(parameterset);
100 
101  parameterset.getIntParamValue("outputLevel", outputLevel_);
102 
103  bool blkGraph = false;
104  bool blkMatrix = false;
105 
106  parameterset.getBoolParamValue("BLOCK_GRAPH", blkGraph);
107  parameterset.getBoolParamValue("BLOCK_MATRIX", blkMatrix);
108 
109  blockEntryMatrix_ = (blkGraph || blkMatrix);
110 
111  const fei::Param* param = parameterset.get("Trilinos_Solver");
112  if (param != 0) {
113  if (param->getType() == fei::Param::STRING) {
114  std::string strval = param->getStringValue();
115  std::string::size_type ii = strval.find("Amesos");
116  if (ii != std::string::npos) {
117  useAmesos_ = true;
118  }
119  ii = strval.find("Belos");
120  if (ii != std::string::npos) {
121  useBelos_ = true;
122  }
123  }
124  }
125 
126  parameterset.getBoolParamValue("LPM_EpetraBasic", use_lpm_epetrabasic_);
127  if (use_lpm_epetrabasic_) {
129  }
130 
131  parameterset.getBoolParamValue("USE_FEI_MATRIX_LOCAL", use_feiMatrixLocal_);
132 
133  parameterset.getBoolParamValue("ORDER_ROWS_WITH_LOCAL_COLS_FIRST",
135 }
136 
140  const char* name)
141 {
142  static fei::MatrixGraph_Impl2::Factory factory2;
143  if (!use_lpm_epetrabasic_) {
144  return(factory2.createMatrixGraph(rowSpace, colSpace, name));
145  }
146 
147  return(factory2.createMatrixGraph(rowSpace, colSpace, name));
148 }
149 
151 Factory_Trilinos::wrapVector(fei::SharedPtr<fei::VectorSpace> vecSpace,
153 {
155  if (vecSpace->getNumIndices_Owned() != multiVec->Map().NumMyPoints()) {
156  //vecSpace isn't compatible with multiVec's map, so return empty vecptr
157  return(vecptr);
158  }
159 
160  vecptr.reset(new fei::Vector_Impl<Epetra_MultiVector>(vecSpace,
161  multiVec.get(),
162  multiVec->Map().NumMyPoints(), false));
163  return(vecptr);
164 }
165 
167 Factory_Trilinos::wrapVector(fei::SharedPtr<fei::MatrixGraph> matGraph,
169 {
171  if (matGraph->getRowSpace()->getNumIndices_Owned() !=
172  multiVec->Map().NumMyPoints()) {
173  //vector-space isn't compatible with multiVec's map, so return empty vecptr
174  return(vecptr);
175  }
176 
178  multiVec.get(),
179  multiVec->Map().NumMyPoints(), false));
180  return(vecptr);
181 }
182 
185  bool isSolutionVector,
186  int numVectors)
187 {
188  if (use_feiMatrixLocal_) {
189  fei::SharedPtr<fei::Vector> feivec(new fei::Vector_Local(vecSpace));
190  return(feivec);
191  }
192 
193 // if (blockEntryMatrix_) {
194 // throw std::runtime_error("Factory_Trilinos: fei ERROR, block-entry matrices/vectors not currently supported.");
195 // }
196 
197  std::vector<int> indices;
198  int err = 0, localSize;
199  if (reducer_.get() != NULL) {
200  indices = reducer_->getLocalReducedEqns();
201  localSize = indices.size();
202  }
203  else {
204  if (blockEntryMatrix_) {
205  localSize = vecSpace->getNumBlkIndices_Owned();
206  indices.resize(localSize*2);
207  err = vecSpace->getBlkIndices_Owned(localSize, &indices[0], &indices[localSize], localSize);
208  }
209  else {
210  localSize = vecSpace->getNumIndices_Owned();
211  indices.resize(localSize);
212  err = vecSpace->getIndices_Owned(indices);
213  }
214  }
215  if (err != 0) {
216  throw std::runtime_error("error in vecSpace->getIndices_Owned");
217  }
218 
219  fei::SharedPtr<fei::Vector> feivec, tmpvec;
220  if (!use_lpm_epetrabasic_) {
221  try {
223  Trilinos_Helpers::create_Epetra_BlockMap(vecSpace) :
224  Trilinos_Helpers::create_Epetra_Map(comm_, indices);
225 
226  Epetra_MultiVector* emvec = new Epetra_MultiVector(emap, numVectors);
227 
228  tmpvec.reset(new fei::Vector_Impl<Epetra_MultiVector>(vecSpace,
229  emvec, localSize,
230  isSolutionVector, true));
231  }
232  catch(std::runtime_error& exc) {
233  fei::console_out() << "Factory_Trilinos::createVector: caught exception '"
234  << exc.what() << "', re-throwing..." << FEI_ENDL;
235  throw exc;
236  }
237  }
238  else {
240  lpm_epetrabasic_.get(), localSize,
241  isSolutionVector));
242  }
243 
244  if (reducer_.get() != NULL) {
245  feivec.reset(new fei::VectorReducer(reducer_,
246  tmpvec, isSolutionVector));
247  }
248  else {
249  feivec = tmpvec;
250  }
251 
252  return(feivec);
253 }
254 
257  int numVectors)
258 {
259  bool isSolnVector = false;
260  return(createVector(vecSpace, isSolnVector, numVectors));
261 }
262 
265  int numVectors)
266 {
267  bool isSolnVector = false;
268  return(createVector(matrixGraph, isSolnVector, numVectors));
269 }
270 
273  bool isSolutionVector,
274  int numVectors)
275 {
276  if (use_feiMatrixLocal_) {
277  fei::SharedPtr<fei::Vector> feivec(new fei::Vector_Local(matrixGraph->getRowSpace()));
278  return(feivec);
279  }
280 
281  int globalNumSlaves = matrixGraph->getGlobalNumSlaveConstraints();
282 
283  if (globalNumSlaves > 0 && reducer_.get()==NULL) {
284  reducer_ = matrixGraph->getReducer();
285  }
286 
287  fei::SharedPtr<fei::Vector> feivec, tmpvec;
288 
289  std::vector<int> indices;
290  int err = 0, localSize;
291  fei::SharedPtr<fei::VectorSpace> vecSpace = matrixGraph->getRowSpace();
292  if (reducer_.get() != NULL) {
293  indices = reducer_->getLocalReducedEqns();
294  localSize = indices.size();
295  }
296  else {
297  localSize = vecSpace->getNumIndices_Owned();
298  indices.resize(localSize);
299  err = vecSpace->getIndices_Owned(indices);
300  }
301  if (err != 0) {
302  throw std::runtime_error("error in vecSpace->getIndices_Owned");
303  }
304 
305  if (!use_lpm_epetrabasic_) {
306 #ifdef HAVE_FEI_EPETRA
307  try {
309  Trilinos_Helpers::create_Epetra_BlockMap(vecSpace) :
310  Trilinos_Helpers::create_Epetra_Map(comm_, indices);
311 
312  Epetra_MultiVector* emvec = new Epetra_MultiVector(emap, numVectors);
313 
314  tmpvec.reset(new fei::Vector_Impl<Epetra_MultiVector>(matrixGraph->getRowSpace(), emvec,
315  localSize, isSolutionVector, true));
316  }
317  catch(std::runtime_error& exc) {
318  fei::console_out() << "Factory_Trilinos::createVector: caught exception '"
319  << exc.what() << "', re-throwing..." << FEI_ENDL;
320  throw exc;
321  }
322 #else
323  fei::console_out() << "fei_Factory_Trilinos::createVector ERROR, HAVE_FEI_EPETRA not defined."
324  << FEI_ENDL;
325 #endif
326  }
327  else {
328 #ifdef HAVE_FEI_EPETRA
329  vecSpace = matrixGraph->getRowSpace();
330 
333  lpm_epetrabasic_.get(), localSize, isSolutionVector));
334 #endif
335  }
336 
337  if (reducer_.get() != NULL) {
338  feivec.reset(new fei::VectorReducer(reducer_, tmpvec, isSolutionVector));
339  }
340  else {
341  feivec = tmpvec;
342  }
343 
344  return(feivec);
345 }
346 
349 {
351 #ifdef HAVE_FEI_EPETRA
352  int globalNumSlaves = matrixGraph->getGlobalNumSlaveConstraints();
353 
354  if (globalNumSlaves > 0 && reducer_.get()==NULL) {
355  reducer_ = matrixGraph->getReducer();
356  }
357 
358  if (use_lpm_epetrabasic_) {
359  return(
360  Trilinos_Helpers::create_from_LPM_EpetraBasic(matrixGraph,
362  reducer_,
364  );
365  }
366 
367  if (use_feiMatrixLocal_) {
369  }
370 
371  return(
372  Trilinos_Helpers::create_from_Epetra_Matrix(matrixGraph, blockEntryMatrix_,
374  );
375 #else
376  fei::console_out() << "fei_Factory_Trilinos::createMatrix ERROR, HAVE_FEI_EPETRA "
377  << "not defined."<<FEI_ENDL;
378  return feimat;
379 #endif
380 }
381 
383 Factory_Trilinos::createSolver(const char* name)
384 {
386 
387  std::string::size_type ii_amesos = std::string::npos;
388  std::string::size_type ii_belos = std::string::npos;
389 
390  if (name != 0) {
391  std::string sname(name);
392  ii_amesos = sname.find("Amesos");
393  ii_belos = sname.find("Belos");
394  }
395 
396  if (useAmesos_ || ii_amesos != std::string::npos) {
397 #ifdef HAVE_FEI_AMESOS
398  solver.reset(new Solver_Amesos);
399  return solver;
400 #else
401  fei::console_out() << "fei_Factory_Trilinos::createSolver: ERROR, Amesos not available (not enabled at compile-time?)"<<FEI_ENDL;
402 #endif
403  }
404 
405  if (useBelos_ || ii_belos != std::string::npos) {
406 #ifdef HAVE_FEI_BELOS
407  solver.reset(new Solver_Belos);
408  return solver;
409 #else
410  fei::console_out() << "fei_Factory_Trilinos::createSolver: ERROR, Belos not available (not enabled at compile-time?)"<<FEI_ENDL;
411 #endif
412  }
413 
414 #ifdef HAVE_FEI_AZTECOO
415  solver.reset(new Solver_AztecOO);
416  return solver;
417 #else
418  fei::console_out() << "fei_Factory_Trilinos::createSolver: ERROR, AztecOO not "
419  << "available." << FEI_ENDL;
420  return solver;
421 #endif
422 }
423 
424 void Factory_Trilinos::create_LinProbMgr(bool replace_if_already_created)
425 {
426  if (!use_lpm_epetrabasic_) {
427  return;
428  }
429 
430  bool need_to_create_lpm = false;
431 
432  if (lpm_epetrabasic_.get() == NULL) {
433  need_to_create_lpm = true;
434  }
435 
436  if (replace_if_already_created) {
437  need_to_create_lpm = true;
438  }
439 
440  if (need_to_create_lpm) {
441 #ifdef HAVE_FEI_EPETRA
443  newlpm(new LinProbMgr_EpetraBasic(comm_));
444 
445  lpm_epetrabasic_ = newlpm;
446 #else
447  fei::console_out() << "fei_Factory_Trilinos::create_LinProbMgr ERROR, HAVE_FEI_EPETRA"
448  <<" not defined."<<FEI_ENDL;
449 #endif
450  }
451 }
452 
453 //HAVE_FEI_EPETRA
454 #endif
455 
void char_ptrs_to_strings(int numStrings, const char *const *charstrings, std::vector< std::string > &stdstrings)
Definition: fei_utils.cpp:164
int getBlkIndices_Owned(int lenBlkIndices, int *globalBlkIndices, int *blkSizes, int &numBlkIndices)
ParamType getType() const
Definition: fei_Param.hpp:98
Factory_Trilinos(MPI_Comm comm)
virtual int getGlobalNumSlaveConstraints() const =0
const Param * get(const char *name) const
virtual void parameters(const fei::ParameterSet &paramset)
Definition: fei_Factory.cpp:38
virtual ~Factory_Trilinos()
virtual fei::SharedPtr< fei::Reducer > getReducer()=0
int getBoolParamValue(const char *name, bool &paramValue) const
fei::SharedPtr< fei::Matrix > createMatrix(fei::SharedPtr< fei::MatrixGraph > matrixGraph)
fei::SharedPtr< fei::Solver > createSolver(const char *name=0)
void reset(T *p=0)
#define MPI_Comm
Definition: fei_mpi.h:56
virtual fei::SharedPtr< fei::VectorSpace > getRowSpace()=0
virtual int parameters(int numParams, const char *const *paramStrings)
virtual const Epetra_BlockMap & Map() const =0
virtual fei::SharedPtr< fei::MatrixGraph > createMatrixGraph(fei::SharedPtr< fei::VectorSpace > rowSpace, fei::SharedPtr< fei::VectorSpace > columnSpace, const char *name)
void create_LinProbMgr(bool replace_if_already_created=false)
fei::SharedPtr< fei::Reducer > reducer_
fei::SharedPtr< fei::Vector > createVector(fei::SharedPtr< fei::VectorSpace > vecSpace, int numVectors=1)
T * get() const
int getIndices_Owned(std::vector< int > &globalIndices) const
#define FEI_ENDL
std::ostream & console_out()
void parse_strings(std::vector< std::string > &stdstrings, const char *separator_string, fei::ParameterSet &paramset)
Definition: fei_utils.cpp:191
fei::SharedPtr< fei::LinearProblemManager > lpm_epetrabasic_
const std::string & getStringValue() const
Definition: fei_Param.hpp:104
int getNumBlkIndices_Owned() const
std::vector< int > & getLocalReducedEqns()
virtual void setRowDistribution(const std::vector< int > &ownedGlobalRows)=0
int NumMyPoints() const
static fei::SharedPtr< fei::Matrix > create_Matrix_Local(fei::SharedPtr< fei::MatrixGraph > matrixGraph, bool blockEntry)
fei::SharedPtr< fei::MatrixGraph > createMatrixGraph(fei::SharedPtr< fei::VectorSpace > rowSpace, fei::SharedPtr< fei::VectorSpace > colSpace, const char *name)
int getNumIndices_Owned() const
int getIntParamValue(const char *name, int &paramValue) const