47 #include "Teko_SIMPLEPreconditionerFactory.hpp" 50 #include "Teko_InverseFactory.hpp" 51 #include "Teko_BlockLowerTriInverseOp.hpp" 52 #include "Teko_BlockUpperTriInverseOp.hpp" 53 #include "Teko_DiagonalPreconditionerFactory.hpp" 55 #include "Teuchos_Time.hpp" 56 #include "Epetra_FECrsGraph.h" 57 #include "Epetra_FECrsMatrix.h" 58 #include "EpetraExt_PointToBlockDiagPermute.h" 59 #include "EpetraExt_MatrixMatrix.h" 60 #include "Thyra_EpetraOperatorWrapper.hpp" 61 #include "Thyra_EpetraLinearOp.hpp" 73 : invVelFactory_(inverse), invPrsFactory_(inverse), alpha_(alpha), fInverseType_(
Diagonal), useMass_(false)
78 const RCP<InverseFactory> & invPFact,
80 : invVelFactory_(invVFact), invPrsFactory_(invPFact), alpha_(alpha), fInverseType_(
Diagonal), useMass_(false)
84 : alpha_(1.0), fInverseType_(
Diagonal), useMass_(false)
92 Teko_DEBUG_SCOPE(
"SIMPLEPreconditionerFactory::buildPreconditionerOperator",10);
93 Teko_DEBUG_EXPR(Teuchos::Time timer(
""));
98 TEUCHOS_ASSERT(rows==2);
99 TEUCHOS_ASSERT(cols==2);
101 bool buildExplicitSchurComplement =
true;
104 const LinearOp F =
getBlock(0,0,blockOp);
105 const LinearOp Bt =
getBlock(0,1,blockOp);
106 const LinearOp B =
getBlock(1,0,blockOp);
107 const LinearOp C =
getBlock(1,1,blockOp);
111 TEUCHOS_ASSERT(massMatrix_!=Teuchos::null);
116 std::string fApproxStr =
"<error>";
120 fApproxStr = customHFactory_->toString();
123 buildExplicitSchurComplement =
false;
125 else if(fInverseType_==
BlkDiag) {
138 buildExplicitSchurComplement =
true;
143 H = getInvDiagonalOp(matF,fInverseType_);
144 fApproxStr = getDiagonalName(fInverseType_);
149 RCP<const Teuchos::ParameterList> pl = state.getParameterList();
151 if(pl->isParameter(
"stepsize")) {
153 double stepsize = pl->get<
double>(
"stepsize");
157 H =
scale(stepsize,H);
164 if(buildExplicitSchurComplement) {
170 mHBt = explicitMultiply(H,Bt,mHBt);
174 BHBt = explicitMultiply(B,HBt,BHBt);
177 mhatS = explicitAdd(C,
scale(-1.0,BHBt),mhatS);
182 HBt = multiply(H,Bt);
184 hatS = add(C,
scale(-1.0,multiply(B,HBt)));
189 if(invF==Teuchos::null)
196 if(invS==Teuchos::null)
201 std::vector<LinearOp> invDiag(2);
204 BlockedLinearOp L = zeroBlockedOp(blockOp);
210 LinearOp invL = createBlockLowerTriInverseOp(L,invDiag);
213 BlockedLinearOp U = zeroBlockedOp(blockOp);
219 LinearOp invU = createBlockUpperTriInverseOp(U,invDiag);
222 return multiply(invU,invL,
"SIMPLE_"+fApproxStr);
232 customHFactory_ = Teuchos::null;
236 std::string invStr=
"", invVStr=
"", invPStr=
"";
240 if(pl.isParameter(
"Inverse Type"))
241 invStr = pl.get<std::string>(
"Inverse Type");
242 if(pl.isParameter(
"Inverse Velocity Type"))
243 invVStr = pl.get<std::string>(
"Inverse Velocity Type");
244 if(pl.isParameter(
"Inverse Pressure Type"))
245 invPStr = pl.get<std::string>(
"Inverse Pressure Type");
246 if(pl.isParameter(
"Alpha"))
247 alpha_ = pl.get<
double>(
"Alpha");
248 if(pl.isParameter(
"Explicit Velocity Inverse Type")) {
249 std::string fInverseStr = pl.get<std::string>(
"Explicit Velocity Inverse Type");
252 fInverseType_ = getDiagonalType(fInverseStr);
254 customHFactory_ = invLib->getInverseFactory(fInverseStr);
258 BlkDiagList_=pl.sublist(
"H options");
260 if(pl.isParameter(
"Use Mass Scaling"))
261 useMass_ = pl.get<
bool>(
"Use Mass Scaling");
263 Teko_DEBUG_MSG_BEGIN(5)
264 DEBUG_STREAM <<
"SIMPLE Parameters: " << std::endl;
265 DEBUG_STREAM <<
" inv type = \"" << invStr <<
"\"" << std::endl;
266 DEBUG_STREAM <<
" inv v type = \"" << invVStr <<
"\"" << std::endl;
267 DEBUG_STREAM <<
" inv p type = \"" << invPStr <<
"\"" << std::endl;
268 DEBUG_STREAM <<
" alpha = " << alpha_ << std::endl;
269 DEBUG_STREAM <<
" use mass = " << useMass_ << std::endl;
270 DEBUG_STREAM <<
" vel scaling = " << getDiagonalName(fInverseType_) << std::endl;
271 DEBUG_STREAM <<
"SIMPLE Parameter list: " << std::endl;
272 pl.print(DEBUG_STREAM);
276 if(invStr==
"") invStr =
"Amesos";
277 if(invVStr==
"") invVStr = invStr;
278 if(invPStr==
"") invPStr = invStr;
281 RCP<InverseFactory> invVFact, invPFact;
284 invVFact = invLib->getInverseFactory(invVStr);
287 invPFact = invLib->getInverseFactory(invPStr);
290 invVelFactory_ = invVFact;
291 invPrsFactory_ = invPFact;
295 rh->preRequest<Teko::LinearOp>(Teko::RequestMesg(
"Velocity Mass Matrix"));
297 = rh->request<Teko::LinearOp>(Teko::RequestMesg(
"Velocity Mass Matrix"));
305 Teuchos::RCP<Teuchos::ParameterList> result;
306 Teuchos::RCP<Teuchos::ParameterList> pl = rcp(
new Teuchos::ParameterList());
309 RCP<Teuchos::ParameterList> vList = invVelFactory_->getRequestedParameters();
310 if(vList!=Teuchos::null) {
311 Teuchos::ParameterList::ConstIterator itr;
312 for(itr=vList->begin();itr!=vList->end();++itr)
313 pl->setEntry(itr->first,itr->second);
318 RCP<Teuchos::ParameterList> pList = invPrsFactory_->getRequestedParameters();
319 if(pList!=Teuchos::null) {
320 Teuchos::ParameterList::ConstIterator itr;
321 for(itr=pList->begin();itr!=pList->end();++itr)
322 pl->setEntry(itr->first,itr->second);
327 if(customHFactory_!=Teuchos::null) {
328 RCP<Teuchos::ParameterList> hList = customHFactory_->getRequestedParameters();
329 if(hList!=Teuchos::null) {
330 Teuchos::ParameterList::ConstIterator itr;
331 for(itr=hList->begin();itr!=hList->end();++itr)
332 pl->setEntry(itr->first,itr->second);
343 Teko_DEBUG_SCOPE(
"InvLSCStrategy::updateRequestedParameters",10);
347 result &= invVelFactory_->updateRequestedParameters(pl);
348 result &= invPrsFactory_->updateRequestedParameters(pl);
349 if(customHFactory_!=Teuchos::null)
350 result &= customHFactory_->updateRequestedParameters(pl);
Specifies that a block diagonal approximation is to be used.
Specifies that just the diagonal is used.
int blockRowCount(const BlockedLinearOp &blo)
Get the row count in a block linear operator.
InverseLinearOp buildInverse(const InverseFactory &factory, const LinearOp &A)
Build an inverse operator using a factory and a linear operator.
void scale(const double alpha, MultiVector &x)
Scale a multivector by a constant.
void endBlockFill(BlockedLinearOp &blo)
Notify the blocked operator that the fill stage is completed.
virtual bool updateRequestedParameters(const Teuchos::ParameterList &pl)
For assisting in construction of the preconditioner.
Preconditioner factory for building explcit inverse of diagonal operators. This includes block operat...
void rebuildInverse(const InverseFactory &factory, const LinearOp &A, InverseLinearOp &invA)
An implementation of a state object for block preconditioners.
VectorSpace rangeSpace(const LinearOp &lo)
Get the range space of a linear operator.
MultiVector getBlock(int i, const BlockedMultiVector &bmv)
Get the ith block from a BlockedMultiVector object.
virtual void setMassMatrix(Teko::LinearOp &mass)
Set the mass matrix for this factory.
For user convenience, if Teko recieves this value, exceptions will be thrown.
virtual void initializeFromParameterList(const Teuchos::ParameterList &pl)
Initialize from a parameter list.
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.
int blockColCount(const BlockedLinearOp &blo)
Get the column count in a block linear operator.
virtual Teuchos::RCP< Teuchos::ParameterList > getRequestedParameters() const
For assisting in construction of the preconditioner.
void setBlock(int i, int j, BlockedLinearOp &blo, const LinearOp &lo)
Set the i,j block in a BlockedLinearOp object.
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.
SIMPLEPreconditionerFactory()
Default constructor.
LinearOp buildPreconditionerOperator(LinearOp &lo, PreconditionerState &state) const
virtual Teko::ModifiableLinearOp & getModifiableOp(const std::string &name)
Add a named operator to the state object.