33 #include "Thyra_AmesosLinearOpWithSolveFactory.hpp" 35 #include "Thyra_AmesosLinearOpWithSolve.hpp" 36 #include "Thyra_EpetraOperatorViewExtractorStd.hpp" 38 #include "Teuchos_dyn_cast.hpp" 39 #include "Teuchos_TimeMonitor.hpp" 41 #ifdef HAVE_AMESOS_KLU 44 #ifdef HAVE_AMESOS_LAPACK 47 #ifdef HAVE_AMESOS_MUMPS 50 #ifdef HAVE_AMESOS_SCALAPACK 53 #ifdef HAVE_AMESOS_UMFPACK 56 #ifdef HAVE_AMESOS_SUPERLUDIST 59 #ifdef HAVE_AMESOS_SUPERLU 62 #ifdef HAVE_AMESOS_DSCPACK 65 #if defined(HAVE_AMESOS_PARDISO) || defined(HAVE_AMESOS_PARDISO_MKL) 68 #ifdef HAVE_AMESOS_TAUCS 71 #ifdef HAVE_AMESOS_PARAKLETE 77 Teuchos::RCP<Teuchos::Time> overallTimer, constructTimer, symbolicTimer, factorTimer;
79 const std::string epetraFwdOp_str =
"epetraFwdOp";
88 const std::string AmesosLinearOpWithSolveFactory::SolverType_name =
"Solver Type";
90 const std::string AmesosLinearOpWithSolveFactory::RefactorizationPolicy_name =
"Refactorization Policy";
92 const std::string AmesosLinearOpWithSolveFactory::ThrowOnPreconditionerInput_name =
"Throw on Preconditioner Input";
94 const std::string AmesosLinearOpWithSolveFactory::Amesos_Settings_name =
"Amesos Settings";
98 AmesosLinearOpWithSolveFactory::~AmesosLinearOpWithSolveFactory()
102 paramList_->validateParameters(
103 *this->getValidParameters(),0
108 AmesosLinearOpWithSolveFactory::AmesosLinearOpWithSolveFactory(
109 const Amesos::ESolverType solverType
110 ,
const Amesos::ERefactorizationPolicy refactorizationPolicy
111 ,
const bool throwOnPrecInput
113 :epetraFwdOpViewExtractor_(Teuchos::rcp(
new EpetraOperatorViewExtractorStd()))
114 ,solverType_(solverType)
115 ,refactorizationPolicy_(refactorizationPolicy)
116 ,throwOnPrecInput_(throwOnPrecInput)
123 bool AmesosLinearOpWithSolveFactory::isCompatible(
124 const LinearOpSourceBase<double> &fwdOpSrc
127 Teuchos::RCP<const LinearOpBase<double> >
128 fwdOp = fwdOpSrc.getOp();
129 Teuchos::RCP<const Epetra_Operator> epetraFwdOp;
130 ETransp epetraFwdOpTransp;
131 EApplyEpetraOpAs epetraFwdOpApplyAs;
132 EAdjointEpetraOp epetraFwdOpAdjointSupport;
133 double epetraFwdOpScalar;
134 epetraFwdOpViewExtractor_->getEpetraOpView(
136 ,&epetraFwdOp,&epetraFwdOpTransp,&epetraFwdOpApplyAs,&epetraFwdOpAdjointSupport,&epetraFwdOpScalar
138 if( !dynamic_cast<const Epetra_RowMatrix*>(&*epetraFwdOp) )
143 Teuchos::RCP<LinearOpWithSolveBase<double> >
144 AmesosLinearOpWithSolveFactory::createOp()
const 146 return Teuchos::rcp(
new AmesosLinearOpWithSolve());
149 void AmesosLinearOpWithSolveFactory::initializeOp(
150 const Teuchos::RCP<
const LinearOpSourceBase<double> > &fwdOpSrc
151 ,LinearOpWithSolveBase<double> *Op
152 ,
const ESupportSolveUse supportSolveUse
155 Teuchos::TimeMonitor overallTimeMonitor(*overallTimer);
157 TEUCHOS_TEST_FOR_EXCEPT(Op==NULL);
159 const Teuchos::RCP<const LinearOpBase<double> >
160 fwdOp = fwdOpSrc->getOp();
164 Teuchos::RCP<const Epetra_Operator> epetraFwdOp;
165 ETransp epetraFwdOpTransp;
166 EApplyEpetraOpAs epetraFwdOpApplyAs;
167 EAdjointEpetraOp epetraFwdOpAdjointSupport;
168 double epetraFwdOpScalar;
169 epetraFwdOpViewExtractor_->getEpetraOpView(
170 fwdOp,&epetraFwdOp,&epetraFwdOpTransp,&epetraFwdOpApplyAs,&epetraFwdOpAdjointSupport,&epetraFwdOpScalar
173 AmesosLinearOpWithSolve
174 *amesosOp = &Teuchos::dyn_cast<AmesosLinearOpWithSolve>(*Op);
178 bool startOver = ( amesosOp->get_amesosSolver()==Teuchos::null );
182 epetraFwdOpTransp != amesosOp->get_amesosSolverTransp() ||
183 epetraFwdOp.get() != amesosOp->get_epetraLP()->GetOperator()
199 Teuchos::RCP<Epetra_LinearProblem>
201 epetraLP->SetOperator(const_cast<Epetra_Operator*>(&*epetraFwdOp));
202 Teuchos::set_extra_data< Teuchos::RCP<const Epetra_Operator> >( epetraFwdOp,
203 epetraFwdOp_str, Teuchos::inOutArg(epetraLP) );
205 Teuchos::RCP<Amesos_BaseSolver>
208 Teuchos::TimeMonitor constructTimeMonitor(*constructTimer);
209 switch(solverType_) {
213 #ifdef HAVE_AMESOS_KLU 215 amesosSolver = Teuchos::rcp(
new Amesos_Klu(*epetraLP));
218 #ifdef HAVE_AMESOS_MUMPS 220 amesosSolver = Teuchos::rcp(
new Amesos_Mumps(*epetraLP));
223 #ifdef HAVE_AMESOS_SCALAPACK 228 #ifdef HAVE_AMESOS_UMFPACK 233 #ifdef HAVE_AMESOS_SUPERLUDIST 238 #ifdef HAVE_AMESOS_SUPERLU 243 #ifdef HAVE_AMESOS_DSCPACK 248 #ifdef HAVE_AMESOS_PARDISO 253 #ifdef HAVE_AMESOS_TAUCS 255 amesosSolver = Teuchos::rcp(
new Amesos_Taucs(*epetraLP));
258 #ifdef HAVE_AMESOS_PARAKLETE 264 TEUCHOS_TEST_FOR_EXCEPTION(
265 true, std::logic_error
266 ,
"Error, the solver type ID = " << solverType_ <<
" is invalid!" 271 if(paramList_.get()) amesosSolver->setParameterList(sublist(paramList_,
"Amesos Settings"));
274 Teuchos::TimeMonitor symbolicTimeMonitor(*symbolicTimer);
275 amesosSolver->SymbolicFactorization();
278 Teuchos::TimeMonitor factorTimeMonitor(*factorTimer);
279 amesosSolver->NumericFactorization();
282 amesosOp->initialize(fwdOp,fwdOpSrc,epetraLP,amesosSolver,epetraFwdOpTransp,epetraFwdOpScalar);
291 Teuchos::RCP<Epetra_LinearProblem>
293 Teuchos::RCP<Amesos_BaseSolver>
294 amesosSolver = amesosOp->get_amesosSolver();
296 epetraLP->
SetOperator(const_cast<Epetra_Operator*>(&*epetraFwdOp));
297 Teuchos::get_extra_data< Teuchos::RCP<const Epetra_Operator> >(epetraLP,epetraFwdOp_str) = epetraFwdOp;
299 if(paramList_.get()) amesosSolver->setParameterList(sublist(paramList_,Amesos_Settings_name));
301 if(refactorizationPolicy_==Amesos::REPIVOT_ON_REFACTORIZATION) {
302 Teuchos::TimeMonitor symbolicTimeMonitor(*symbolicTimer);
303 amesosSolver->SymbolicFactorization();
306 Teuchos::TimeMonitor factorTimeMonitor(*factorTimer);
307 amesosSolver->NumericFactorization();
311 amesosOp->initialize(fwdOp,fwdOpSrc,epetraLP,amesosSolver,epetraFwdOpTransp,epetraFwdOpScalar);
313 amesosOp->setOStream(this->getOStream());
314 amesosOp->setVerbLevel(this->getVerbLevel());
317 bool AmesosLinearOpWithSolveFactory::supportsPreconditionerInputType(
const EPreconditionerInputType precOpType)
const 322 void AmesosLinearOpWithSolveFactory::initializePreconditionedOp(
323 const Teuchos::RCP<
const LinearOpSourceBase<double> > &fwdOpSrc
324 ,
const Teuchos::RCP<
const PreconditionerBase<double> > &prec
325 ,LinearOpWithSolveBase<double> *Op
326 ,
const ESupportSolveUse supportSolveUse
329 TEUCHOS_TEST_FOR_EXCEPTION(
330 this->throwOnPrecInput_, std::logic_error
331 ,
"Error, the concrete implementation described as \'"<<this->description()<<
"\' does not support preconditioners " 332 "and has been configured to throw this exception when the initializePreconditionedOp(...) function is called!" 334 this->initializeOp(fwdOpSrc,Op,supportSolveUse);
337 void AmesosLinearOpWithSolveFactory::initializePreconditionedOp(
338 const Teuchos::RCP<
const LinearOpSourceBase<double> > &fwdOpSrc
339 ,
const Teuchos::RCP<
const LinearOpSourceBase<double> > &approxFwdOpSrc
340 ,LinearOpWithSolveBase<double> *Op
341 ,
const ESupportSolveUse supportSolveUse
344 TEUCHOS_TEST_FOR_EXCEPTION(
345 this->throwOnPrecInput_, std::logic_error
346 ,
"Error, the concrete implementation described as \'"<<this->description()<<
"\' does not support preconditioners " 347 "and has been configured to throw this exception when the initializePreconditionedOp(...) function is called!" 349 this->initializeOp(fwdOpSrc,Op,supportSolveUse);
352 void AmesosLinearOpWithSolveFactory::uninitializeOp(
353 LinearOpWithSolveBase<double> *Op
354 ,Teuchos::RCP<
const LinearOpSourceBase<double> > *fwdOpSrc
355 ,Teuchos::RCP<
const PreconditionerBase<double> > *prec
356 ,Teuchos::RCP<
const LinearOpSourceBase<double> > *approxFwdOpSrc
357 ,ESupportSolveUse *supportSolveUse
361 TEUCHOS_TEST_FOR_EXCEPT(Op==NULL);
363 AmesosLinearOpWithSolve
364 *amesosOp = &Teuchos::dyn_cast<AmesosLinearOpWithSolve>(*Op);
365 Teuchos::RCP<const LinearOpSourceBase<double> >
366 _fwdOpSrc = amesosOp->extract_fwdOpSrc();
367 if(_fwdOpSrc.get()) {
369 Teuchos::RCP<Epetra_LinearProblem> epetraLP = amesosOp->get_epetraLP();
370 Teuchos::get_extra_data< Teuchos::RCP<const Epetra_Operator> >(
371 epetraLP,epetraFwdOp_str
379 if(fwdOpSrc) *fwdOpSrc = _fwdOpSrc;
380 if(prec) *prec = Teuchos::null;
381 if(approxFwdOpSrc) *approxFwdOpSrc = Teuchos::null;
386 void AmesosLinearOpWithSolveFactory::setParameterList(
387 Teuchos::RCP<Teuchos::ParameterList>
const& paramList
390 TEUCHOS_TEST_FOR_EXCEPT(paramList.get()==NULL);
391 paramList->validateParameters(*this->getValidParameters(),0);
392 paramList_ = paramList;
394 Amesos::solverTypeNameToEnumMap.get<Amesos::ESolverType>(
397 ,Amesos::toString(solverType_)
399 ,paramList_->name()+
"->"+SolverType_name
401 refactorizationPolicy_ =
402 Amesos::refactorizationPolicyNameToEnumMap.get<Amesos::ERefactorizationPolicy>(
404 RefactorizationPolicy_name
405 ,Amesos::toString(refactorizationPolicy_)
407 ,paramList_->name()+
"->"+RefactorizationPolicy_name
409 throwOnPrecInput_ = paramList_->get(ThrowOnPreconditionerInput_name,throwOnPrecInput_);
412 Teuchos::RCP<Teuchos::ParameterList>
413 AmesosLinearOpWithSolveFactory::getNonconstParameterList()
418 Teuchos::RCP<Teuchos::ParameterList>
419 AmesosLinearOpWithSolveFactory::unsetParameterList()
421 Teuchos::RCP<Teuchos::ParameterList> _paramList = paramList_;
422 paramList_ = Teuchos::null;
426 Teuchos::RCP<const Teuchos::ParameterList>
427 AmesosLinearOpWithSolveFactory::getParameterList()
const 432 Teuchos::RCP<const Teuchos::ParameterList>
433 AmesosLinearOpWithSolveFactory::getValidParameters()
const 435 return generateAndGetValidParameters();
440 std::string AmesosLinearOpWithSolveFactory::description()
const 442 std::ostringstream oss;
443 oss <<
"Thyra::AmesosLinearOpWithSolveFactory{";
444 oss <<
"solverType=" <<
toString(solverType_);
452 void AmesosLinearOpWithSolveFactory::initializeTimers()
454 if(!overallTimer.get()) {
455 overallTimer = Teuchos::TimeMonitor::getNewTimer(
"AmesosLOWSF");
456 constructTimer = Teuchos::TimeMonitor::getNewTimer(
"AmesosLOWSF:InitConstruct");
457 symbolicTimer = Teuchos::TimeMonitor::getNewTimer(
"AmesosLOWSF:Symbolic");
458 factorTimer = Teuchos::TimeMonitor::getNewTimer(
"AmesosLOWSF:Factor");
462 Teuchos::RCP<const Teuchos::ParameterList>
463 AmesosLinearOpWithSolveFactory::generateAndGetValidParameters()
465 static Teuchos::RCP<Teuchos::ParameterList> validParamList;
466 if(validParamList.get()==NULL) {
467 validParamList = Teuchos::rcp(
new Teuchos::ParameterList(
"Amesos"));
470 #ifdef HAVE_AMESOS_KLU
476 validParamList->set(RefactorizationPolicy_name,Amesos::toString(Amesos::REPIVOT_ON_REFACTORIZATION));
477 validParamList->set(ThrowOnPreconditionerInput_name,
bool(
true));
480 return validParamList;
Amesos_Klu: A serial, unblocked code ideal for getting started and for very sparse matrices...
void SetOperator(Epetra_RowMatrix *A)
Amesos_Paraklete: A serial, unblocked code ideal for getting started and for very sparse matrices...
Amesos_Superludist: An object-oriented wrapper for Superludist.
Amesos_Mumps: An object-oriented wrapper for the double precision version of MUMPS.
Amesos_Dscpack: An object-oriented wrapper for Dscpack.
static Teuchos::ParameterList GetValidParameters()
Get the list of valid parameters.
Amesos_Pardiso: Interface to the PARDISO package.
Amesos_Superlu: Amesos interface to Xioye Li's SuperLU 3.0 serial code.
std::string toString(ModelEvaluator::EDerivativeMultiVectorOrientation orientation)
Amesos_Lapack: an interface to LAPACK.
Class Amesos_Umfpack: An object-oriented wrapper for UMFPACK.
Amesos_Scalapack: A serial and parallel dense solver. For now, we implement only the unsymmetric ScaL...
Amesos_Taucs: An interface to the TAUCS package.