]> 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:12 +0000 (11:39 +0000)
committerjanssen <janssen@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 11 Dec 2011 11:39:12 +0000 (11:39 +0000)
git-svn-id: https://svn.dealii.org/trunk@24820 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/multigrid/mg_constrained_dofs.h

index ed4db82c04df8e300505b560d21ea56b5cd0694b..52a739bbda45d2087478f6e74ecc16f726c20a22 100644 (file)
@@ -74,6 +74,14 @@ class MGConstrainedDoFs : public Subscriptor
     bool is_boundary_index (const unsigned int level,
                            const unsigned int index) const;
 
+                                    /**
+                                     * Determine whether a dof index
+                                     * is at an edge that is not
+                                      * a refinement edge.
+                                     */
+    bool non_refinement_edge_index (const unsigned int level,
+                            const unsigned int index) const;
+
                                     /**
                                      * Determine whether a dof index
                                      * is at the refinement edge.
@@ -98,6 +106,14 @@ class MGConstrainedDoFs : public Subscriptor
     const std::vector<std::set<unsigned int> > &
       get_boundary_indices () const;
 
+                                    /**
+                                     * Return the indices of dofs for each
+                                     * level that lie on the boundary of the
+                                     * domain.
+                                     */
+    const std::vector<std::set<unsigned int> > &
+      get_non_refinement_edge_indices () const;
+
                                     /**
                                      * Return the indices of dofs for each
                                      * level that lie on the refinement edge
@@ -138,6 +154,14 @@ class MGConstrainedDoFs : public Subscriptor
                                      */
     std::vector<std::set<unsigned int> > boundary_indices;
 
+                                    /**
+                                     * The degrees of freedom on egdges
+                                      * that are not a
+                                     * refinement edge between a
+                                     * level and coarser cells.
+                                     */
+    std::vector<std::set<unsigned int> > non_refinement_edge_indices;
+
                                     /**
                                      * The degrees of freedom on the
                                      * refinement edge between a
@@ -166,6 +190,7 @@ MGConstrainedDoFs::initialize(const MGDoFHandler<dim,spacedim>& dof)
   const unsigned int nlevels = dof.get_tria().n_levels();
   refinement_edge_indices.resize(nlevels);
   refinement_edge_boundary_indices.resize(nlevels);
+  non_refinement_edge_indices.resize(nlevels);
   for(unsigned int l=0; l<nlevels; ++l)
     {
       refinement_edge_indices[l].resize(dof.n_dofs(l));
@@ -173,6 +198,8 @@ MGConstrainedDoFs::initialize(const MGDoFHandler<dim,spacedim>& dof)
     }
   MGTools::extract_inner_interface_dofs (dof, refinement_edge_indices,
                                         refinement_edge_boundary_indices);
+
+  MGTools::extract_non_interface_dofs (dof, non_refinement_edge_indices);
 }
 
 
@@ -188,6 +215,7 @@ MGConstrainedDoFs::initialize(
   boundary_indices.resize(nlevels);
   refinement_edge_indices.resize(nlevels);
   refinement_edge_boundary_indices.resize(nlevels);
+  non_refinement_edge_indices.resize(nlevels);
 
   for(unsigned int l=0; l<nlevels; ++l)
     {
@@ -199,6 +227,7 @@ MGConstrainedDoFs::initialize(
   MGTools::make_boundary_list (dof, function_map, boundary_indices, component_mask);
   MGTools::extract_inner_interface_dofs (dof, refinement_edge_indices,
                                         refinement_edge_boundary_indices);
+  MGTools::extract_non_interface_dofs (dof, non_refinement_edge_indices);
 }
 
 
@@ -214,6 +243,9 @@ MGConstrainedDoFs::clear()
 
   for(unsigned int l=0; l<refinement_edge_boundary_indices.size(); ++l)
     refinement_edge_boundary_indices[l].clear();
+
+  for(unsigned int l=0; l<non_refinement_edge_indices.size(); ++l)
+    non_refinement_edge_indices[l].clear();
 }
 
 
@@ -229,6 +261,18 @@ MGConstrainedDoFs::is_boundary_index (const unsigned int level,
     return false;
 }
 
+inline
+bool
+MGConstrainedDoFs::non_refinement_edge_index (const unsigned int level,
+                            const unsigned int index) const
+{
+  AssertIndexRange(level, non_refinement_edge_indices.size());
+
+  if(non_refinement_edge_indices[level].find(index) != non_refinement_edge_indices[level].end())
+    return true;
+  else
+    return false;
+}
 
 inline
 bool
@@ -260,6 +304,13 @@ MGConstrainedDoFs::get_boundary_indices () const
   return boundary_indices;
 }
 
+inline
+const std::vector<std::set<unsigned int> > &
+MGConstrainedDoFs::get_non_refinement_edge_indices () const
+{
+  return non_refinement_edge_indices;
+}
+
 inline
 const std::vector<std::vector<bool> > &
 MGConstrainedDoFs::get_refinement_edge_indices () const

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.