MueLu  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MueLu_MultiPhys_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // MueLu: A package for multigrid based preconditioning
4 //
5 // Copyright 2012 NTESS and the MueLu contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef MUELU_MULTIPHYS_DECL_HPP
11 #define MUELU_MULTIPHYS_DECL_HPP
12 
13 #include "MueLu_ConfigDefs.hpp"
14 #include "MueLu_BaseClass.hpp"
15 
16 #include "MueLu_SaPFactory_fwd.hpp"
17 
19 #include "MueLu_TrilinosSmoother.hpp"
20 #include "MueLu_Utilities_fwd.hpp"
21 #include "MueLu_Level_fwd.hpp"
22 #include "MueLu_Hierarchy_fwd.hpp"
23 #include "MueLu_RAPFactory_fwd.hpp"
24 #include "MueLu_PerfUtils_fwd.hpp"
25 #include "MueLu_SmootherBase.hpp"
26 
27 #include "Xpetra_Map_fwd.hpp"
28 #include "Xpetra_Matrix_fwd.hpp"
33 
34 namespace MueLu {
35 
42 template <class Scalar,
43  class LocalOrdinal,
44  class GlobalOrdinal,
45  class Node>
46 class MultiPhys : public VerboseObject, public Xpetra::Operator<Scalar, LocalOrdinal, GlobalOrdinal, Node> {
47 #undef MUELU_MULTIPHYS_SHORT
48 #include "MueLu_UseShortNames.hpp"
49 
50  public:
54 
57  : AmatMultiphysics_(Teuchos::null)
58  , hierarchyMultiphysics_(Teuchos::null)
59  , nBlks_(0) {
60  }
61 
74  MultiPhys(const Teuchos::RCP<Matrix>& AmatMultiPhysics,
75  const Teuchos::ArrayRCP<RCP<Matrix>> arrayOfAuxMatrices,
76  const Teuchos::ArrayRCP<Teuchos::RCP<MultiVector>> arrayOfNullspaces,
78  const int nBlks,
80  bool ComputePrec = true,
81  const Teuchos::ArrayRCP<Teuchos::RCP<MultiVector>> arrayOfMaterials = Teuchos::null,
82  bool OmitSubblockSmoother = true)
83  : AmatMultiphysics_(AmatMultiPhysics)
84  , arrayOfAuxMatrices_(arrayOfAuxMatrices)
85  , arrayOfNullspaces_(arrayOfNullspaces)
86  , arrayOfCoords_(arrayOfCoords)
87  , arrayOfMaterials_(arrayOfMaterials)
88  , OmitSubblockSmoother_(OmitSubblockSmoother)
89  , nBlks_(nBlks) {
90  initialize(AmatMultiPhysics, arrayOfAuxMatrices, arrayOfNullspaces, arrayOfCoords, nBlks, List, arrayOfMaterials);
91  compute(false);
92  }
93 
95  virtual ~MultiPhys() {}
96 
99 
101  const Teuchos::RCP<const Map> getRangeMap() const;
102 
105 
107  void compute(bool reuse = false);
108 
110  void resetMatrix(Teuchos::RCP<Matrix> SM_Matrix_new, bool ComputePrec = true);
111 
115  void apply(const MultiVector& X, MultiVector& Y,
119 
121  bool hasTransposeApply() const;
122 
124 
126  void residual(const MultiVector& X,
127  const MultiVector& B,
128  MultiVector& R) const {
129  using STS = Teuchos::ScalarTraits<Scalar>;
130  R.update(STS::one(), B, STS::zero());
131  this->apply(X, R, Teuchos::NO_TRANS, -STS::one(), STS::one());
132  }
133 
136  return hierarchyMultiphysics_;
137  }
138 
141  return arrayOfHierarchies_;
142  }
143 
144  private:
155  void initialize(const Teuchos::RCP<Matrix>& AmatMultiPhysics,
156  const Teuchos::ArrayRCP<RCP<Matrix>> arrayOfAuxMatrices,
157  const Teuchos::ArrayRCP<Teuchos::RCP<MultiVector>> arrayOfNullspaces,
159  const int nBlks,
161  const Teuchos::ArrayRCP<Teuchos::RCP<MultiVector>> arrayOfMaterials);
162 
164  void applyInverse(const MultiVector& RHS, MultiVector& X) const;
165 
167  void allocateMemory(int numVectors) const;
168 
170  Teuchos::RCP<Teuchos::TimeMonitor> getTimer(std::string name, RCP<const Teuchos::Comm<int>> comm = Teuchos::null) const;
171 
174 
176 
177  Teuchos::RCP<Matrix> AmatMultiphysics_; // multiphysics discretization matrix
178  Teuchos::RCP<Teuchos::ParameterList> paramListMultiphysics_; // array of parameter lists directing MueLu's construct of subblock P operators
179  Teuchos::RCP<Hierarchy> hierarchyMultiphysics_; // multiphysics discretization matrix
180 
181  Teuchos::ArrayRCP<Teuchos::RCP<Teuchos::ParameterList>> arrayOfParamLists_; // array of parameter lists directing MueLu's construct of subblock P operators
183  Teuchos::ArrayRCP<Teuchos::RCP<Matrix>> arrayOfAuxMatrices_; // array of discretization/auxiliary matrices used to generate subblock prolongators
184  Teuchos::ArrayRCP<Teuchos::RCP<MultiVector>> arrayOfNullspaces_; // array of nullspaces for smoothed aggregation.
185  Teuchos::ArrayRCP<Teuchos::RCP<RealValuedMultiVector>> arrayOfCoords_; // array of coordinates for smoothed aggregation/rebalancing.
186  Teuchos::ArrayRCP<Teuchos::RCP<MultiVector>> arrayOfMaterials_; // array of materials for smoothed aggregation.
187 
189 
190  int nBlks_; // number of PDE sub-systems within multiphysics system
192 };
193 
194 } // namespace MueLu
195 
196 #define MUELU_MULTIPHYS_SHORT
197 #endif // MUELU_MULTIPHYS_DECL_HPP
MultiPhys(const Teuchos::RCP< Matrix > &AmatMultiPhysics, const Teuchos::ArrayRCP< RCP< Matrix >> arrayOfAuxMatrices, const Teuchos::ArrayRCP< Teuchos::RCP< MultiVector >> arrayOfNullspaces, const Teuchos::ArrayRCP< Teuchos::RCP< RealValuedMultiVector >> arrayOfCoords, const int nBlks, Teuchos::ParameterList &List, bool ComputePrec=true, const Teuchos::ArrayRCP< Teuchos::RCP< MultiVector >> arrayOfMaterials=Teuchos::null, bool OmitSubblockSmoother=true)
MueLu::DefaultLocalOrdinal LocalOrdinal
Teuchos::RCP< Hierarchy > multiphysicsHierarchy()
Return the multiphysics hierarchy.
Preconditioner (wrapped as a Xpetra::Operator) for solving MultiPhysics PDEs.
Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > MultiVector
Xpetra::MultiVector< coordinateType, LO, GO, NO > RealValuedMultiVector
Teuchos::ArrayRCP< Teuchos::RCP< Matrix > > arrayOfAuxMatrices_
Teuchos::ScalarTraits< Scalar >::magnitudeType magnitudeType
MueLu::DefaultNode Node
virtual ~MultiPhys()
Destructor.
Verbose class for MueLu classes.
void setParameters(Teuchos::ParameterList &list)
Set parameters.
virtual void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)=0
Teuchos::ScalarTraits< Scalar >::coordinateType coordinateType
Teuchos::ArrayRCP< Teuchos::RCP< RealValuedMultiVector > > arrayOfCoords_
void applyInverse(const MultiVector &RHS, MultiVector &X) const
apply standard MultiPhys cycle
const Teuchos::RCP< const Map > getDomainMap() const
Returns the Xpetra::Map object associated with the domain of this operator.
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
const Teuchos::RCP< const Map > getRangeMap() const
Returns the Xpetra::Map object associated with the range of this operator.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::VERB_HIGH) const
Teuchos::RCP< Matrix > AmatMultiphysics_
Hierarchies: used to define P for (0,0)-block, .... (nBlks_-1,nBlks_-1) block.
MultiPhys()
Constructor.
bool hasTransposeApply() const
Indicates whether this operator supports applying the adjoint operator.
Teuchos::ArrayRCP< Teuchos::RCP< Teuchos::ParameterList > > arrayOfParamLists_
Teuchos::RCP< Teuchos::TimeMonitor > getTimer(std::string name, RCP< const Teuchos::Comm< int >> comm=Teuchos::null) const
get a (synced) timer
Teuchos::ParameterList parameterList_
ParameterLists.
Teuchos::RCP< Hierarchy > hierarchyMultiphysics_
Teuchos::ArrayRCP< Teuchos::RCP< Hierarchy > > arrayOfHierarchies_
void initialize(const Teuchos::RCP< Matrix > &AmatMultiPhysics, const Teuchos::ArrayRCP< RCP< Matrix >> arrayOfAuxMatrices, const Teuchos::ArrayRCP< Teuchos::RCP< MultiVector >> arrayOfNullspaces, const Teuchos::ArrayRCP< Teuchos::RCP< RealValuedMultiVector >> arrayOfCoords, const int nBlks, Teuchos::ParameterList &List, const Teuchos::ArrayRCP< Teuchos::RCP< MultiVector >> arrayOfMaterials)
void apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const
Teuchos::ArrayRCP< Teuchos::RCP< Hierarchy > > subblockHierarchies()
Return an array of hierarchies corresponding to each diagonal subblock.
Teuchos::RCP< Teuchos::ParameterList > paramListMultiphysics_
void residual(const MultiVector &X, const MultiVector &B, MultiVector &R) const
Compute a residual R = B - (*this) * X.
void resetMatrix(Teuchos::RCP< Matrix > SM_Matrix_new, bool ComputePrec=true)
Reset system matrix.
Teuchos::ArrayRCP< Teuchos::RCP< MultiVector > > arrayOfMaterials_
void compute(bool reuse=false)
Setup the preconditioner.
Teuchos::ArrayRCP< Teuchos::RCP< MultiVector > > arrayOfNullspaces_
void allocateMemory(int numVectors) const
allocate multivectors for solve