47 #include "Teko_TimingsSIMPLEPreconditionerFactory.hpp"
50 #include "Teko_InverseFactory.hpp"
51 #include "Teko_BlockLowerTriInverseOp.hpp"
52 #include "Teko_BlockUpperTriInverseOp.hpp"
53 #ifdef TEKO_HAVE_EPETRA
54 #include "Teko_DiagonalPreconditionerFactory.hpp"
57 #include "Teuchos_Time.hpp"
66 const RCP<InverseFactory>& inverse,
double alpha)
67 : invVelFactory_(inverse),
68 invPrsFactory_(inverse),
70 fInverseType_(Diagonal),
72 constrTotal_(
"SIMPLE Constr: Total"),
73 subTotal_(
"SIMPLE Constr: Subs"),
76 TimingsSIMPLEPreconditionerFactory ::TimingsSIMPLEPreconditionerFactory(
77 const RCP<InverseFactory>& invVFact,
const RCP<InverseFactory>& invPFact,
double alpha)
78 : invVelFactory_(invVFact),
79 invPrsFactory_(invPFact),
81 fInverseType_(Diagonal),
83 constrTotal_(
"SIMPLE Constr: Total"),
84 subTotal_(
"SIMPLE Constr: Subs"),
87 TimingsSIMPLEPreconditionerFactory::TimingsSIMPLEPreconditionerFactory()
89 fInverseType_(Diagonal),
91 constrTotal_(
"SIMPLE Constr: Total"),
92 subTotal_(
"SIMPLE Constr: Subs"),
96 if (constrTotal_.totalElapsedTime() > 0.0) {
97 Teuchos::FancyOStream out(Teuchos::rcpFromRef(std::cout));
98 out.setOutputToRootOnly(0);
100 out <<
"==========================================================================="
103 out <<
"SIMPLE Construction Count = " << constrCount_ << std::endl;
104 out <<
"SIMPLE Construction Total = " << constrTotal_.totalElapsedTime() << std::endl;
105 out <<
"SIMPLE Sub Components Total = " << subTotal_.totalElapsedTime() << std::endl;
107 out <<
"==========================================================================="
115 constrTotal_.start();
117 Teko_DEBUG_SCOPE(
"TimingsSIMPLEPreconditionerFactory::buildPreconditionerOperator", 10);
118 Teko_DEBUG_EXPR(Teuchos::Time timer(
""));
120 int rows = blockRowCount(blockOp);
121 int cols = blockColCount(blockOp);
123 TEUCHOS_ASSERT(rows == 2);
124 TEUCHOS_ASSERT(cols == 2);
126 bool buildExplicitSchurComplement =
true;
129 const LinearOp F = getBlock(0, 0, blockOp);
130 const LinearOp Bt = getBlock(0, 1, blockOp);
131 const LinearOp B = getBlock(1, 0, blockOp);
132 const LinearOp C = getBlock(1, 1, blockOp);
136 TEUCHOS_ASSERT(massMatrix_ != Teuchos::null);
141 std::string fApproxStr =
"<error>";
143 if (fInverseType_ == NotDiag) {
145 fApproxStr = customHFactory_->toString();
148 buildExplicitSchurComplement =
false;
149 }
else if (fInverseType_ == BlkDiag) {
150 #ifdef TEKO_HAVE_EPETRA
157 buildExplicitSchurComplement =
true;
160 throw std::logic_error(
161 "TimingsSIMPLEPreconditionerFactory fInverseType_ == "
162 "BlkDiag but EPETRA is turned off!");
167 H = getInvDiagonalOp(matF, fInverseType_);
169 fApproxStr = getDiagonalName(fInverseType_);
174 RCP<const Teuchos::ParameterList> pl = state.getParameterList();
176 if (pl->isParameter(
"stepsize")) {
178 double stepsize = pl->get<
double>(
"stepsize");
181 if (stepsize > 0.0) H = scale(stepsize, H);
188 if (buildExplicitSchurComplement) {
195 mHBt = explicitMultiply(H, Bt, mHBt);
201 BHBt = explicitMultiply(B, HBt, BHBt);
206 mhatS = explicitAdd(C, scale(-1.0, BHBt), mhatS);
211 HBt = multiply(H, Bt);
213 hatS = add(C, scale(-1.0, multiply(B, HBt)));
217 if (timed_HBt_ == Teuchos::null) {
220 timed_HBt_->setLinearOp(HBt);
224 if (timed_B_ == Teuchos::null) {
227 timed_B_->setLinearOp(B);
233 if (invF == Teuchos::null) {
240 timed_invF_->setLinearOp(invF);
247 if (invS == Teuchos::null) {
254 timed_invS_->setLinearOp(invS);
258 std::vector<LinearOp> invDiag(2);
261 BlockedLinearOp L = zeroBlockedOp(blockOp);
262 setBlock(1, 0, L, timed_B_);
265 invDiag[0] = timed_invF_;
266 invDiag[1] = timed_invS_;
267 LinearOp invL = createBlockLowerTriInverseOp(L, invDiag);
270 BlockedLinearOp U = zeroBlockedOp(blockOp);
271 setBlock(0, 1, U, scale<double>(1.0 / alpha_, timed_HBt_.getConst()));
274 invDiag[0] = identity(rangeSpace(invF));
275 invDiag[1] = scale(alpha_, identity(rangeSpace(invS)));
276 LinearOp invU = createBlockUpperTriInverseOp(U, invDiag);
279 Teko::LinearOp iU_t_iL = multiply(invU, invL,
"SIMPLE_" + fApproxStr);
282 if (timed_iU_t_iL_ == Teuchos::null)
283 timed_iU_t_iL_ = Teuchos::rcp(
new DiagnosticLinearOp(getOutputStream(), iU_t_iL,
"iU_t_iL"));
285 timed_iU_t_iL_->setLinearOp(iU_t_iL);
291 return timed_iU_t_iL_;
296 const Teuchos::ParameterList& pl) {
301 customHFactory_ = Teuchos::null;
302 fInverseType_ = Diagonal;
305 std::string invStr =
"", invVStr =
"", invPStr =
"";
309 if (pl.isParameter(
"Inverse Type")) invStr = pl.get<std::string>(
"Inverse Type");
310 if (pl.isParameter(
"Inverse Velocity Type"))
311 invVStr = pl.get<std::string>(
"Inverse Velocity Type");
312 if (pl.isParameter(
"Inverse Pressure Type"))
313 invPStr = pl.get<std::string>(
"Inverse Pressure Type");
314 if (pl.isParameter(
"Alpha")) alpha_ = pl.get<
double>(
"Alpha");
315 if (pl.isParameter(
"Explicit Velocity Inverse Type")) {
316 std::string fInverseStr = pl.get<std::string>(
"Explicit Velocity Inverse Type");
319 fInverseType_ = getDiagonalType(fInverseStr);
320 if (fInverseType_ == NotDiag) customHFactory_ = invLib->getInverseFactory(fInverseStr);
323 if (fInverseType_ == BlkDiag) BlkDiagList_ = pl.sublist(
"H options");
325 if (pl.isParameter(
"Use Mass Scaling")) useMass_ = pl.get<
bool>(
"Use Mass Scaling");
327 Teko_DEBUG_MSG_BEGIN(5) DEBUG_STREAM <<
"SIMPLE Parameters: " << std::endl;
328 DEBUG_STREAM <<
" inv type = \"" << invStr <<
"\"" << std::endl;
329 DEBUG_STREAM <<
" inv v type = \"" << invVStr <<
"\"" << std::endl;
330 DEBUG_STREAM <<
" inv p type = \"" << invPStr <<
"\"" << std::endl;
331 DEBUG_STREAM <<
" alpha = " << alpha_ << std::endl;
332 DEBUG_STREAM <<
" use mass = " << useMass_ << std::endl;
333 DEBUG_STREAM <<
" vel scaling = " << getDiagonalName(fInverseType_) << std::endl;
334 DEBUG_STREAM <<
"SIMPLE Parameter list: " << std::endl;
335 pl.print(DEBUG_STREAM);
339 #if defined(Teko_ENABLE_Amesos)
340 if (invStr ==
"") invStr =
"Amesos";
341 #elif defined(Teko_ENABLE_Amesos2)
342 if (invStr ==
"") invStr =
"Amesos2";
344 if (invVStr ==
"") invVStr = invStr;
345 if (invPStr ==
"") invPStr = invStr;
348 RCP<InverseFactory> invVFact, invPFact;
351 invVFact = invLib->getInverseFactory(invVStr);
353 if (invVStr != invPStr)
354 invPFact = invLib->getInverseFactory(invPStr);
357 invVelFactory_ = invVFact;
358 invPrsFactory_ = invPFact;
362 rh->preRequest<Teko::LinearOp>(Teko::RequestMesg(
"Velocity Mass Matrix"));
363 Teko::LinearOp mass = rh->request<Teko::LinearOp>(Teko::RequestMesg(
"Velocity Mass Matrix"));
371 Teuchos::RCP<Teuchos::ParameterList> result;
372 Teuchos::RCP<Teuchos::ParameterList> pl = rcp(
new Teuchos::ParameterList());
375 RCP<Teuchos::ParameterList> vList = invVelFactory_->getRequestedParameters();
376 if (vList != Teuchos::null) {
377 Teuchos::ParameterList::ConstIterator itr;
378 for (itr = vList->begin(); itr != vList->end(); ++itr) pl->setEntry(itr->first, itr->second);
383 RCP<Teuchos::ParameterList> pList = invPrsFactory_->getRequestedParameters();
384 if (pList != Teuchos::null) {
385 Teuchos::ParameterList::ConstIterator itr;
386 for (itr = pList->begin(); itr != pList->end(); ++itr) pl->setEntry(itr->first, itr->second);
391 if (customHFactory_ != Teuchos::null) {
392 RCP<Teuchos::ParameterList> hList = customHFactory_->getRequestedParameters();
393 if (hList != Teuchos::null) {
394 Teuchos::ParameterList::ConstIterator itr;
395 for (itr = hList->begin(); itr != hList->end(); ++itr) pl->setEntry(itr->first, itr->second);
405 const Teuchos::ParameterList& pl) {
406 Teko_DEBUG_SCOPE(
"InvLSCStrategy::updateRequestedParameters", 10);
410 result &= invVelFactory_->updateRequestedParameters(pl);
411 result &= invPrsFactory_->updateRequestedParameters(pl);
412 if (customHFactory_ != Teuchos::null) result &= customHFactory_->updateRequestedParameters(pl);
void rebuildInverse(const InverseFactory &factory, const LinearOp &A, InverseLinearOp &invA)
virtual bool updateRequestedParameters(const Teuchos::ParameterList &pl)
For assisting in construction of the preconditioner.
virtual ~TimingsSIMPLEPreconditionerFactory()
Destructor that outputs construction timings.
TimingsSIMPLEPreconditionerFactory()
Default constructor.
virtual void setMassMatrix(Teko::LinearOp &mass)
Set the mass matrix for this factory.
virtual Teuchos::RCP< Teuchos::ParameterList > getRequestedParameters() const
For assisting in construction of the preconditioner.
virtual void initializeFromParameterList(const Teuchos::ParameterList &pl)
Initialize from a parameter list.
Preconditioner factory for building explcit inverse of diagonal operators. This includes block operat...
An implementation of a state object for block preconditioners.
This linear operator prints diagnostics about operator application and creation times. It is useful for debugging problems and determining bottle necks.
InverseLinearOp buildInverse(const InverseFactory &factory, const LinearOp &A)
Build an inverse operator using a factory and a linear operator.
virtual void initializeFromParameterList(const Teuchos::ParameterList &pl)
Initialize from a parameter list.
Teuchos::RCP< const InverseLibrary > getInverseLibrary() const
Get the inverse library used by this preconditioner factory.
LinearOp buildPreconditionerOperator(BlockedLinearOp &blo, BlockPreconditionerState &state) const
Function that is called to build the preconditioner for the linear operator that is passed in...
Teuchos::RCP< RequestHandler > getRequestHandler() const
Get the request handler with pointers to the appropriate callbacks.
LinearOp buildPreconditionerOperator(LinearOp &lo, PreconditionerState &state) const
virtual Teko::ModifiableLinearOp & getModifiableOp(const std::string &name)
Add a named operator to the state object.