]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Provide the make_sparsity_pattern functions with a subdomain_id argument on which...
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 29 Apr 2009 14:18:10 +0000 (14:18 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 29 Apr 2009 14:18:10 +0000 (14:18 +0000)
git-svn-id: https://svn.dealii.org/trunk@18786 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/dofs/dof_tools.inst.in

index d62ecf4e6f96213867d40cc09682712c0bb42f2b..cdc2086ed137b57ecfd85cf0b6b8e726268055e5 100644 (file)
@@ -428,6 +428,24 @@ class DoFTools
                                      * <tt>keep_constrained_dofs</tt>
                                      * needs to be set to
                                      * <tt>false</tt>.
+                                     *
+                                     * If the @p subdomain_id parameter is
+                                     * given, the sparsity pattern is built
+                                     * only on cells that have a subdomain_id
+                                     * equal to the given argument. This is
+                                     * useful in parallel contexts where the
+                                     * matrix and sparsity pattern (for
+                                     * example a
+                                     * TrilinosWrappers::SparsityPattern) may
+                                     * be distributed and not every MPI
+                                     * process needs to build the entire
+                                     * sparsity pattern; in that case, it is
+                                     * sufficient if every process only
+                                     * builds that part of the sparsity
+                                     * pattern that corresponds to the
+                                     * subdomain_id for which it is
+                                     * responsible. This feature is
+                                     * @ref step_32 "step-32".
                                      */
     template <class DH, class SparsityPattern>
     static
@@ -435,7 +453,8 @@ class DoFTools
     make_sparsity_pattern (const DH               &dof,
                           SparsityPattern        &sparsity_pattern,
                           const ConstraintMatrix &constraints = ConstraintMatrix(),
-                          const bool              keep_constrained_dofs = true);
+                          const bool              keep_constrained_dofs = true,
+                          const unsigned int      subdomain_id = numbers::invalid_unsigned_int);
 
                                     /**
                                      * Locate non-zero entries for
@@ -565,6 +584,24 @@ class DoFTools
                                      * <tt>keep_constrained_dofs</tt>
                                      * needs to be set to
                                      * <tt>false</tt>.
+                                     *
+                                     * If the @p subdomain_id parameter is
+                                     * given, the sparsity pattern is built
+                                     * only on cells that have a subdomain_id
+                                     * equal to the given argument. This is
+                                     * useful in parallel contexts where the
+                                     * matrix and sparsity pattern (for
+                                     * example a
+                                     * TrilinosWrappers::SparsityPattern) may
+                                     * be distributed and not every MPI
+                                     * process needs to build the entire
+                                     * sparsity pattern; in that case, it is
+                                     * sufficient if every process only
+                                     * builds that part of the sparsity
+                                     * pattern that corresponds to the
+                                     * subdomain_id for which it is
+                                     * responsible. This feature is
+                                     * @ref step_32 "step-32".
                                      */
     template <class DH, class SparsityPattern>
     static
@@ -573,7 +610,8 @@ class DoFTools
                           const Table<2, Coupling> &coupling,
                           SparsityPattern          &sparsity_pattern,
                           const ConstraintMatrix   &constraints = ConstraintMatrix(),
