]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
multigrid adapted to level-less faces and edges
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 6 Jul 2006 14:58:35 +0000 (14:58 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 6 Jul 2006 14:58:35 +0000 (14:58 +0000)
git-svn-id: https://svn.dealii.org/trunk@13348 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/mg_dof_accessor.h
deal.II/deal.II/source/multigrid/mg_dof_accessor.cc
deal.II/deal.II/source/multigrid/mg_dof_handler.cc
deal.II/deal.II/source/multigrid/mg_dof_tools.cc
deal.II/deal.II/source/multigrid/mg_smoother.cc

index 6b52bcc4ac58c2d77b87be95a071cc0d96eabfa3..5c6427ef11e73db04bfb05f17f5ff6da1b9699a8 100644 (file)
@@ -153,7 +153,8 @@ class MGDoFAccessor : public MGDoFObjectAccessor_Inheritance<structdim, dim>::Ba
                                      * vertex for the level this
                                      * object lives on.
                                      */
-    unsigned int mg_vertex_dof_index (const unsigned int vertex,
+    unsigned int mg_vertex_dof_index (const int level,
+                                     const unsigned int vertex,
                                      const unsigned int i) const;
 
                                     /**
@@ -161,7 +162,8 @@ class MGDoFAccessor : public MGDoFObjectAccessor_Inheritance<structdim, dim>::Ba
                                      * on the @p vertexth vertex to @p index
                                      * for the level this object lives on.
                                      */
-    void set_mg_vertex_dof_index (const unsigned int vertex,
+    void set_mg_vertex_dof_index (const int level,
+                                 const unsigned int vertex,
                                  const unsigned int i,
                                  const unsigned int index) const;
 
@@ -171,14 +173,16 @@ class MGDoFAccessor : public MGDoFObjectAccessor_Inheritance<structdim, dim>::Ba
                                      * on the level this line lives
                                      * on.
                                      */
-    unsigned int mg_dof_index (const unsigned int i) const;
+    unsigned int mg_dof_index (const int level,
+                              const unsigned int i) const;
 
                                     /**
                                      * Set the index of the @p ith degree
                                      * of freedom of this line on the
                                      * level this line lives on to @p index.
                                      */
-    void set_mg_dof_index (const unsigned int i,
+    void set_mg_dof_index (const int level,
+                          const unsigned int i,
                           const unsigned int index) const;
 
                                     /**
@@ -287,7 +291,8 @@ class MGDoFObjectAccessor<1, dim> :  public MGDoFAccessor<1,dim>
                                      * indices refer to the local numbering
                                      * for the level this line lives on.
                                      */
-    void get_mg_dof_indices (std::vector<unsigned int> &dof_indices) const;
+    void get_mg_dof_indices (const int level,
+                            std::vector<unsigned int> &dof_indices) const;
 
                                     /**
                                      * Return the value of the given vector
@@ -307,7 +312,8 @@ class MGDoFObjectAccessor<1, dim> :  public MGDoFAccessor<1,dim>
                                      * freedom on this level.
                                      */
     template <typename number>
-    void get_mg_dof_values (const Vector<number> &values,
+    void get_mg_dof_values (const int level,
+                           const Vector<number> &values,
                            Vector<number>       &dof_values) const;
 };
 
@@ -370,7 +376,8 @@ class MGDoFObjectAccessor<2, dim> :  public MGDoFAccessor<2,dim>
                                      * indices refer to the local numbering
                                      * for the level this quad lives on.
                                      */
-    void get_mg_dof_indices (std::vector<unsigned int> &dof_indices) const;
+    void get_mg_dof_indices (const int level,
+                            std::vector<unsigned int> &dof_indices) const;
  
                                     /**
                                      * Return the value of the given vector
@@ -390,7 +397,8 @@ class MGDoFObjectAccessor<2, dim> :  public MGDoFAccessor<2,dim>
                                      * freedom on this level.
                                      */
     template <typename number>
-    void get_mg_dof_values (const Vector<number> &values,
+    void get_mg_dof_values (const int level,
+                           const Vector<number> &values,
                            Vector<number>       &dof_values) const;
 
                                     /**
@@ -461,7 +469,8 @@ class MGDoFObjectAccessor<3, dim> :  public MGDoFAccessor<3,dim>
                                      * indices refer to the local numbering
                                      * for the level this hex lives on.
                                      */
-    void get_mg_dof_indices (std::vector<unsigned int> &dof_indices) const;
+    void get_mg_dof_indices (const int level,
+                            std::vector<unsigned int> &dof_indices) const;
 
                                     /**
                                      * Return the value of the given vector
@@ -481,7 +490,8 @@ class MGDoFObjectAccessor<3, dim> :  public MGDoFAccessor<3,dim>
                                      * freedom on this level.
                                      */
     template <typename number>
-    void get_mg_dof_values (const Vector<number> &values,
+    void get_mg_dof_values (const int level,
+                           const Vector<number> &values,
                            Vector<number>       &dof_values) const;
 
                                     /**
@@ -540,6 +550,41 @@ class MGDoFCellAccessor :  public MGDoFObjectAccessor<dim, dim>
                       const int                 index,
                       const AccessorData       *local_data);
 
+                                    /**
+                                     * Return the indices of the dofs of this
+                                     * hex in the standard ordering: dofs
+                                     * on vertex 0, dofs on vertex 1, etc,
+                                     * dofs on line 0, dofs on line 1, etc,
+                                     * dofs on quad 0, etc.
+                                     *
+                                     * It is assumed that the vector already
+                                     * has the right size beforehand. The
+                                     * indices refer to the local numbering
+                                     * for the level this hex lives on.
+                                     */
+    void get_mg_dof_indices (std::vector<unsigned int> &dof_indices) const;
+
+                                    /**
+                                     * Return the value of the given vector
+                                     * restricted to the dofs of this
+                                     * cell in the standard ordering: dofs
+                                     * on vertex 0, dofs on vertex 1, etc,
+                                     * dofs on line 0, dofs on line 1, etc,
+                                     * dofs on quad 0, etc.
+                                     *
+                                     * It is assumed that the vector already
+                                     * has the right size beforehand. The
+                                     * indices refer to the multilevel
+                                     * numbering local to the present
+                                     * level of this cell. The vector shall
+                                     * therefore have the same number of
+                                     * entries as there are degrees of
+                                     * freedom on this level.
+                                     */
+    template <typename number>
+    void get_mg_dof_values (const Vector<number> &values,
+                           Vector<number>       &dof_values) const;
+
                                     /**
                                      * Return the @p ith neighbor as a MGDoF cell
                                      * iterator. This function is needed since
index 943ab81ee90fb4e6ef60557ed083136e004ec46a..6500466da94f8ea6e10a997195bf5d45e5103bdf 100644 (file)
@@ -80,7 +80,8 @@ MGDoFAccessor<structdim, dim>::operator = (const MGDoFAccessor &da)
 
 template <int structdim,int dim>
 unsigned int
-MGDoFAccessor<structdim, dim>::mg_vertex_dof_index (const unsigned int vertex,
+MGDoFAccessor<structdim, dim>::mg_vertex_dof_index (const int level,
+                                                   const unsigned int vertex,
                                                    const unsigned int i) const
 {
   Assert (this->dof_handler != 0,
@@ -95,13 +96,14 @@ MGDoFAccessor<structdim, dim>::mg_vertex_dof_index (const unsigned int vertex,
          ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_vertex));
 
   return (this->mg_dof_handler->mg_vertex_dofs[this->vertex_index(vertex)]
-         .get_index (this->present_level, i, this->dof_handler->get_fe().dofs_per_vertex));
+         .get_index (level, i, this->dof_handler->get_fe().dofs_per_vertex));
 }
 
 
 template <int structdim, int dim>
 void
-MGDoFAccessor<structdim, dim>::set_mg_vertex_dof_index (const unsigned int vertex,
+MGDoFAccessor<structdim, dim>::set_mg_vertex_dof_index (const int level,
+                                                       const unsigned int vertex,
                                                        const unsigned int i,
                                                        const unsigned int index) const
 {
@@ -117,17 +119,18 @@ MGDoFAccessor<structdim, dim>::set_mg_vertex_dof_index (const unsigned int verte
          ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_vertex));
 
   this->mg_dof_handler->mg_vertex_dofs[this->vertex_index(vertex)]
-    .set_index (this->present_level, i, this->dof_handler->get_fe().dofs_per_vertex, index);
+    .set_index (level, i, this->dof_handler->get_fe().dofs_per_vertex, index);
 }
 
 
 
 template <int structdim, int dim>
 unsigned int
-MGDoFAccessor<structdim, dim>::mg_dof_index (const unsigned int i) const
+MGDoFAccessor<structdim, dim>::mg_dof_index (const int level,
+                                            const unsigned int i) const
 {
   return this->mg_dof_handler
-    ->template get_dof_index<structdim>(this->present_level,
+    ->template get_dof_index<structdim>(level,
                                        this->present_index,
                                        0,
                                        i);
@@ -135,11 +138,12 @@ MGDoFAccessor<structdim, dim>::mg_dof_index (const unsigned int i) const
 
 
 template <int structdim, int dim>
-void MGDoFAccessor<structdim, dim>::set_mg_dof_index (const unsigned int i,
-                                                   const unsigned int index) const
+void MGDoFAccessor<structdim, dim>::set_mg_dof_index (const int level,
+                                                     const unsigned int i,
+                                                     const unsigned int index) const
 {
   this->mg_dof_handler
-    ->template set_dof_index<structdim>(this->present_level,
+    ->template set_dof_index<structdim>(level,
                                        this->present_index,
                                        0,
                                        i,
@@ -215,7 +219,8 @@ MGDoFObjectAccessor<1, dim>::MGDoFObjectAccessor (const Triangulation<dim> *tria
 
 template <int dim>
 void
-MGDoFObjectAccessor<1, dim>::get_mg_dof_indices (std::vector<unsigned int> &dof_indices) const
+MGDoFObjectAccessor<1, dim>::get_mg_dof_indices (const int level,
+                                                std::vector<unsigned int> &dof_indices) const
 {
   Assert (this->dof_handler != 0,
          typename BaseClass::ExcInvalidObject());
@@ -232,9 +237,9 @@ MGDoFObjectAccessor<1, dim>::get_mg_dof_indices (std::vector<unsigned int> &dof_
   std::vector<unsigned int>::iterator next = dof_indices.begin();
   for (unsigned int vertex=0; vertex<2; ++vertex)
     for (unsigned int d=0; d<dofs_per_vertex; ++d)
-      *next++ = this->mg_vertex_dof_index(vertex,d);
+      *next++ = this->mg_vertex_dof_index(level,vertex,d);
   for (unsigned int d=0; d<dofs_per_line; ++d)
-    *next++ = this->mg_dof_index(d);
+    *next++ = this->mg_dof_index(level,d);
   
   Assert (next == dof_indices.end(),
          ExcInternalError());
@@ -244,7 +249,8 @@ MGDoFObjectAccessor<1, dim>::get_mg_dof_indices (std::vector<unsigned int> &dof_
 template <int dim>
 template <typename number>
 void
-MGDoFObjectAccessor<1,dim>::get_mg_dof_values (const Vector<number> &values,
+MGDoFObjectAccessor<1,dim>::get_mg_dof_values (const int level,
+                                              const Vector<number> &values,
                                               Vector<number>       &dof_values) const
 {
   Assert (this->dof_handler != 0,
@@ -263,9 +269,9 @@ MGDoFObjectAccessor<1,dim>::get_mg_dof_values (const Vector<number> &values,
   typename Vector<number>::iterator next_dof_value=dof_values.begin();
   for (unsigned int vertex=0; vertex<2; ++vertex)
     for (unsigned int d=0; d<dofs_per_vertex; ++d)
-      *next_dof_value++ = values(this->mg_vertex_dof_index(vertex,d));
+      *next_dof_value++ = values(this->mg_vertex_dof_index(level,vertex,d));
   for (unsigned int d=0; d<dofs_per_line; ++d)
-    *next_dof_value++ = values(this->mg_dof_index(d));
+    *next_dof_value++ = values(this->mg_dof_index(level,d));
   
   Assert (next_dof_value == dof_values.end(),
          ExcInternalError());
@@ -289,7 +295,8 @@ MGDoFObjectAccessor<2, dim>::MGDoFObjectAccessor (const Triangulation<dim> *tria
 
 template <int dim>
 void
-MGDoFObjectAccessor<2, dim>::get_mg_dof_indices (std::vector<unsigned int> &dof_indices) const
+MGDoFObjectAccessor<2, dim>::get_mg_dof_indices (const int level,
+                                                std::vector<unsigned int> &dof_indices) const
 {
   Assert (this->dof_handler != 0,
          typename BaseClass::ExcInvalidObject());
@@ -308,12 +315,12 @@ MGDoFObjectAccessor<2, dim>::get_mg_dof_indices (std::vector<unsigned int> &dof_
   std::vector<unsigned int>::iterator next = dof_indices.begin();
   for (unsigned int vertex=0; vertex<4; ++vertex)
     for (unsigned int d=0; d<dofs_per_vertex; ++d)
-      *next++ = this->mg_vertex_dof_index(vertex,d);
+      *next++ = this->mg_vertex_dof_index(level,vertex,d);
   for (unsigned int line=0; line<4; ++line)
     for (unsigned int d=0; d<dofs_per_line; ++d)
-      *next++ = this->line(line)->mg_dof_index(d);
+      *next++ = this->line(line)->mg_dof_index(level,d);
   for (unsigned int d=0; d<dofs_per_quad; ++d)
-    *next++ = this->mg_dof_index(d);
+    *next++ = this->mg_dof_index(level,d);
   
   Assert (next == dof_indices.end(),
          ExcInternalError());
@@ -323,7 +330,8 @@ MGDoFObjectAccessor<2, dim>::get_mg_dof_indices (std::vector<unsigned int> &dof_
 template <int dim>
 template <typename number>
 void
-MGDoFObjectAccessor<2,dim>::get_mg_dof_values (const Vector<number> &values,
+MGDoFObjectAccessor<2,dim>::get_mg_dof_values (const int level,
+                                              const Vector<number> &values,
                                               Vector<number>       &dof_values) const
 {
   Assert (this->dof_handler != 0,
@@ -334,7 +342,7 @@ MGDoFObjectAccessor<2,dim>::get_mg_dof_values (const Vector<number> &values,
          typename BaseClass::ExcInvalidObject());
   Assert (dof_values.size() == this->dof_handler->get_fe().dofs_per_cell,
          typename BaseClass::ExcVectorDoesNotMatch());
-  Assert (values.size() == this->mg_dof_handler->n_dofs(this->present_level),
+  Assert (values.size() == this->mg_dof_handler->n_dofs(level),
          typename BaseClass::ExcVectorDoesNotMatch());
 
   const unsigned int dofs_per_vertex = this->dof_handler->get_fe().dofs_per_vertex,
@@ -343,12 +351,12 @@ MGDoFObjectAccessor<2,dim>::get_mg_dof_values (const Vector<number> &values,
   typename Vector<number>::iterator next_dof_value=dof_values.begin();
   for (unsigned int vertex=0; vertex<4; ++vertex)
     for (unsigned int d=0; d<dofs_per_vertex; ++d)
-      *next_dof_value++ = values(this->mg_vertex_dof_index(vertex,d));
+      *next_dof_value++ = values(this->mg_vertex_dof_index(level,vertex,d));
   for (unsigned int line=0; line<4; ++line)
     for (unsigned int d=0; d<dofs_per_line; ++d)
-      *next_dof_value++ = values(this->line(line)->mg_dof_index(d));
+      *next_dof_value++ = values(this->line(line)->mg_dof_index(level,d));
   for (unsigned int d=0; d<dofs_per_quad; ++d)
-    *next_dof_value++ = values(this->mg_dof_index(d));
+    *next_dof_value++ = values(this->mg_dof_index(level,d));
   
   Assert (next_dof_value == dof_values.end(),
          ExcInternalError());
@@ -387,7 +395,8 @@ MGDoFObjectAccessor<3, dim>::MGDoFObjectAccessor (const Triangulation<dim> *tria
 
 template <int dim>
 void
-MGDoFObjectAccessor<3, dim>::get_mg_dof_indices (std::vector<unsigned int> &dof_indices) const
+MGDoFObjectAccessor<3, dim>::get_mg_dof_indices (const int level,
+                                                std::vector<unsigned int> &dof_indices) const
 {
   Assert (this->dof_handler != 0,
          typename BaseClass::ExcInvalidObject());
@@ -408,15 +417,15 @@ MGDoFObjectAccessor<3, dim>::get_mg_dof_indices (std::vector<unsigned int> &dof_
   std::vector<unsigned int>::iterator next = dof_indices.begin();
   for (unsigned int vertex=0; vertex<8; ++vertex)
     for (unsigned int d=0; d<dofs_per_vertex; ++d)
-      *next++ = this->mg_vertex_dof_index(vertex,d);
+      *next++ = this->mg_vertex_dof_index(level,vertex,d);
   for (unsigned int line=0; line<12; ++line)
     for (unsigned int d=0; d<dofs_per_line; ++d)
-      *next++ = this->line(line)->mg_dof_index(d);
+      *next++ = this->line(line)->mg_dof_index(level,d);
   for (unsigned int quad=0; quad<6; ++quad)
     for (unsigned int d=0; d<dofs_per_quad; ++d)
-      *next++ = this->quad(quad)->mg_dof_index(d);
+      *next++ = this->quad(quad)->mg_dof_index(level,d);
   for (unsigned int d=0; d<dofs_per_hex; ++d)
-    *next++ = this->mg_dof_index(d);
+    *next++ = this->mg_dof_index(level,d);
   
   Assert (next == dof_indices.end(),
          ExcInternalError());
@@ -426,7 +435,8 @@ MGDoFObjectAccessor<3, dim>::get_mg_dof_indices (std::vector<unsigned int> &dof_
 template <int dim>
 template <typename number>
 void
-MGDoFObjectAccessor<3,dim>::get_mg_dof_values (const Vector<number> &values,
+MGDoFObjectAccessor<3,dim>::get_mg_dof_values (const int level,
+                                              const Vector<number> &values,
                                               Vector<number>       &dof_values) const
 {
   Assert (this->dof_handler != 0,
@@ -437,7 +447,7 @@ MGDoFObjectAccessor<3,dim>::get_mg_dof_values (const Vector<number> &values,
          typename BaseClass::ExcInvalidObject());
   Assert (dof_values.size() == this->dof_handler->get_fe().dofs_per_cell,
          typename BaseClass::ExcVectorDoesNotMatch());
-  Assert (values.size() == this->mg_dof_handler->n_dofs(this->present_level),
+  Assert (values.size() == this->mg_dof_handler->n_dofs(level),
          typename BaseClass::ExcVectorDoesNotMatch());
 
   const unsigned int dofs_per_vertex = this->dof_handler->get_fe().dofs_per_vertex,
@@ -447,13 +457,13 @@ MGDoFObjectAccessor<3,dim>::get_mg_dof_values (const Vector<number> &values,
   typename Vector<number>::iterator next_dof_value=dof_values.begin();
   for (unsigned int vertex=0; vertex<8; ++vertex)
     for (unsigned int d=0; d<dofs_per_vertex; ++d)
-      *next_dof_value++ = values(this->mg_vertex_dof_index(vertex,d));
+      *next_dof_value++ = values(this->mg_vertex_dof_index(level,vertex,d));
   for (unsigned int line=0; line<12; ++line)
     for (unsigned int d=0; d<dofs_per_line; ++d)
-      *next_dof_value++ = values(this->line(line)->mg_dof_index(d));
+      *next_dof_value++ = values(this->line(line)->mg_dof_index(level,d));
   for (unsigned int quad=0; quad<6; ++quad)
     for (unsigned int d=0; d<dofs_per_quad; ++d)
-      *next_dof_value++ = values(this->quad(quad)->mg_dof_index(d));
+      *next_dof_value++ = values(this->quad(quad)->mg_dof_index(level,d));
   for (unsigned int d=0; d<dofs_per_hex; ++d)
     *next_dof_value++ = values(this->dof_index(d));
   
@@ -507,6 +517,25 @@ MGDoFCellAccessor<dim>::MGDoFCellAccessor (const Triangulation<dim> *tria,
 {}
 
 
+template <int dim>
+void
+MGDoFCellAccessor<dim>::get_mg_dof_indices (std::vector<unsigned int> &dof_indices) const
+{
+  MGDoFObjectAccessor<dim, dim>::get_mg_dof_indices (this->level(), dof_indices);
+}
+
+
+template <int dim>
+template <typename number>
+void
+MGDoFCellAccessor<dim>::get_mg_dof_values (const Vector<number> &values,
+                                          Vector<number>       &dof_values) const
+{
+  MGDoFObjectAccessor<dim, dim>::get_mg_dof_values (this->level(), values, dof_values);
+}
+
+
+
 template <int dim>
 TriaIterator<dim,MGDoFCellAccessor<dim> >
 MGDoFCellAccessor<dim>::neighbor (const unsigned int i) const
@@ -601,13 +630,15 @@ neighbor_child_on_subface (const unsigned int face,
 template
 void
 MGDoFObjectAccessor<1,deal_II_dimension>::
-get_mg_dof_values (const Vector<double> &values,
+get_mg_dof_values (const int level,
+                  const Vector<double> &values,
                   Vector<double>       &dof_values) const;
 
 template
 void
 MGDoFObjectAccessor<1,deal_II_dimension>::
-get_mg_dof_values (const Vector<float> &values,
+get_mg_dof_values (const int level,
+                  const Vector<float> &values,
                   Vector<float>       &dof_values) const;
 
 
@@ -616,13 +647,15 @@ get_mg_dof_values (const Vector<float> &values,
 template
 void
 MGDoFObjectAccessor<2,deal_II_dimension>::
-get_mg_dof_values (const Vector<double> &values,
+get_mg_dof_values (const int level,
+                  const Vector<double> &values,
                   Vector<double>       &dof_values) const;
 
 template
 void
 MGDoFObjectAccessor<2,deal_II_dimension>::
-get_mg_dof_values (const Vector<float> &values,
+get_mg_dof_values (const int level,
+                  const Vector<float> &values,
                   Vector<float>       &dof_values) const;
 
 #endif
@@ -633,16 +666,30 @@ get_mg_dof_values (const Vector<float> &values,
 template
 void
 MGDoFObjectAccessor<3,deal_II_dimension>::
-get_mg_dof_values (const Vector<double> &values,
+get_mg_dof_values (const int level,
+                  const Vector<double> &values,
                   Vector<double>       &dof_values) const;
 
 template
 void
 MGDoFObjectAccessor<3,deal_II_dimension>::
-get_mg_dof_values (const Vector<float> &values,
+get_mg_dof_values (const int level,
+                  const Vector<float> &values,
                   Vector<float>       &dof_values) const;
 
 #endif
+template
+void
+MGDoFCellAccessor<deal_II_dimension>::
+get_mg_dof_values (const Vector<double> &values,
+                  Vector<double>       &dof_values) const;
+
+template
+void
+MGDoFCellAccessor<deal_II_dimension>::
+get_mg_dof_values (const Vector<float> &values,
+                  Vector<float>       &dof_values) const;
+
 
 
 #if deal_II_dimension == 1
index bbd5a8b1c24ed517429408e37f57ae8fa1007899..d4c699166cab45634365bcf01423d416b07985e3 100644 (file)
@@ -1733,12 +1733,12 @@ MGDoFHandler<1>::distribute_dofs_on_cell (cell_iterator &cell,
              {
                if (v==0) 
                  for (unsigned int d=0; d<this->selected_fe->dofs_per_vertex; ++d)
-                   cell->set_mg_vertex_dof_index (0, d,
-                                                  neighbor->mg_vertex_dof_index (1, d));
+                   cell->set_mg_vertex_dof_index (cell->level(), 0, d,
+                                                  neighbor->mg_vertex_dof_index (cell->level(), 1, d));
                else
                  for (unsigned int d=0; d<this->selected_fe->dofs_per_vertex; ++d)
-                   cell->set_mg_vertex_dof_index (1, d,
-                                                  neighbor->mg_vertex_dof_index (0, d));
+                   cell->set_mg_vertex_dof_index (cell->level(), 1, d,
+                                                  neighbor->mg_vertex_dof_index (cell->level(), 0, d));
                
                                                 // next neighbor
                continue;
@@ -1747,13 +1747,13 @@ MGDoFHandler<1>::distribute_dofs_on_cell (cell_iterator &cell,
        
                                         // otherwise: create dofs newly
        for (unsigned int d=0; d<this->selected_fe->dofs_per_vertex; ++d)
-         cell->set_mg_vertex_dof_index (v, d, next_free_dof++);
+         cell->set_mg_vertex_dof_index (cell->level(), v, d, next_free_dof++);
       };
   
                                   // dofs of line
   if (this->selected_fe->dofs_per_line > 0)
     for (unsigned int d=0; d<this->selected_fe->dofs_per_line; ++d)
-      cell->set_mg_dof_index (d, next_free_dof++);
+      cell->set_mg_dof_index (cell->level(), d, next_free_dof++);
 
                                   // note that this cell has been processed
   cell->set_user_flag ();
@@ -1776,9 +1776,9 @@ MGDoFHandler<2>::distribute_dofs_on_cell (cell_iterator &cell,
                                       // check whether dofs for this
                                       // vertex have been distributed
                                       // (only check the first dof)
-      if (cell->mg_vertex_dof_index(vertex, 0) == DoFHandler<2>::invalid_dof_index)
+      if (cell->mg_vertex_dof_index(cell->level(), vertex, 0) == DoFHandler<2>::invalid_dof_index)
        for (unsigned int d=0; d<this->selected_fe->dofs_per_vertex; ++d)
-         cell->set_mg_vertex_dof_index (vertex, d, next_free_dof++);
+         cell->set_mg_vertex_dof_index (cell->level(), vertex, d, next_free_dof++);
     
                                   // for the four sides
   if (this->selected_fe->dofs_per_line > 0)
@@ -1789,17 +1789,17 @@ MGDoFHandler<2>::distribute_dofs_on_cell (cell_iterator &cell,
                                         // distribute dofs if necessary:
                                         // check whether line dof is already
                                         // numbered (check only first dof)
-       if (line->mg_dof_index(0) == DoFHandler<2>::invalid_dof_index)
+       if (line->mg_dof_index(cell->level(), 0) == DoFHandler<2>::invalid_dof_index)
                                           // if not: distribute dofs
          for (unsigned int d=0; d<this->selected_fe->dofs_per_line; ++d)
-           line->set_mg_dof_index (d, next_free_dof++);            
+           line->set_mg_dof_index (cell->level(), d, next_free_dof++);     
       };
 
 
                                   // dofs of quad
   if (this->selected_fe->dofs_per_quad > 0)
     for (unsigned int d=0; d<this->selected_fe->dofs_per_quad; ++d)
-      cell->set_mg_dof_index (d, next_free_dof++);
+      cell->set_mg_dof_index (cell->level(), d, next_free_dof++);
 
 
                                   // note that this cell has been processed
@@ -1823,9 +1823,9 @@ MGDoFHandler<3>::distribute_dofs_on_cell (cell_iterator &cell,
                                       // check whether dofs for this
                                       // vertex have been distributed
                                       // (only check the first dof)
-      if (cell->mg_vertex_dof_index(vertex, 0) == DoFHandler<3>::invalid_dof_index)
+      if (cell->mg_vertex_dof_index(cell->level(), vertex, 0) == DoFHandler<3>::invalid_dof_index)
        for (unsigned int d=0; d<this->selected_fe->dofs_per_vertex; ++d)
-         cell->set_mg_vertex_dof_index (vertex, d, next_free_dof++);
+         cell->set_mg_vertex_dof_index (cell->level(), vertex, d, next_free_dof++);
     
                                   // for the lines
   if (this->selected_fe->dofs_per_line > 0)
@@ -1836,10 +1836,10 @@ MGDoFHandler<3>::distribute_dofs_on_cell (cell_iterator &cell,
                                         // distribute dofs if necessary:
                                         // check whether line dof is already
                                         // numbered (check only first dof)
-       if (line->mg_dof_index(0) == DoFHandler<3>::invalid_dof_index)
+       if (line->mg_dof_index(cell->level(), 0) == DoFHandler<3>::invalid_dof_index)
                                           // if not: distribute dofs
          for (unsigned int d=0; d<this->selected_fe->dofs_per_line; ++d)
-           line->set_mg_dof_index (d, next_free_dof++);            
+           line->set_mg_dof_index (cell->level(), d, next_free_dof++);     
       };
 
                                   // for the quads
@@ -1851,17 +1851,17 @@ MGDoFHandler<3>::distribute_dofs_on_cell (cell_iterator &cell,
                                         // distribute dofs if necessary:
                                         // check whether line dof is already
                                         // numbered (check only first dof)
-       if (quad->mg_dof_index(0) == DoFHandler<3>::invalid_dof_index)
+       if (quad->mg_dof_index(cell->level(), 0) == DoFHandler<3>::invalid_dof_index)
                                           // if not: distribute dofs
          for (unsigned int d=0; d<this->selected_fe->dofs_per_quad; ++d)
-           quad->set_mg_dof_index (d, next_free_dof++);            
+           quad->set_mg_dof_index (cell->level(), d, next_free_dof++);     
       };
 
 
                                   // dofs of cell
   if (this->selected_fe->dofs_per_hex > 0)
     for (unsigned int d=0; d<this->selected_fe->dofs_per_hex; ++d)
-      cell->set_mg_dof_index (d, next_free_dof++);
+      cell->set_mg_dof_index (cell->level(), d, next_free_dof++);
 
 
                                   // note that this cell has been processed
@@ -1975,7 +1975,7 @@ void MGDoFHandler<2>::renumber_dofs (const unsigned int  level,
        if (line->user_flag_set())
          {
            for (unsigned int d=0; d<this->selected_fe->dofs_per_line; ++d)
-             line->set_mg_dof_index (d, new_numbers[line->mg_dof_index(d)]);
+             line->set_mg_dof_index (level, d, new_numbers[line->mg_dof_index(level, d)]);
            line->clear_user_flag();
          }
                                       // finally, restore user flags
@@ -2042,7 +2042,7 @@ void MGDoFHandler<3>::renumber_dofs (const unsigned int  level,
        if (line->user_flag_set())
          {
            for (unsigned int d=0; d<this->selected_fe->dofs_per_line; ++d)
-             line->set_mg_dof_index (d, new_numbers[line->mg_dof_index(d)]);
+             line->set_mg_dof_index (level, d, new_numbers[line->mg_dof_index(level, d)]);
            line->clear_user_flag();
          }
                                       // finally, restore user flags
@@ -2074,7 +2074,7 @@ void MGDoFHandler<3>::renumber_dofs (const unsigned int  level,
        if (quad->user_flag_set())
          {
            for (unsigned int d=0; d<this->selected_fe->dofs_per_quad; ++d)
-             quad->set_mg_dof_index (d, new_numbers[quad->mg_dof_index(d)]);
+             quad->set_mg_dof_index (level, d, new_numbers[quad->mg_dof_index(level, d)]);
            quad->clear_user_flag();
          }
                                       // finally, restore user flags
index 124e7b3450a5436b199ede4cd52faefe7b6c752e..b10aca91d43a8fdfbdf7cfca9ae88aea44bd5af5 100644 (file)
@@ -1111,7 +1111,7 @@ MGTools::make_boundary_list(
                                              // boundary values of dofs on this
                                              // face
             face_dofs.resize (fe.dofs_per_face);
-            face->get_mg_dof_indices (face_dofs);
+            face->get_mg_dof_indices (level, face_dofs);
             if (fe_is_system)
               {
                                                  // enter those dofs
index f71d9e98750fbaf96982450c64558b2d723e91a8..e916067704c6bec6d3237b35974a0aad7ca4cdde 100644 (file)
@@ -91,7 +91,7 @@ MGSmootherContinuous<VECTOR>::MGSmootherContinuous (
               == level-1))
            {
                                               // get indices of this face
-             cell->face(face)->get_mg_dof_indices (dofs_on_face);
+             cell->face(face)->get_mg_dof_indices (level, dofs_on_face);
                                               // append them to the levelwise
                                               // list
              boundary_dofs.insert (boundary_dofs.end(),

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.