]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Allow AffineConstraints argument to MGTools::make_sparsity_pattern
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 11 Aug 2020 14:47:18 +0000 (16:47 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 11 Aug 2020 14:47:18 +0000 (16:47 +0200)
include/deal.II/multigrid/mg_tools.h
source/multigrid/mg_tools.cc
source/multigrid/mg_tools.inst.in

index b93635f2aaef707401af651b2f88e197baae5279..c5cf619e3eddc4b79254d1880364d3498fe1af6c 100644 (file)
@@ -75,16 +75,25 @@ namespace MGTools
    * Write the sparsity structure of the matrix belonging to the specified @p
    * level. The sparsity pattern is not compressed, so before creating the
    * actual matrix you have to compress the matrix yourself, using
-   * <tt>SparseMatrixStruct::compress()</tt>.
+   * <tt>SparsityPatternType::compress()</tt>.
    *
-   * There is no need to consider hanging nodes here, since only one level is
-   * considered.
+   * The optional AffineConstraints argument allows to define constraints of
+   * the level matrices like Dirichlet boundary conditions. Note that there is
+   * need to consider hanging nodes on the typical level matrices, since only
+   * one level is considered. See DoFTools::make_sparsity_pattern() for more
+   * details about the arguments.
    */
-  template <int dim, int spacedim, typename SparsityPatternType>
+  template <int dim,
+            int spacedim,
+            typename SparsityPatternType,
+            typename number = double>
   void
-  make_sparsity_pattern(const DoFHandler<dim, spacedim> &dof_handler,
-                        SparsityPatternType &            sparsity,
-                        const unsigned int               level);
+  make_sparsity_pattern(
+    const DoFHandler<dim, spacedim> &dof_handler,
+    SparsityPatternType &            sparsity,
+    const unsigned int               level,
+    const AffineConstraints<number> &constraints = AffineConstraints<number>(),
+    const bool                       keep_constrained_dofs = true);
 
   /**
    * Make a sparsity pattern including fluxes of discontinuous Galerkin
index bde170821e826f1dabeba1db1fbb96578c58ce78..a1f65565175749f103203efd53c8577f8924733e 100644 (file)
@@ -567,11 +567,16 @@ namespace MGTools
 
 
 
-  template <int dim, int spacedim, typename SparsityPatternType>
+  template <int dim,
+            int spacedim,
+            typename SparsityPatternType,
+            typename number>
   void
   make_sparsity_pattern(const DoFHandler<dim, spacedim> &dof,
                         SparsityPatternType &            sparsity,
-                        const unsigned int               level)
+                        const unsigned int               level,
+                        const AffineConstraints<number> &constraints,
+                        const bool                       keep_constrained_dofs)
   {
     const types::global_dof_index n_dofs = dof.n_dofs(level);
     (void)n_dofs;
@@ -592,10 +597,9 @@ namespace MGTools
             dof.get_triangulation().locally_owned_subdomain())
         {
           cell->get_mg_dof_indices(dofs_on_this_cell);
-          // make sparsity pattern for this cell
-          for (unsigned int i = 0; i < dofs_per_cell; ++i)
-            for (unsigned int j = 0; j < dofs_per_cell; ++j)
-              sparsity.add(dofs_on_this_cell[i], dofs_on_this_cell[j]);
+          constraints.add_entries_local_to_global(dofs_on_this_cell,
+                                                  sparsity,
+                                                  keep_constrained_dofs);
         }
   }
 
index 28c8c89ffa58cded79964ec5443407bd891792b2..0938499d6e01dbd54ac33e785f6c3614d05a6ca0 100644 (file)
@@ -25,10 +25,24 @@ for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
       template void
       make_sparsity_pattern<deal_II_dimension,
                             deal_II_space_dimension,
-                            PATTERN>(
+                            PATTERN,
+                            double>(
         const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
         PATTERN &,
-        const unsigned int);
+        const unsigned int,
+        const AffineConstraints<double> &,
+        const bool);
+
+      template void
+      make_sparsity_pattern<deal_II_dimension,
+                            deal_II_space_dimension,
+                            PATTERN,
+                            float>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        PATTERN &,
+        const unsigned int,
+        const AffineConstraints<float> &,
+        const bool);
 #endif
     \}
   }

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.