]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Provide an alternative count_dofs_with_subdomain_association function that splits...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 18 Oct 2008 03:10:44 +0000 (03:10 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 18 Oct 2008 03:10:44 +0000 (03:10 +0000)
git-svn-id: https://svn.dealii.org/trunk@17271 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/doc/news/changes.h

index 2634881624ed8cd475f95a2c8ae91cca6804a9d5..e2a3962092397d01b79b14e75099c6678f8d8989 100644 (file)
@@ -1160,6 +1160,28 @@ class DoFTools
     count_dofs_with_subdomain_association (const DH           &dof_handler,
                                            const unsigned int  subdomain);
 
+                                     /**
+                                      * Count how many degrees of freedom are
+                                      * uniquely associated with the given
+                                      * @p subdomain index.
+                                      *
+                                      * This function does what the previous
+                                      * one does except that it splits the
+                                      * result among the vector components of
+                                      * the finite element in use by the
+                                      * DoFHandler object. The last argument
+                                      * (which must have a length equal to the
+                                      * number of vector components) will
+                                      * therefore store how many degrees of
+                                      * freedom of each vector component are
+                                      * associated with the given subdomain.
+                                      */
+    template <class DH>
+    static void
+    count_dofs_with_subdomain_association (const DH           &dof_handler,
+                                           const unsigned int  subdomain,
+                                          std::vector<unsigned int> &n_dofs_on_subdomain);
+    
                                      /**
                                      * Count how many degrees of
                                      * freedom out of the total
index fbbac19936ef31494d6f908c9ebd2a905b72eaeb..559c480adc9baa4bbd579060093b9e5678216989 100644 (file)
@@ -3738,9 +3738,8 @@ DoFTools::get_subdomain_association (const DH                  &dof_handler,
 
 template <class DH>
 unsigned int
-DoFTools::count_dofs_with_subdomain_association (
-  const DH           &dof_handler,
-  const unsigned int  subdomain)
+DoFTools::count_dofs_with_subdomain_association (const DH           &dof_handler,
+                                                const unsigned int  subdomain)
 {
                                    // in debug mode, make sure that there are
                                    // some cells at least with this subdomain
@@ -3771,6 +3770,55 @@ DoFTools::count_dofs_with_subdomain_association (
 
 
 
+template <class DH>
+void
+DoFTools::count_dofs_with_subdomain_association (const DH           &dof_handler,
+                                                const unsigned int  subdomain,
+                                                std::vector<unsigned int> &n_dofs_on_subdomain)
+{
+  Assert (n_dofs_on_subdomain.size() == dof_handler.get_fe().n_components(),
+         ExcDimensionMismatch (n_dofs_on_subdomain.size(),
+                               dof_handler.get_fe().n_components()));
+  std::fill (n_dofs_on_subdomain.begin(), n_dofs_on_subdomain.end(), 0);
+  
+                                   // in debug mode, make sure that there are
+                                   // some cells at least with this subdomain
+                                   // id
+#ifdef DEBUG
+  {
+    bool found = false;
+    for (typename Triangulation<DH::dimension>::active_cell_iterator
+           cell=dof_handler.get_tria().begin_active();
+         cell!=dof_handler.get_tria().end(); ++cell)
+      if (cell->subdomain_id() == subdomain)
+        {
+          found = true;
+          break;
+        }
+    Assert (found == true,
+            ExcMessage ("There are no cells for the given subdomain!"));
+  } 
+#endif
+
+  std::vector<bool> subdomain_association (dof_handler.n_dofs());
+  extract_subdomain_dofs (dof_handler, subdomain, subdomain_association);
+
+  std::vector<bool> component_association (dof_handler.n_dofs());
+  for (unsigned int c=0; c<dof_handler.get_fe().n_components(); ++c)
+    {
+      std::vector<bool> component_mask (dof_handler.get_fe().n_components(), false);
+      component_mask[c] = true;
+      extract_dofs (dof_handler, component_mask, component_association);
+
+      for (unsigned int i=0; i<dof_handler.n_dofs(); ++i)
+       if ((subdomain_association[i] == true) &&
+           (component_association[i] == true))
+         ++n_dofs_on_subdomain[c];
+    }
+}
+
+
+
 template <int dim>
 void
 DoFTools::count_dofs_per_component (
@@ -5788,15 +5836,33 @@ DoFTools::count_dofs_with_subdomain_association<DoFHandler<deal_II_dimension> >
 (const DoFHandler<deal_II_dimension> &,
  const unsigned int);
 template
+void
+DoFTools::count_dofs_with_subdomain_association<DoFHandler<deal_II_dimension> >
+(const DoFHandler<deal_II_dimension> &,
+ const unsigned int,
+ std::vector<unsigned int> &);
+template
 unsigned int
 DoFTools::count_dofs_with_subdomain_association<hp::DoFHandler<deal_II_dimension> >
 (const hp::DoFHandler<deal_II_dimension> &,
  const unsigned int);
 template
+void
+DoFTools::count_dofs_with_subdomain_association<hp::DoFHandler<deal_II_dimension> >
+(const hp::DoFHandler<deal_II_dimension> &,
+ const unsigned int,
+ std::vector<unsigned int> &);
+template
 unsigned int
 DoFTools::count_dofs_with_subdomain_association<MGDoFHandler<deal_II_dimension> >
 (const MGDoFHandler<deal_II_dimension> &,
  const unsigned int);
+template
+void
+DoFTools::count_dofs_with_subdomain_association<MGDoFHandler<deal_II_dimension> >
+(const MGDoFHandler<deal_II_dimension> &,
+ const unsigned int,
+ std::vector<unsigned int> &);
 
 
 template
index 00a7296b5ffc4d726c2f4c4c5b5277ddf04cc197..94005cd21631088bb4374d267cb9b8e9a4c4a7b1 100644 (file)
@@ -353,6 +353,16 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
+  <li>
+  <p>
+  New: There is now a second DoFTools::count_dofs_with_subdomain_association function that
+  calculates the number of degrees of freedom associated with a certain subdomain and
+  splits the result up according to the vector component of each degree of freedom. This
+  function is needed when splitting block matrices in parallel computations.
+  <br>
+  (WB 2008/10/07)
+  </p> 
+
   <li>
   <p>
   Fixed: The GridOut::write_gnuplot function had a bug that made it output only the

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.