v0.9.1
PCMGSetUpViaApproxOrders.hpp
Go to the documentation of this file.
1 /** \file PCMGSetUpViaApproxOrders.hpp
2  * \brief header of multi-grid solver for p- adaptivity
3  *
4  * MoFEM is free software: you can redistribute it and/or modify it under
5  * the terms of the GNU Lesser General Public License as published by the
6  * Free Software Foundation, either version 3 of the License, or (at your
7  * option) any later version.
8  *
9  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
11  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
12  * License for more details.
13  *
14  * You should have received a copy of the GNU Lesser General Public
15  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
16 */
17 
18 #ifndef __PCMGSETUP_VIA_APPROX_ORDERS_HPP__
19 #define __PCMGSETUP_VIA_APPROX_ORDERS_HPP__
20 
21 static const int DMMGVIAAPPROXORDERSCTX_INTERFACE = 1<<2;
23 
25  Mat A;
26  Vec X,F;
27  IS iS;
28  VecScatter sCat;
29  PCMGSubMatrixCtx(Mat a,IS is);
30  virtual ~PCMGSubMatrixCtx();
31 };
32 
33 /**
34  * \brief Structure for DM for multi-grid via approximation orders
35  * \ingroup dm
36  */
38 
39  MoFEMErrorCode query_interface(const MOFEMuuid& uuid,MoFEM::UnknownInterface** iface) const;
40 
42  virtual ~DMMGViaApproxOrdersCtx();
43 
45 
46  AO aO;
47  std::vector<IS> coarseningIS; ///< Coarsening IS
48  std::vector<Mat> kspOperators; ///< Get KSP operators
49  boost::ptr_vector<PCMGSubMatrixCtx> shellMatrixCtxPtr; ///< Shell sub-matrix context
50 
51 };
52 
53 /**
54  * Get DM Ctx
55  */
57 
58 /**
59  * \brief Set DM ordering
60  *
61  * IS can be given is some other ordering, AO will transform indices from coarseningIS ordering
62  * to ordering used to construct fine matrix.
63  *
64  * @param dm [description]
65  * @param ao [description]
66  * @return [description]
67  */
69 
70 /**
71  * \brief Gets size of coarseningIS in internal data struture DMMGViaApproxOrdersCtx
72  * @param dm DM
73  * @param size size of coarseningIS
74  * @return Error code
75  */
77 
78 /**
79  * \brief Push back coarsening level to MG via approximation orders
80  *
81  * @param DM discrete manager
82  * @param is Push back IS used for coarsening
83  * @param A Get sub-matrix of A using is (that sets operators for coarsening levels)
84  * @param subA Returning pointer to created sub matrix
85  * @param subA If true create sub matrix, otherwise in subA has to be valid pointer to subA
86  * @return Error code
87  *
88  * \ingroup dm
89  */
91  DM,IS is,Mat A,Mat *subA,bool create_sub_matrix,bool shell_sub_a
92 );
93 
94 /**
95  * \brief Pop is form MG via approximation orders
96  * @param DM dm
97  * @param is pop back IS
98  * @return error code
99  *
100  * \ingroup dm
101  */
103 
104 /**
105  * \brief Clear approximation orders
106  * @param DM dm
107  * @return Error code
108  *
109  * \ingroup dm
110  */
112 
113 /**
114  * \brief Replace coarsening IS in DM via approximation orders
115  * @param dm dm
116  * @param is_vec Pointer to vector of is
117  * @param nb_elems Number of elements
118  * @param A Fine matrix
119  * @return Error code
120  */
122  DM dm,IS *is_vec,int nb_elems,Mat A,int verb = 0
123 );
124 
125 /**
126  * \brief Get context for DM via approximation orders
127  * @param dm the DM object
128  * @param ctx data context
129  * @return error code
130  */
132 
133 /**
134  * \brief Register DM for Multi-Grid via approximation orders
135  * @param sname problem/dm registered name
136  * @return error code
137  * \ingroup dm
138  */
139 MoFEMErrorCode DMRegister_MGViaApproxOrders(const char sname[]);
140 
141 /**
142  * \brief Create DM data structure for Multi-Grid via approximation orders
143  *
144  * It set data structure and operators needed
145  *
146  * @param dm Discrete manager
147  * @return Error code
148  */
150 
151 /**
152  * \brief Create matrix for Multi-Grid via approximation orders
153  *
154  * Not used directly by user
155  *
156  * @param dm Distributed mesh data structure
157  * @param M Matrix
158  * @return Error code
159  * \ingroup dm
160  */
162 
163 /**
164  * \brief Coarsen DM
165  *
166  * Not used directly by user
167  *
168  * @param dm Distributed mesh data structure
169  * @param comm Communicator
170  * @param dmc Coarse distributed mesh data structure
171  * @return Error code
172  *
173  * \ingroup dm
174  */
175 MoFEMErrorCode DMCoarsen_MGViaApproxOrders(DM dm, MPI_Comm comm, DM *dmc);
176 
177 /**
178  * \brief Create interpolation matrix between data managers dm1 and dm2
179  * @param dm1 Distributed mesh data structure
180  * @param dm2 Distributed mesh data structure
181  * @param mat Pointer to returned interpolation matrix
182  * @param vec Pointer to scaling vector here returned NULL
183  * @return Error code
184  */
185 MoFEMErrorCode DMCreateInterpolation_MGViaApproxOrders(DM dm1,DM dm2,Mat *mat,Vec *vec);
186 
187 /**
188  * \brief Create global vector for DMGViaApproxOrders
189  * @param dm Distributed mesh data structure
190  * @param g returned pointer to vector
191  * @return Error code
192  */
194 
195 /**
196  * \brief Set data structures of MG pre-conditioner via approximation orders
197  */
199 
200  // Interface *mFieldPtr; ///< MoFEM interface
201  // string problemName; ///< Problem name
202 
203  DM dM; ///< Distributed mesh manager
204  Mat A; ///< Matrix at fine level
205 
207  DM dm,Mat a,bool shell_sub_a
208  ):
209  // mFieldPtr(mfield_ptr),
210  // problemName(problem_name),
211  dM(dm),
212  A(a),
213  nbLevels(2),
214  coarseOrder(2),
215  orderAtLastLevel(1000),
216  shellSubA(shell_sub_a),
217  verboseLevel(0) {
218  }
219 
221  }
222 
223  int nbLevels; ///< number of multi-grid levels
224  int coarseOrder; ///< approximation order of coarse level
225  int orderAtLastLevel; ///< set maximal evaluated order
226 
227  bool shellSubA;
229 
230  /**
231  * \brief get options from line command
232  * @return error code
233  */
234  virtual MoFEMErrorCode getOptions();
235 
236  /**
237  * \brief Set IS for levels
238  * @param kk level
239  * @param is pointer to IS
240  * @return error code
241  */
242  virtual MoFEMErrorCode createIsAtLevel(int kk,IS *is);
243 
244  /**
245  * \brief Destroy IS if internally created
246  * @param kk level
247  * @param is pointer to is
248  * @return error code
249  */
250  virtual MoFEMErrorCode destroyIsAtLevel(int kk,IS *is);
251 
252  /**
253  * \brief Set up data structures for MG
254  * @param pc MG pre-conditioner <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/PC/PCMG.html>
255  * @param verb verbosity level
256  * @return error code
257  */
258  virtual MoFEMErrorCode buildProlongationOperator(bool use_mat_a,int verb = 0);
259 
260 
261  // DEPRECATED virtual MoFEMErrorCode buildProlongationOperator(PC pc,int verb = 0) {
262  // return buildProlongationOperator(false,verb);
263  // }
264 
265 };
266 
267 /**
268  * \brief Function build MG structure
269  * @param pc MG pre-conditioner <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/PC/PCMG.html>
270  * @param ctx MoFEM data structure for MG
271  * @param verb verbosity level
272  * @return error code
273  */
275  PC pc,PCMGSetUpViaApproxOrdersCtx *ctx,int verb = 0
276 );
277 
278 
279 #endif //__PCMGSETUP_VIA_APPROX_ORDERS_HPP__
MoFEMErrorCode DMMGViaApproxOrdersSetAO(DM dm, AO ao)
Set DM ordering.
PCMGSetUpViaApproxOrdersCtx(DM dm, Mat a, bool shell_sub_a)
MoFEMErrorCode DMCoarsen_MGViaApproxOrders(DM dm, MPI_Comm comm, DM *dmc)
Coarsen DM.
DM dM
Distributed mesh manager.
MoFEMErrorCode DMMGViaApproxOrdersClearCoarseningIS(DM)
Clear approximation orders.
MoFEMErrorCode DMMGViaApproxOrdersPushBackCoarseningIS(DM, IS is, Mat A, Mat *subA, bool create_sub_matrix, bool shell_sub_a)
Push back coarsening level to MG via approximation orders.
base class for all interface classes
MoFEMErrorCode DMMGViaApproxOrdersGetCoarseningISSize(DM dm, int *size)
Gets size of coarseningIS in internal data struture DMMGViaApproxOrdersCtx.
static const MOFEMuuid IDD_DMMGVIAAPPROXORDERSCTX
boost::ptr_vector< PCMGSubMatrixCtx > shellMatrixCtxPtr
Shell sub-matrix context.
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, MoFEM::UnknownInterface **iface) const
int nbLevels
number of multi-grid levels
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:55
MoFEMErrorCode DMCreateInterpolation_MGViaApproxOrders(DM dm1, DM dm2, Mat *mat, Vec *vec)
Create interpolation matrix between data managers dm1 and dm2.
Structure for DM for multi-grid via approximation orders.
MoFEMErrorCode DMCreateMatrix_MGViaApproxOrders(DM dm, Mat *M)
Create matrix for Multi-Grid via approximation orders.
MoFEMErrorCode PCMGSetUpViaApproxOrders(PC pc, PCMGSetUpViaApproxOrdersCtx *ctx, int verb=0)
Function build MG structure.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
static const int DMMGVIAAPPROXORDERSCTX_INTERFACE
int orderAtLastLevel
set maximal evaluated order
Set data structures of MG pre-conditioner via approximation orders.
virtual MoFEMErrorCode buildProlongationOperator(bool use_mat_a, int verb=0)
Set up data structures for MG.
MoFEMErrorCode DMRegister_MGViaApproxOrders(const char sname[])
Register DM for Multi-Grid via approximation orders.
MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS(DM)
Pop is form MG via approximation orders.
FTensor::Index< 'M', 3 > M
Definition: PlasticOps.hpp:72
int coarseOrder
approximation order of coarse level
virtual MoFEMErrorCode createIsAtLevel(int kk, IS *is)
Set IS for levels.
MoFEMErrorCode DMMGViaApproxOrdersReplaceCoarseningIS(DM dm, IS *is_vec, int nb_elems, Mat A, int verb=0)
Replace coarsening IS in DM via approximation orders.
PETSc Discrete Manager data structure.
Definition: DMMoFEM.hpp:877
virtual MoFEMErrorCode destroyIsAtLevel(int kk, IS *is)
Destroy IS if internally created.
MoFEMErrorCode DMCreate_MGViaApproxOrders(DM dm)
Create DM data structure for Multi-Grid via approximation orders.
std::vector< Mat > kspOperators
Get KSP operators.
std::vector< IS > coarseningIS
Coarsening IS.
virtual MoFEMErrorCode getOptions()
get options from line command
MoFEMErrorCode DMCreateGlobalVector_MGViaApproxOrders(DM dm, Vec *g)
Create global vector for DMGViaApproxOrders.
MoFEMErrorCode DMMGViaApproxOrdersGetCtx(DM dm, DMMGViaApproxOrdersCtx **ctx)