Teko Version of the Day
Loading...
Searching...
No Matches
Teko_ModALPreconditionerFactory.cpp
1/*
2 * Author: Zhen Wang
3 * Email: wangz@ornl.gov
4 * zhen.wang@alum.emory.edu
5 */
6
7#include "Teko_ModALPreconditionerFactory.hpp"
8
9#include "Thyra_DefaultMultipliedLinearOp.hpp"
10#include "Thyra_DefaultAddedLinearOp.hpp"
11#include "Thyra_DefaultIdentityLinearOp.hpp"
12#include "Thyra_DefaultZeroLinearOp.hpp"
13
15#include "Teko_Utilities.hpp"
16#include "Teko_BlockLowerTriInverseOp.hpp"
17#include "Teko_BlockUpperTriInverseOp.hpp"
18#include "Teko_StaticLSCStrategy.hpp"
19#include "Teko_InvLSCStrategy.hpp"
20#include "Teko_PresLaplaceLSCStrategy.hpp"
21
22#include "Teuchos_Time.hpp"
23
24namespace Teko
25{
26
27namespace NS
28{
29
30ModALPrecondState::ModALPrecondState():
31pressureMassMatrix_(Teuchos::null), invPressureMassMatrix_(Teuchos::null)
32{
33}
34
35ModALPreconditionerFactory::ModALPreconditionerFactory(const Teuchos::RCP<InverseFactory> & factory) :
36 invOpsStrategy_(Teuchos::rcp(new InvModALStrategy(factory))),
37 isSymmetric_(true)
38{
39}
40
41ModALPreconditionerFactory::ModALPreconditionerFactory(const Teuchos::RCP<InverseFactory> & invFactoryA,
42 const Teuchos::RCP<InverseFactory> & invFactoryS) :
43 invOpsStrategy_(Teuchos::rcp(new InvModALStrategy(invFactoryA, invFactoryS))),
44 isSymmetric_(true)
45{
46}
47
48ModALPreconditionerFactory::ModALPreconditionerFactory(const Teuchos::RCP<InverseFactory> & factory,
49 LinearOp & pressureMassMatrix) :
50 invOpsStrategy_(Teuchos::rcp(new InvModALStrategy(factory, pressureMassMatrix))), isSymmetric_(true)
51{
52}
53
54ModALPreconditionerFactory::ModALPreconditionerFactory(const Teuchos::RCP<InverseFactory> & invFactoryA,
55 const Teuchos::RCP<InverseFactory> & invFactoryS,
56 LinearOp & pressureMassMatrix) :
57 invOpsStrategy_(Teuchos::rcp(new InvModALStrategy(invFactoryA, invFactoryS, pressureMassMatrix))),
58 isSymmetric_(true)
59{
60}
61
62// Construct from a strategy.
63ModALPreconditionerFactory::ModALPreconditionerFactory(const Teuchos::RCP<InvModALStrategy> & strategy) :
64 invOpsStrategy_(strategy), isSymmetric_(true)
65{
66}
67
68LinearOp ModALPreconditionerFactory::buildPreconditionerOperator(BlockedLinearOp & alOp,
69 BlockPreconditionerState & state) const
70{
71 Teko_DEBUG_SCOPE("ModALPreconditionerFactory::buildPreconditionerOperator()", 10);
72 Teko_DEBUG_EXPR(Teuchos::Time timer(""));
73 Teko_DEBUG_EXPR(Teuchos::Time totalTimer(""));
74 Teko_DEBUG_EXPR(totalTimer.start());
75
76 // Only for 2D or 3D problems.
77 int dim = blockRowCount(alOp) - 1;
78 TEUCHOS_ASSERT(dim == 2 || dim == 3);
79
80 // Build what is necessary for the state object.
81 Teko_DEBUG_EXPR(timer.start(true));
82 invOpsStrategy_->buildState(alOp, state);
83 Teko_DEBUG_EXPR(timer.stop());
84 Teko_DEBUG_MSG("ModALPreconditionerFactory::buildPreconditionerOperator():BuildStateTime = "
85 << timer.totalElapsedTime(), 2);
86
87 // Extract inverse operators from strategy
88 Teko_DEBUG_EXPR(timer.start(true));
89 LinearOp invA11p = invOpsStrategy_->getInvA11p(state);
90 LinearOp invA22p = invOpsStrategy_->getInvA22p(state);
91 LinearOp invA33p;
92 if(dim == 3)
93 {
94 invA33p = invOpsStrategy_->getInvA33p(state);
95 }
96
97 // The inverse of S can be built from strategy,
98 // or just a diagonal matrix.
99 ModALPrecondState * modALState = dynamic_cast<ModALPrecondState*>(&state);
100 TEUCHOS_ASSERT(modALState != NULL);
101 LinearOp invS;
102 if(modALState->isStabilized_)
103 {
104 invS = invOpsStrategy_->getInvS(state);
105 }
106 else
107 {
108 invS = scale(modALState->gamma_, modALState->invPressureMassMatrix_);
109 }
110
111 Teko_DEBUG_EXPR(timer.stop());
112 Teko_DEBUG_MSG("ModALPrecFact::buildPreconditionerOperator(): GetInvTime = "
113 << timer.totalElapsedTime(), 2);
114
115 // Build diagonal operations.
116 std::vector<LinearOp> invDiag;
117 invDiag.resize(dim + 1);
118 invDiag[0] = invA11p;
119 invDiag[1] = invA22p;
120 if(dim == 2)
121 {
122 invDiag[2] = scale(-1.0, invS);
123 }
124 else if(dim == 3)
125 {
126 invDiag[2] = invA33p;
127 invDiag[3] = scale(-1.0, invS);
128 }
129
130 // Get the upper triangular matrix.
131 BlockedLinearOp U = getUpperTriBlocks(alOp);
132
133 Teko_DEBUG_EXPR(totalTimer.stop());
134 Teko_DEBUG_MSG("ModALPrecFact::buildPreconditionerOperator TotalTime = "
135 << totalTimer.totalElapsedTime(), 2);
136
137 // Create the preconditioner
138 return createBlockUpperTriInverseOp(U, invDiag, "Modified AL preconditioner-Upper");
139}
140
141} // end namespace NS
142
143} // end namespace Teko
int blockRowCount(const BlockedLinearOp &blo)
Get the row count in a block linear operator.
An implementation of a state object for block preconditioners.
Class for saving state variables for ModALPreconditionerFactory.