Teko  Version of the Day
 All Classes Files Functions Variables Pages
Teko_GaussSeidelPreconditionerFactory.cpp
1 /*
2 // @HEADER
3 //
4 // ***********************************************************************
5 //
6 // Teko: A package for block and physics based preconditioning
7 // Copyright 2010 Sandia Corporation
8 //
9 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
10 // the U.S. Government retains certain rights in this software.
11 //
12 // Redistribution and use in source and binary forms, with or without
13 // modification, are permitted provided that the following conditions are
14 // met:
15 //
16 // 1. Redistributions of source code must retain the above copyright
17 // notice, this list of conditions and the following disclaimer.
18 //
19 // 2. Redistributions in binary form must reproduce the above copyright
20 // notice, this list of conditions and the following disclaimer in the
21 // documentation and/or other materials provided with the distribution.
22 //
23 // 3. Neither the name of the Corporation nor the names of the
24 // contributors may be used to endorse or promote products derived from
25 // this software without specific prior written permission.
26 //
27 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
28 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
29 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
30 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
31 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
32 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
33 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
34 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
35 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
36 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
37 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
38 //
39 // Questions? Contact Eric C. Cyr (eccyr@sandia.gov)
40 //
41 // ***********************************************************************
42 //
43 // @HEADER
44 
45 */
46 
47 #include "Teko_GaussSeidelPreconditionerFactory.hpp"
48 
49 #include "Teko_BlockUpperTriInverseOp.hpp"
50 #include "Teko_BlockLowerTriInverseOp.hpp"
51 
52 using Teuchos::rcp;
53 using Teuchos::RCP;
54 
55 namespace Teko {
56 
57 GaussSeidelPreconditionerFactory::GaussSeidelPreconditionerFactory(TriSolveType solveType,const LinearOp & invD0,const LinearOp & invD1)
58  : invOpsStrategy_(rcp(new StaticInvDiagStrategy(invD0,invD1))), solveType_(solveType)
59 { }
60 
61 GaussSeidelPreconditionerFactory::GaussSeidelPreconditionerFactory(TriSolveType solveType,const RCP<const BlockInvDiagonalStrategy> & strategy)
62  : invOpsStrategy_(strategy), solveType_(solveType)
63 { }
64 
66  : solveType_(GS_UseLowerTriangle)
67 { }
68 
70 {
71  int rows = blockRowCount(blo);
72  int cols = blockColCount(blo);
73 
74  TEUCHOS_ASSERT(rows==cols);
75 
76  // get diagonal blocks
77  std::vector<LinearOp> invDiag;
78  invOpsStrategy_->getInvD(blo,state,invDiag);
79  TEUCHOS_ASSERT(rows==(int) invDiag.size());
80 
81  if(solveType_==GS_UseUpperTriangle) {
82  // create a blocked linear operator
83  BlockedLinearOp U = getUpperTriBlocks(blo);
84 
85  return createBlockUpperTriInverseOp(U,invDiag,"Gauss Seidel");
86  }
87  else if(solveType_==GS_UseLowerTriangle) {
88  // create a blocked linear operator
89  BlockedLinearOp L = getLowerTriBlocks(blo);
90 
91  return createBlockLowerTriInverseOp(L,invDiag,"Gauss Seidel");
92  }
93 
94  TEUCHOS_ASSERT(false); // we should never make it here!
95 }
96 
99 {
100  Teko_DEBUG_SCOPE("GaussSeidelPreconditionerFactory::initializeFromParameterList",10);
101  Teko_DEBUG_MSG_BEGIN(9);
102  DEBUG_STREAM << "Parameter list: " << std::endl;
103  pl.print(DEBUG_STREAM);
104  Teko_DEBUG_MSG_END();
105 
106  const std::string inverse_type = "Inverse Type";
107  std::vector<RCP<InverseFactory> > inverses;
108 
109  RCP<const InverseLibrary> invLib = getInverseLibrary();
110 
111  // get string specifying default inverse
112  std::string invStr ="Amesos";
113  if(pl.isParameter(inverse_type))
114  invStr = pl.get<std::string>(inverse_type);
115  if(pl.isParameter("Use Upper Triangle"))
116  solveType_ = pl.get<bool>("Use Upper Triangle") ? GS_UseUpperTriangle : GS_UseLowerTriangle;
117 
118  Teko_DEBUG_MSG("GSPrecFact: Building default inverse \"" << invStr << "\"",5);
119  RCP<InverseFactory> defaultInverse = invLib->getInverseFactory(invStr);
120 
121  // now check individual solvers
122  Teuchos::ParameterList::ConstIterator itr;
123  for(itr=pl.begin();itr!=pl.end();++itr) {
124  std::string fieldName = itr->first;
125  Teko_DEBUG_MSG("GSPrecFact: checking fieldName = \"" << fieldName << "\"",9);
126 
127  // figure out what the integer is
128  if(fieldName.compare(0,inverse_type.length(),inverse_type)==0 && fieldName!=inverse_type) {
129  int position = -1;
130  std::string inverse,type;
131 
132  // figure out position
133  std::stringstream ss(fieldName);
134  ss >> inverse >> type >> position;
135 
136  if(position<=0) {
137  Teko_DEBUG_MSG("Gauss-Seidel \"Inverse Type\" must be a (strictly) positive integer",1);
138  }
139 
140  // inserting inverse factory into vector
141  std::string invStr = pl.get<std::string>(fieldName);
142  Teko_DEBUG_MSG("GSPrecFact: Building inverse " << position << " \"" << invStr << "\"",5);
143  if(position>(int) inverses.size()) {
144  inverses.resize(position,defaultInverse);
145  inverses[position-1] = invLib->getInverseFactory(invStr);
146  }
147  else
148  inverses[position-1] = invLib->getInverseFactory(invStr);
149  }
150  }
151 
152  // use default inverse
153  if(inverses.size()==0)
154  inverses.push_back(defaultInverse);
155 
156  // based on parameter type build a strategy
157  invOpsStrategy_ = rcp(new InvFactoryDiagStrategy(inverses,defaultInverse));
158 }
159 
160 } // end namspace Teko
An implementation of a state object for block preconditioners.
Teuchos::RCP< const BlockInvDiagonalStrategy > invOpsStrategy_
some members
LinearOp buildPreconditionerOperator(BlockedLinearOp &blo, BlockPreconditionerState &state) const
Create the Gauss-Seidel preconditioner operator.
Teuchos::RCP< const InverseLibrary > getInverseLibrary() const
Get the inverse library used by this preconditioner factory.
virtual void initializeFromParameterList(const Teuchos::ParameterList &pl)
Initialize from a parameter list.