Line data Source code
1 : // The libMesh Finite Element Library.
2 : // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 :
4 : // This library is free software; you can redistribute it and/or
5 : // modify it under the terms of the GNU Lesser General Public
6 : // License as published by the Free Software Foundation; either
7 : // version 2.1 of the License, or (at your option) any later version.
8 :
9 : // This library is distributed in the hope that it will be useful,
10 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 : // Lesser General Public License for more details.
13 :
14 : // You should have received a copy of the GNU Lesser General Public
15 : // License along with this library; if not, write to the Free Software
16 : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 :
18 : #ifndef LIBMESH_PETSC_DM_WRAPPER_H
19 : #define LIBMESH_PETSC_DM_WRAPPER_H
20 :
21 : #include "libmesh/petsc_macro.h"
22 :
23 : #ifdef LIBMESH_HAVE_PETSC
24 : #if !PETSC_VERSION_LESS_THAN(3,7,3)
25 : #if defined(LIBMESH_ENABLE_AMR) && defined(LIBMESH_HAVE_METAPHYSICL)
26 :
27 : // libMesh includes
28 : #include "libmesh/petsc_macro.h"
29 : #include "libmesh/petsc_matrix_base.h"
30 : #include "libmesh/petsc_vector.h"
31 : #include "libmesh/wrapped_petsc.h"
32 :
33 : // PETSc includes
34 : #ifdef I
35 : # define LIBMESH_SAW_I
36 : #endif
37 : #include <petsc.h>
38 : #ifndef LIBMESH_SAW_I
39 : # undef I // Avoid complex.h contamination
40 : #endif
41 :
42 : // C++ includes
43 : #include <vector>
44 : #include <memory>
45 : #include <unordered_map>
46 : #include <map>
47 :
48 : namespace libMesh
49 : {
50 : // Forward declarations
51 : class System;
52 : class DofObject;
53 :
54 : //! Struct to house data regarding where in the mesh hierarchy we are located.
55 528 : struct PetscDMContext
56 : {
57 : int n_dofs;
58 : int mesh_dim;
59 : DM * coarser_dm;
60 : DM * finer_dm;
61 : DM * global_dm;
62 : PetscMatrixBase<libMesh::Number> * K_interp_ptr;
63 : PetscMatrixBase<libMesh::Number> * K_sub_interp_ptr;
64 : PetscMatrixBase<libMesh::Number> * K_restrict_ptr;
65 : PetscVector<libMesh::Number> * current_vec;
66 :
67 : //! Stores local dofs for each var for use in subprojection matrixes
68 : std::vector<std::vector<numeric_index_type>> dof_vec;
69 :
70 : //! Stores subfield ids for use in subprojection matrixes on coarser DMs
71 : std::vector<PetscInt> subfields;
72 :
73 560 : PetscDMContext() :
74 528 : n_dofs(-12345),
75 528 : mesh_dim(-12345),
76 528 : coarser_dm(nullptr),
77 528 : finer_dm(nullptr),
78 528 : global_dm(nullptr),
79 528 : K_interp_ptr(nullptr),
80 528 : K_sub_interp_ptr(nullptr),
81 528 : K_restrict_ptr(nullptr),
82 560 : current_vec(nullptr)
83 16 : {}
84 :
85 : };
86 :
87 : /**
88 : * This class defines a wrapper around the PETSc DM infrastructure.
89 : * By coordinating DM data structures with libMesh, we can use libMesh
90 : * mesh hierarchies for geometric multigrid. Additionally, by setting the
91 : * DM data, we can additionally (with or without multigrid) define recursive
92 : * fieldsplits of our variables.
93 : *
94 : * \author Paul T. Bauman, Boris Boutkov
95 : * \date 2018
96 : */
97 216 : class PetscDMWrapper
98 : {
99 : public:
100 :
101 1836 : PetscDMWrapper() = default;
102 :
103 : ~PetscDMWrapper();
104 :
105 : //! Destroys and clears all build DM-related data
106 : void clear();
107 :
108 : void init_and_attach_petscdm(System & system, SNES snes);
109 : void init_and_attach_petscdm(System & system, KSP ksp);
110 :
111 : private:
112 : /**
113 : * Initialize the PETSc DM and return the number of geometric multigrid levels
114 : */
115 : unsigned int init_petscdm(System & system);
116 :
117 : /**
118 : * Vector of DMs for all grid levels. These are PETSc objects
119 : * created by calling DMShellCreate(), so we are responsible for
120 : * cleaning them up.
121 : */
122 : std::vector<WrappedPetsc<DM>> _dms;
123 :
124 : /**
125 : * Vector of PETScSections for all grid levels. These are PETSc
126 : * objects which are attached to the DM by calling DMSetLocalSection().
127 : * While the DM takes care of destroying existing PetscSections in
128 : * calls to DMSetLocalSection(), it does not appear to clean up
129 : * PetscSections it holds when it is destroyed, so we manage their
130 : * lifetimes using WrappedPetsc objects.
131 : */
132 : std::vector<WrappedPetsc<PetscSection>> _sections;
133 :
134 : /**
135 : * Vector of star forests for all grid levels. These are PETSc
136 : * objects which are attached to the DM by calling
137 : * DMSetSectionSF(). The DM seems to take care of cleaning these
138 : * up itself as far as I can tell, so we do not try to manage
139 : * their lifetime in any way.
140 : */
141 : std::vector<PetscSF> _star_forests;
142 :
143 : /**
144 : * Vector of projection matrixes for all grid levels. These are
145 : * C++ objects, they are cleaned up automatically by their
146 : * destructors.
147 : */
148 : std::vector<std::unique_ptr<PetscMatrixBase<Number>>> _pmtx_vec;
149 :
150 : /**
151 : * Vector of sub projection matrixes for all grid levels for
152 : * fieldsplit. These are C++ objects, they are cleaned up
153 : * automatically by their destructors.
154 : */
155 : std::vector<std::unique_ptr<PetscMatrixBase<Number>>> _subpmtx_vec;
156 :
157 : /**
158 : * Vector of internal PetscDM context structs for all grid levels
159 : * Pointers to these C++ objects are passed to DMShellSetContext(),
160 : * they are cleaned up automatically by their destructors.
161 : */
162 : std::vector<PetscDMContext> _ctx_vec;
163 :
164 : /**
165 : * Vector of solution vectors for all grid levels. These are C++
166 : * objects, they are cleaned up automatically by their
167 : * destructors.
168 : */
169 : std::vector<std::unique_ptr<PetscVector<Number>>> _vec_vec;
170 :
171 : //! Stores n_dofs for each grid level, to be used for projection matrix sizing
172 : std::vector<unsigned int> _mesh_dof_sizes;
173 :
174 : //! Stores n_local_dofs for each grid level, to be used for projection vector sizing
175 : std::vector<unsigned int> _mesh_dof_loc_sizes;
176 :
177 : //! Init all the n_mesh_level dependent data structures
178 : void init_dm_data(unsigned int n_levels, const Parallel::Communicator & comm);
179 :
180 : /**
181 : * Get reference to DM for the given mesh level.
182 : * init_dm_data() should be called before this function.
183 : */
184 56 : DM & get_dm(unsigned int level)
185 56 : { libmesh_assert_less(level, _dms.size());
186 1432 : return *_dms[level]; }
187 :
188 : /**
189 : * Get reference to PetscSection for the given mesh level.
190 : * init_dm_data() should be called before this function.
191 : */
192 16 : PetscSection & get_section(unsigned int level)
193 16 : { libmesh_assert_less(level, _sections.size());
194 32 : return *_sections[level]; }
195 :
196 : /**
197 : * Get reference to PetscSF for the given mesh level.
198 : * init_dm_data() should be called before this function.
199 : */
200 16 : PetscSF & get_star_forest(unsigned int level)
201 16 : { libmesh_assert_less(level, _star_forests.size());
202 32 : return _star_forests[level]; }
203 :
204 : /**
205 : * Takes System, empty PetscSection and populates the PetscSection.
206 : * Take the System in its current state and an empty PetscSection and then
207 : * populate the PetscSection. The PetscSection is comprised of global "point"
208 : * numbers, where a "point" in PetscDM parlance is a geometric entity: node, edge,
209 : * face, or element. Then, we also add the DoF numbering for each variable
210 : * for each of the "points". The PetscSection, together the with PetscSF
211 : * will allow for recursive fieldsplits from the command line using PETSc.
212 : */
213 : void build_section(const System & system, PetscSection & section);
214 :
215 : /**
216 : * Takes System, empty PetscSF and populates the PetscSF.
217 : * The PetscSF (star forest) is a cousin of PetscSection. PetscSection
218 : * has the DoF info, and PetscSF gives the parallel distribution of the
219 : * DoF info. So PetscSF should only be necessary when we have more than
220 : * one MPI rank. Essentially, we are copying the DofMap.send_list(): we
221 : * are specifying the local dofs, what rank communicates that dof info
222 : * (for off-processor dofs that are communicated) and the dofs local
223 : * index on that rank.
224 : *
225 : * https://jedbrown.org/files/StarForest.pdf
226 : */
227 : void build_sf( const System & system, PetscSF & star_forest );
228 :
229 : /**
230 : * Helper function for build_section.
231 : * This function will count how many "points" on the current processor have
232 : * DoFs associated with them and give that count to PETSc. We need to cache
233 : * a mapping between the global node id and our local count that we do in this
234 : * function because we will need the local number again in the add_dofs_to_section
235 : * function.
236 : */
237 : void set_point_range_in_section( const System & system,
238 : PetscSection & section,
239 : std::unordered_map<dof_id_type,dof_id_type> & node_map,
240 : std::unordered_map<dof_id_type,dof_id_type> & elem_map,
241 : std::map<dof_id_type,unsigned int> & scalar_map);
242 :
243 : /**
244 : * Helper function for build_section.
245 : * This function will set the DoF info for each "point" in the PetscSection.
246 : */
247 : void add_dofs_to_section (const System & system,
248 : PetscSection & section,
249 : const std::unordered_map<dof_id_type,dof_id_type> & node_map,
250 : const std::unordered_map<dof_id_type,dof_id_type> & elem_map,
251 : const std::map<dof_id_type,unsigned int> & scalar_map);
252 :
253 : /**
254 : * Helper function to sanity check PetscSection construction
255 : * The PetscSection contains local dof information. This helper function just facilitates
256 : * sanity checking that in fact it only has n_local_dofs.
257 : */
258 : dof_id_type check_section_n_dofs( PetscSection & section );
259 :
260 : // Helper function to reduce code duplication when setting dofs in section
261 : void add_dofs_helper (const System & system,
262 : const DofObject & dof_object,
263 : dof_id_type local_id,
264 : PetscSection & section);
265 : };
266 :
267 : }
268 :
269 : #endif // #if LIBMESH_ENABLE_AMR && LIBMESH_HAVE_METAPHYSICL
270 : #endif // #if PETSC_VERSION
271 : #endif // #ifdef LIBMESH_HAVE_PETSC
272 :
273 : #endif // LIBMESH_PETSC_DM_WRAPPER_H
|
