]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
extract also non interface dofs
authorjanssen <janssen@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 11 Dec 2011 11:39:53 +0000 (11:39 +0000)
committerjanssen <janssen@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 11 Dec 2011 11:39:53 +0000 (11:39 +0000)
git-svn-id: https://svn.dealii.org/trunk@24821 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/multigrid/mg_tools.h
deal.II/source/multigrid/mg_tools.cc
deal.II/source/multigrid/mg_tools.inst.in

index 5fb1420f5cbcabb8b4356867b0931f6eacc5fb9a..1734a4ea679a711512330d543807c7f327c887a2 100644 (file)
@@ -297,6 +297,11 @@ namespace MGTools
   void
   extract_inner_interface_dofs (const MGDoFHandler<dim,spacedim> &mg_dof_handler,
                                std::vector<std::vector<bool> >  &interface_dofs);
+
+  template <int dim, int spacedim>
+  void
+  extract_non_interface_dofs (const MGDoFHandler<dim,spacedim> &mg_dof_handler,
+                               std::vector<std::set<unsigned int> >  &non_interface_dofs);
 }
 
 /* @} */
index f021ab74b177e867ad5dca77b4e048017c508247..76b6fc0085996facb93a80d4dd049b2ca594bd96 100644 (file)
@@ -1460,6 +1460,72 @@ namespace MGTools
   }
 
 
+template <int dim, int spacedim>
+void
+extract_non_interface_dofs (const MGDoFHandler<dim,spacedim> &mg_dof_handler,
+                           std::vector<std::set<unsigned int> >  &non_interface_dofs)
+{
+  Assert (non_interface_dofs.size() == mg_dof_handler.get_tria().n_levels(),
+         ExcDimensionMismatch (non_interface_dofs.size(),
+                               mg_dof_handler.get_tria().n_levels()));
+
+  const unsigned int nlevels = mg_dof_handler.get_tria().n_levels();
+  const FiniteElement<dim,spacedim> &fe = mg_dof_handler.get_fe();
+
+  const unsigned int   dofs_per_cell   = fe.dofs_per_cell;
+  const unsigned int   dofs_per_face   = fe.dofs_per_face;
+
+  std::vector<unsigned int> local_dof_indices (dofs_per_cell);
+  std::vector<bool> cell_dofs(dofs_per_cell, false);
+  std::vector<bool> cell_dofs_interface(dofs_per_cell, false);
+
+  typename MGDoFHandler<dim>::cell_iterator cell = mg_dof_handler.begin(),
+                                           endc = mg_dof_handler.end();
+
+
+  bool any_dof_set = false;
+  for (; cell!=endc; ++cell)
+    {
+      std::fill (cell_dofs.begin(), cell_dofs.end(), false);
+      std::fill (cell_dofs_interface.begin(), cell_dofs_interface.end(), false);
+
+      for (unsigned int face_nr=0; face_nr<GeometryInfo<dim>::faces_per_cell; ++face_nr)
+       {
+         const typename DoFHandler<dim,spacedim>::face_iterator face = cell->face(face_nr);
+         if (!face->at_boundary())
+           {
+                                              //interior face
+             const typename MGDoFHandler<dim>::cell_iterator
+               neighbor = cell->neighbor(face_nr);
+
+             if ((neighbor->level() < cell->level()))
+               {
+                 for (unsigned int j=0; j<dofs_per_face; ++j)
+                   cell_dofs_interface[fe.face_to_cell_index(j,face_nr)] = true;
+               }
+              else
+              {
+                for (unsigned int j=0; j<dofs_per_face; ++j)
+                   cell_dofs[fe.face_to_cell_index(j,face_nr)] = true;
+              }
+           }
+          else
+          {
+                                              //boundary face
+                for (unsigned int j=0; j<dofs_per_face; ++j)
+                   cell_dofs[fe.face_to_cell_index(j,face_nr)] = true;
+          }
+       }
+
+      const unsigned int level = cell->level();
+      cell->get_mg_dof_indices (local_dof_indices);
+
+      for(unsigned int i=0; i<dofs_per_cell; ++i)
+         if (cell_dofs[i] && !cell_dofs_interface[i])
+           non_interface_dofs[level].insert(local_dof_indices[i]);
+    }
+}
+
 
   template <int dim, int spacedim>
   void
index 5d4b07132b604958ddd72ca8d162a61ad8a50183..18b01e595850f5d713dd2506c1cbf5ea8f6e0ced 100644 (file)
@@ -106,6 +106,10 @@ for (deal_II_dimension : DIMENSIONS)
        void
        extract_inner_interface_dofs (const MGDoFHandler<deal_II_dimension> &mg_dof_handler,
                                      std::vector<std::vector<bool> >  &interface_dofs);
+      template
+        void
+        extract_non_interface_dofs (const MGDoFHandler<deal_II_dimension> & mg_dof_handler,
+                                    std::vector<std::set<unsigned int> > &non_interface_dofs);
 
 #if deal_II_dimension < 3
       template void count_dofs_per_block<deal_II_dimension,deal_II_dimension+1> (

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.