88const std::string AmesosLinearOpWithSolveFactory::SolverType_name =
"Solver Type";
90const std::string AmesosLinearOpWithSolveFactory::RefactorizationPolicy_name =
"Refactorization Policy";
92const std::string AmesosLinearOpWithSolveFactory::ThrowOnPreconditionerInput_name =
"Throw on Preconditioner Input";
94const std::string AmesosLinearOpWithSolveFactory::Amesos_Settings_name =
"Amesos Settings";
98AmesosLinearOpWithSolveFactory::~AmesosLinearOpWithSolveFactory()
102 paramList_->validateParameters(
103 *this->getValidParameters(),0
108AmesosLinearOpWithSolveFactory::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)
123bool 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) )
143Teuchos::RCP<LinearOpWithSolveBase<double> >
144AmesosLinearOpWithSolveFactory::createOp()
const
146 return Teuchos::rcp(
new AmesosLinearOpWithSolve());
149void 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>
200 epetraLP = Teuchos::rcp(
new 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_) {
210 case Thyra::Amesos::LAPACK :
213#ifdef HAVE_AMESOS_KLU
214 case Thyra::Amesos::KLU :
215 amesosSolver = Teuchos::rcp(
new Amesos_Klu(*epetraLP));
218#ifdef HAVE_AMESOS_MUMPS
219 case Thyra::Amesos::MUMPS :
220 amesosSolver = Teuchos::rcp(
new Amesos_Mumps(*epetraLP));
223#ifdef HAVE_AMESOS_SCALAPACK
224 case Thyra::Amesos::SCALAPACK :
228#ifdef HAVE_AMESOS_UMFPACK
229 case Thyra::Amesos::UMFPACK :
233#ifdef HAVE_AMESOS_SUPERLUDIST
234 case Thyra::Amesos::SUPERLUDIST :
238#ifdef HAVE_AMESOS_SUPERLU
239 case Thyra::Amesos::SUPERLU :
243#ifdef HAVE_AMESOS_DSCPACK
244 case Thyra::Amesos::DSCPACK :
248#ifdef HAVE_AMESOS_PARDISO
249 case Thyra::Amesos::PARDISO :
253#ifdef HAVE_AMESOS_TAUCS
254 case Thyra::Amesos::TAUCS :
255 amesosSolver = Teuchos::rcp(
new Amesos_Taucs(*epetraLP));
258#ifdef HAVE_AMESOS_PARAKLETE
259 case Thyra::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>
292 epetraLP = Teuchos::rcp_const_cast<Epetra_LinearProblem>(amesosOp->get_epetraLP());
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());
317bool AmesosLinearOpWithSolveFactory::supportsPreconditionerInputType(
const EPreconditionerInputType precOpType)
const
322void 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);
337void 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);
352void 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;
386void 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_);
412Teuchos::RCP<Teuchos::ParameterList>
413AmesosLinearOpWithSolveFactory::getNonconstParameterList()
418Teuchos::RCP<Teuchos::ParameterList>
419AmesosLinearOpWithSolveFactory::unsetParameterList()
421 Teuchos::RCP<Teuchos::ParameterList> _paramList = paramList_;
422 paramList_ = Teuchos::null;
426Teuchos::RCP<const Teuchos::ParameterList>
427AmesosLinearOpWithSolveFactory::getParameterList()
const
432Teuchos::RCP<const Teuchos::ParameterList>
433AmesosLinearOpWithSolveFactory::getValidParameters()
const
435 return generateAndGetValidParameters();
440std::string AmesosLinearOpWithSolveFactory::description()
const
442 std::ostringstream oss;
443 oss <<
"Thyra::AmesosLinearOpWithSolveFactory{";
444 oss <<
"solverType=" << toString(solverType_);
452void 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");
462Teuchos::RCP<const Teuchos::ParameterList>
463AmesosLinearOpWithSolveFactory::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
471 ,Amesos::toString(Amesos::KLU)
473 ,Amesos::toString(Amesos::LAPACK)
476 validParamList->set(RefactorizationPolicy_name,Amesos::toString(Amesos::REPIVOT_ON_REFACTORIZATION));
477 validParamList->set(ThrowOnPreconditionerInput_name,
bool(
true));
480 return validParamList;