-                          const bool                keep_constrained_dofs = true);
+                          const bool                keep_constrained_dofs = true,
+                          const unsigned int        subdomain_id = numbers::invalid_unsigned_int);
     
                                     /**
                                      * @deprecated This is the old
index 3cd9d351e0f240a2772002bcf7fc4088aefc4a9a..d75a185883b7ed10cbb0ada42e2ef938d269c0eb 100644 (file)
@@ -58,7 +58,8 @@ void
 DoFTools::make_sparsity_pattern (const DH               &dof,
                                 SparsityPattern        &sparsity,
                                 const ConstraintMatrix &constraints,
-                                const bool              keep_constrained_dofs)
+                                const bool              keep_constrained_dofs,
+                                const unsigned int      subdomain_id)
 {
   const unsigned int n_dofs = dof.n_dofs();
 
@@ -77,48 +78,35 @@ DoFTools::make_sparsity_pattern (const DH               &dof,
                                   // only have to do the work if the
                                   // current cell is owned by the calling
                                   // processor. Otherwise, just continue.
-                                  // 
-                                  // TODO: here, we should actually check
-                                  // the communicator of the sparsity
-                                  // pattern, not MPI_COMM_WORLD as the
-                                  // Utilities function does!
-  for (; cell!=endc; ++cell) 
-#ifdef DEAL_II_USE_TRILINOS
-    if ((types_are_equal<SparsityPattern,TrilinosWrappers::SparsityPattern>::value
-        ||
-        types_are_equal<SparsityPattern,TrilinosWrappers::BlockSparsityPattern>::value)
-       &&
-       cell->subdomain_id() !=
-       Utilities::Trilinos::get_this_mpi_process(Utilities::Trilinos::comm_world()))
-      continue;
-    else
-        
-#endif
-    {
-      const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell;
-      dofs_on_this_cell.resize (dofs_per_cell);
-      cell->get_dof_indices (dofs_on_this_cell);
-
-                                      // make sparsity pattern for this
-                                      // cell. if no constraints pattern was
-                                      // given, then the following call acts
-                                      // as if simply no constraints existed
-      constraints.add_entries_local_to_global (dofs_on_this_cell,
-                                              sparsity,
-                                              keep_constrained_dofs);
-    }
+  for (; cell!=endc; ++cell)
+    if ((subdomain_id == numbers::invalid_unsigned_int)
+       ||
+       (subdomain_id == cell->subdomain_id()))
+      {
+       const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell;
+       dofs_on_this_cell.resize (dofs_per_cell);
+       cell->get_dof_indices (dofs_on_this_cell);
+
+                                        // make sparsity pattern for this
+                                        // cell. if no constraints pattern was
+                                        // given, then the following call acts
+                                        // as if simply no constraints existed
+       constraints.add_entries_local_to_global (dofs_on_this_cell,
+                                                sparsity,
+                                                keep_constrained_dofs);
+      }
 }
 
 
 
 template <class DH, class SparsityPattern>
 void
-DoFTools::make_sparsity_pattern (
-  const DH                &dof,
-  const Table<2,Coupling> &couplings,
-  SparsityPattern         &sparsity,
-  const ConstraintMatrix  &constraints,
-  const bool               keep_constrained_dofs)
+DoFTools::make_sparsity_pattern (const DH                &dof,
+                                const Table<2,Coupling> &couplings,
+                                SparsityPattern         &sparsity,
+                                const ConstraintMatrix  &constraints,
+                                const bool               keep_constrained_dofs,
+                                const unsigned int       subdomain_id)
 {
   const unsigned int n_dofs = dof.n_dofs();
 
@@ -202,37 +190,29 @@ DoFTools::make_sparsity_pattern (
                                   // current cell is owned by the calling
                                   // processor. Otherwise, just continue.
   for (; cell!=endc; ++cell)
-#ifdef DEAL_II_USE_TRILINOS
-    if ((types_are_equal<SparsityPattern,TrilinosWrappers::SparsityPattern>::value
-        ||
-        types_are_equal<SparsityPattern,TrilinosWrappers::BlockSparsityPattern>::value)
-       &&
-       cell->subdomain_id() != 
-       Utilities::Trilinos::get_this_mpi_process(Utilities::Trilinos::comm_world()))
-      continue;
-    else
-        
-#endif
-    {
-      const unsigned int fe_index
-       = cell->active_fe_index();
+    if ((subdomain_id == numbers::invalid_unsigned_int)
+       ||
+       (subdomain_id == cell->subdomain_id()))
+      {
+       const unsigned int fe_index
+         = cell->active_fe_index();
 
-      const unsigned int dofs_per_cell
-       = fe_collection[fe_index].dofs_per_cell;
+       const unsigned int dofs_per_cell
+         = fe_collection[fe_index].dofs_per_cell;
 
-      dofs_on_this_cell.resize (dofs_per_cell);
-      cell->get_dof_indices (dofs_on_this_cell);
+       dofs_on_this_cell.resize (dofs_per_cell);
+       cell->get_dof_indices (dofs_on_this_cell);
 
 
-                                      // make sparsity pattern for this
-                                      // cell. if no constraints pattern was
-                                      // given, then the following call acts
-                                      // as if simply no constraints existed
-      constraints.add_entries_local_to_global (dofs_on_this_cell,
-                                              sparsity,
-                                              keep_constrained_dofs,
-                                              dof_mask[fe_index]);
-    }
+                                        // make sparsity pattern for this
+                                        // cell. if no constraints pattern was
+                                        // given, then the following call acts
+                                        // as if simply no constraints existed
+       constraints.add_entries_local_to_global (dofs_on_this_cell,
+                                                sparsity,
+                                                keep_constrained_dofs,
+                                                dof_mask[fe_index]);
+      }
 }
 
 
index 8eb98b11e8c23f09d306b6b7fff4ca5369c832bd..efaa23c50b424092926bdc82c6f25cdd40536f7d 100644 (file)
@@ -19,14 +19,16 @@ for (SP : SPARSITY_PATTERNS)
     (const DoFHandler<deal_II_dimension,deal_II_dimension> &dof,
      SP    &sparsity,
      const ConstraintMatrix &,
-     const bool);
+     const bool,
+     const unsigned int);
 
     template void
     DoFTools::make_sparsity_pattern<hp::DoFHandler<deal_II_dimension,deal_II_dimension>, SP>
     (const hp::DoFHandler<deal_II_dimension,deal_II_dimension> &dof,
      SP    &sparsity,
      const ConstraintMatrix &,
-     const bool);
+     const bool,
+     const unsigned int);
 
     template void 
     DoFTools::make_sparsity_pattern<DoFHandler<deal_II_dimension,deal_II_dimension>, SP>
@@ -34,7 +36,8 @@ for (SP : SPARSITY_PATTERNS)
      const Table<2,Coupling>&,
      SP &,
      const ConstraintMatrix &,
-     const bool);
+     const bool,
+     const unsigned int);
 
     template void 
     DoFTools::make_sparsity_pattern<hp::DoFHandler<deal_II_dimension,deal_II_dimension>, SP>
@@ -42,7 +45,8 @@ for (SP : SPARSITY_PATTERNS)
      const Table<2,Coupling>&,
      SP &,
      const ConstraintMatrix &,
-     const bool);
+     const bool,
+     const unsigned int);
 
     template void
     DoFTools::make_sparsity_pattern<DoFHandler<deal_II_dimension,deal_II_dimension>, SP>
@@ -85,14 +89,16 @@ for (SP : SPARSITY_PATTERNS)
     (const DoFHandler<deal_II_dimension,deal_II_dimension+1> &dof,
      SP    &sparsity,
      const ConstraintMatrix &,
-     const bool);
+     const bool,
+     const unsigned int);
 
     template void
     DoFTools::make_sparsity_pattern<hp::DoFHandler<deal_II_dimension,deal_II_dimension+1>, SP>
     (const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> &dof,
      SP    &sparsity,
      const ConstraintMatrix &,
-     const bool);
+     const bool,
+     const unsigned int);
 
     template void 
     DoFTools::make_sparsity_pattern<DoFHandler<deal_II_dimension,deal_II_dimension+1>, SP>
@@ -100,7 +106,8 @@ for (SP : SPARSITY_PATTERNS)
      const Table<2,Coupling>&,
      SP &,
      const ConstraintMatrix &,
-     const bool);
+     const bool,
+     const unsigned int);
 
     template void 
     DoFTools::make_sparsity_pattern<hp::DoFHandler<deal_II_dimension,deal_II_dimension+1>, SP>
@@ -108,7 +115,8 @@ for (SP : SPARSITY_PATTERNS)
      const Table<2,Coupling>&,
      SP &,
      const ConstraintMatrix &,
-     const bool);
+     const bool,
+     const unsigned int);
 
     template void
     DoFTools::make_sparsity_pattern<DoFHandler<deal_II_dimension,deal_II_dimension+1>, SP>

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.