]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Have a DoFAccessor::get_fe() function, and use it instead of dof_handler->selected_fe...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 10 Jan 2003 16:37:09 +0000 (16:37 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 10 Jan 2003 16:37:09 +0000 (16:37 +0000)
git-svn-id: https://svn.dealii.org/trunk@6911 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_accessor.h
deal.II/deal.II/include/dofs/dof_accessor.templates.h
deal.II/deal.II/include/dofs/dof_handler.h
deal.II/deal.II/source/dofs/dof_accessor.cc
deal.II/doc/news/2002/c-3-4.html

index 15800badec733965f02e37ae8ed2f07baca31586..a41663021eba77ef97fcc0418ad03300f64308cc 100644 (file)
@@ -26,6 +26,7 @@ template <typename number> class SparseMatrix;
 template <typename number> class Vector;
 
 template <int dim> class DoFHandler;
+template <int dim> class FiniteElement;
 
 template <int celldim, int dim> class DoFObjectAccessor;
 template <int dim>              class DoFObjectAccessor<0, dim>;
@@ -100,6 +101,14 @@ class DoFAccessor
     const DoFHandler<dim> &
     get_dof_handler () const;
 
+                                     /**
+                                      * Return the finite element that
+                                      * is used on the cell pointed to
+                                      * by this iterator.
+                                      */
+    const FiniteElement<dim> &
+    get_fe () const;
+    
                                     /**
                                      * Copy operator.
                                      */
index 81e62c226c8d180cabb144032256130493412527..d4a7da8e9209223ecf87be962373ffe4ca6f8d29 100644 (file)
@@ -55,6 +55,7 @@ DoFAccessor<dim>::set_dof_handler (DoFHandler<dim> *dh)
 
 
 template <int dim>
+inline
 const DoFHandler<dim> &
 DoFAccessor<dim>::get_dof_handler () const
 {
@@ -63,6 +64,16 @@ DoFAccessor<dim>::get_dof_handler () const
 
 
 
+template <int dim>
+inline
+const FiniteElement<dim> &
+DoFAccessor<dim>::get_fe () const
+{
+  return *dof_handler->selected_fe;
+}
+
+
+
 template <int dim>
 inline
 DoFAccessor<dim> &
@@ -96,12 +107,12 @@ DoFObjectAccessor<1,dim>::dof_index (const unsigned int i) const
   Assert (this->dof_handler != 0, typename BaseClass::ExcInvalidObject());
                                   // make sure a FE has been selected
                                   // and enough room was reserved
-  Assert (this->dof_handler->selected_fe != 0, typename BaseClass::ExcInvalidObject());
-  Assert (i<this->dof_handler->selected_fe->dofs_per_line,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_line));
+  Assert (&this->get_fe() != 0, typename BaseClass::ExcInvalidObject());
+  Assert (i<this->get_fe().dofs_per_line,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_line));
 
   return this->dof_handler->levels[this->present_level]
-    ->line_dofs[this->present_index*this->dof_handler->selected_fe->dofs_per_line+i];
+    ->line_dofs[this->present_index*this->get_fe().dofs_per_line+i];
 }
 
 
@@ -123,13 +134,13 @@ DoFObjectAccessor<1,dim>::vertex_dof_index (const unsigned int vertex,
   typedef DoFAccessor<dim> BaseClass;
   
   Assert (this->dof_handler != 0, typename BaseClass::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0, typename BaseClass::ExcInvalidObject());
+  Assert (&this->get_fe() != 0, typename BaseClass::ExcInvalidObject());
   Assert (vertex<2, ExcIndexRange (i,0,2));
-  Assert (i<this->dof_handler->selected_fe->dofs_per_vertex,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_vertex));
+  Assert (i<this->get_fe().dofs_per_vertex,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_vertex));
 
   const unsigned int dof_number = (this->vertex_index(vertex) *
-                                  this->dof_handler->selected_fe->dofs_per_vertex +
+                                  this->get_fe().dofs_per_vertex +
                                   i);
   return this->dof_handler->vertex_dofs[dof_number];
 }
@@ -152,7 +163,7 @@ DoFObjectAccessor<1,dim>::get_dof_indices (std::vector<unsigned int> &dof_indice
   typedef DoFAccessor<dim> BaseClass;
   
   Assert (this->dof_handler != 0, typename BaseClass::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0, typename BaseClass::ExcInvalidObject());
+  Assert (&this->get_fe() != 0, typename BaseClass::ExcInvalidObject());
   Assert (dof_indices.size() == (2*this->dof_handler->get_fe().dofs_per_vertex +
                                 this->dof_handler->get_fe().dofs_per_line),
          typename BaseClass::ExcVectorDoesNotMatch());
@@ -219,13 +230,13 @@ unsigned int DoFObjectAccessor<2,dim>::dof_index (const unsigned int i) const
          typename DoFAccessor<dim>::ExcInvalidObject());
                                   // make sure a FE has been selected
                                   // and enough room was reserved
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (i<this->dof_handler->selected_fe->dofs_per_quad,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_quad));
+  Assert (i<this->get_fe().dofs_per_quad,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_quad));
 
   return this->dof_handler->levels[this->present_level]
-    ->quad_dofs[this->present_index*this->dof_handler->selected_fe->dofs_per_quad+i];
+    ->quad_dofs[this->present_index*this->get_fe().dofs_per_quad+i];
 }
 
 
@@ -237,14 +248,14 @@ DoFObjectAccessor<2,dim>::vertex_dof_index (const unsigned int vertex,
 {
   Assert (this->dof_handler != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
   Assert (vertex<4, ExcIndexRange (i,0,4));
-  Assert (i<this->dof_handler->selected_fe->dofs_per_vertex,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_vertex));
+  Assert (i<this->get_fe().dofs_per_vertex,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_vertex));
 
   const unsigned int dof_number = (this->vertex_index(vertex) *
-                                  this->dof_handler->selected_fe->dofs_per_vertex +
+                                  this->get_fe().dofs_per_vertex +
                                   i);
   return this->dof_handler->vertex_dofs[dof_number];
 }
@@ -263,7 +274,7 @@ DoFObjectAccessor<2,dim>::get_dof_indices (std::vector<unsigned int> &dof_indice
 {
   Assert (this->dof_handler != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
   Assert (dof_indices.size() == (4*this->dof_handler->get_fe().dofs_per_vertex +
                                 4*this->dof_handler->get_fe().dofs_per_line +
@@ -356,13 +367,13 @@ DoFObjectAccessor<3,dim>::dof_index (const unsigned int i) const
          typename DoFAccessor<dim>::ExcInvalidObject());
                                   // make sure a FE has been selected
                                   // and enough room was reserved
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (i<this->dof_handler->selected_fe->dofs_per_hex,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_hex));
+  Assert (i<this->get_fe().dofs_per_hex,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_hex));
 
   return this->dof_handler->levels[this->present_level]
-    ->hex_dofs[this->present_index*this->dof_handler->selected_fe->dofs_per_hex+i];
+    ->hex_dofs[this->present_index*this->get_fe().dofs_per_hex+i];
 }
 
 
@@ -374,14 +385,14 @@ DoFObjectAccessor<3,dim>::vertex_dof_index (const unsigned int vertex,
 {
   Assert (this->dof_handler != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
   Assert (vertex<8, ExcIndexRange (i,0,8));
-  Assert (i<this->dof_handler->selected_fe->dofs_per_vertex,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_vertex));
+  Assert (i<this->get_fe().dofs_per_vertex,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_vertex));
 
   const unsigned int dof_number = (this->vertex_index(vertex) *
-                                  this->dof_handler->selected_fe->dofs_per_vertex +
+                                  this->get_fe().dofs_per_vertex +
                                   i);
   return this->dof_handler->vertex_dofs[dof_number];
 }
@@ -394,7 +405,7 @@ DoFObjectAccessor<3,dim>::get_dof_indices (std::vector<unsigned int> &dof_indice
 {
   Assert (this->dof_handler != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
   Assert (dof_indices.size() == (8*this->dof_handler->get_fe().dofs_per_vertex +
                                 12*this->dof_handler->get_fe().dofs_per_line +
index ee5c0864a0e6ef6c3450650608445d791fdb4f68..1f1b4f6f33592317b5b18be48f6728a32c859c86 100644 (file)
@@ -23,6 +23,7 @@
 #include <map>
 #include <set>
 
+template <int dim> class DoFAccessor;
 template <int dim> class DoFCellAccessor;
 template <int dim> class DoFLevel;
 template <int celldim, int dim> class DoFObjectAccessor;
@@ -1143,6 +1144,11 @@ class DoFHandler  :  public Subscriptor,
                                      */
     std::vector<unsigned int>      vertex_dofs;
 
+                                    /**
+                                     * Make accessor objects friends.
+                                     */
+    template <int dim1> friend class DoFAccessor;
+
                                     /**
                                      * Make accessor objects friends.
                                      */
index 2c35389d5b02ac3833f3bcb92acbcbde7b383828..7470d91d56317ca66a552592941348b52e6efd95 100644 (file)
@@ -38,13 +38,13 @@ void DoFObjectAccessor<1, dim>::set_dof_index (const unsigned int i,
          typename DoFAccessor<dim>::ExcInvalidObject());
                                   // make sure a FE has been selected
                                   // and enough room was reserved
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (i<this->dof_handler->selected_fe->dofs_per_line,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_line));
+  Assert (i<this->get_fe().dofs_per_line,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_line));
 
   this->dof_handler->levels[this->present_level]
-    ->line_dofs[this->present_index*this->dof_handler->selected_fe->dofs_per_line+i] = index;
+    ->line_dofs[this->present_index*this->get_fe().dofs_per_line+i] = index;
 }
 
 
@@ -68,14 +68,14 @@ void DoFObjectAccessor<1, dim>::set_vertex_dof_index (const unsigned int vertex,
   
   Assert (this->dof_handler != 0,
          typename BaseClass::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename BaseClass::ExcInvalidObject());
   Assert (vertex<2, ExcIndexRange (i,0,2));
-  Assert (i<this->dof_handler->selected_fe->dofs_per_vertex,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_vertex));
+  Assert (i<this->get_fe().dofs_per_vertex,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_vertex));
 
   const unsigned int dof_number = (this->vertex_index(vertex) *
-                                  this->dof_handler->selected_fe->dofs_per_vertex +
+                                  this->get_fe().dofs_per_vertex +
                                   i);
   this->dof_handler->vertex_dofs[dof_number] = index;
 }
@@ -100,7 +100,7 @@ distribute_local_to_global (const Vector<double> &local_source,
   
   Assert (this->dof_handler != 0,
          typename BaseClass::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename BaseClass::ExcInvalidObject());
   Assert (local_source.size() == (2*this->dof_handler->get_fe().dofs_per_vertex +
                                  this->dof_handler->get_fe().dofs_per_line),
@@ -139,7 +139,7 @@ distribute_local_to_global (const FullMatrix<double> &local_source,
   
   Assert (this->dof_handler != 0,
          typename BaseClass::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename BaseClass::ExcInvalidObject());
   Assert (local_source.m() == (2*this->dof_handler->get_fe().dofs_per_vertex +
                               this->dof_handler->get_fe().dofs_per_line),
@@ -237,13 +237,13 @@ void DoFObjectAccessor<2, dim>::set_dof_index (const unsigned int i,
          typename DoFAccessor<dim>::ExcInvalidObject());
                                   // make sure a FE has been selected
                                   // and enough room was reserved
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (i<this->dof_handler->selected_fe->dofs_per_quad,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_quad));
+  Assert (i<this->get_fe().dofs_per_quad,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_quad));
 
   this->dof_handler->levels[this->present_level]
-    ->quad_dofs[this->present_index*this->dof_handler->selected_fe->dofs_per_quad+i] = index;
+    ->quad_dofs[this->present_index*this->get_fe().dofs_per_quad+i] = index;
 }
 
 
@@ -256,14 +256,14 @@ DoFObjectAccessor<2, dim>::set_vertex_dof_index (const unsigned int vertex,
 {
   Assert (this->dof_handler != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
   Assert (vertex<4, ExcIndexRange (i,0,4));
-  Assert (i<this->dof_handler->selected_fe->dofs_per_vertex,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_vertex));
+  Assert (i<this->get_fe().dofs_per_vertex,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_vertex));
 
   const unsigned int dof_number = (this->vertex_index(vertex) *
-                                  this->dof_handler->selected_fe->dofs_per_vertex +
+                                  this->get_fe().dofs_per_vertex +
                                   i);
   this->dof_handler->vertex_dofs[dof_number] = index;
 }
@@ -276,7 +276,7 @@ distribute_local_to_global (const Vector<double> &local_source,
                            Vector<double>       &global_destination) const {
   Assert (this->dof_handler != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
   Assert (local_source.size() == (4*this->dof_handler->get_fe().dofs_per_vertex +
                                  4*this->dof_handler->get_fe().dofs_per_line +
@@ -304,7 +304,7 @@ distribute_local_to_global (const FullMatrix<double> &local_source,
                            SparseMatrix<double>     &global_destination) const {
   Assert (this->dof_handler != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
   Assert (local_source.m() == (4*this->dof_handler->get_fe().dofs_per_vertex +
                               4*this->dof_handler->get_fe().dofs_per_line +
@@ -417,13 +417,13 @@ void DoFObjectAccessor<3, dim>::set_dof_index (const unsigned int i,
          typename DoFAccessor<dim>::ExcInvalidObject());
                                   // make sure a FE has been selected
                                   // and enough room was reserved
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (i<this->dof_handler->selected_fe->dofs_per_hex,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_hex));
+  Assert (i<this->get_fe().dofs_per_hex,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_hex));
 
   this->dof_handler->levels[this->present_level]
-    ->hex_dofs[this->present_index*this->dof_handler->selected_fe->dofs_per_hex+i] = index;
+    ->hex_dofs[this->present_index*this->get_fe().dofs_per_hex+i] = index;
 }
 
 
@@ -435,15 +435,15 @@ void DoFObjectAccessor<3, dim>::set_vertex_dof_index (const unsigned int vertex,
 {
   Assert (this->dof_handler != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
   Assert (vertex<8,
          ExcIndexRange (i,0,8));
-  Assert (i<this->dof_handler->selected_fe->dofs_per_vertex,
-         ExcIndexRange (i, 0, this->dof_handler->selected_fe->dofs_per_vertex));
+  Assert (i<this->get_fe().dofs_per_vertex,
+         ExcIndexRange (i, 0, this->get_fe().dofs_per_vertex));
 
   const unsigned int dof_number = (this->vertex_index(vertex) *
-                                  this->dof_handler->selected_fe->dofs_per_vertex +
+                                  this->get_fe().dofs_per_vertex +
                                   i);
   this->dof_handler->vertex_dofs[dof_number] = index;
 }
@@ -457,7 +457,7 @@ distribute_local_to_global (const Vector<double> &local_source,
 {
   Assert (this->dof_handler != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
   Assert (local_source.size() == (8*this->dof_handler->get_fe().dofs_per_vertex +
                                  12*this->dof_handler->get_fe().dofs_per_line +
@@ -487,7 +487,7 @@ distribute_local_to_global (const FullMatrix<double> &local_source,
 {
   Assert (this->dof_handler != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
-  Assert (this->dof_handler->selected_fe != 0,
+  Assert (&this->get_fe() != 0,
          typename DoFAccessor<dim>::ExcInvalidObject());
   Assert (local_source.m() == (8*this->dof_handler->get_fe().dofs_per_vertex +
                               12*this->dof_handler->get_fe().dofs_per_line +
index 759cdf78b3d4d0f2243092226b29d4f9f8e2a10c..6c3d7ef38702fde66f7aa3712db7fc1e91ee1ec9 100644 (file)
@@ -539,6 +539,17 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p> 
+       New: The DoF accessor classes now have a function <code
+       class="member">get_fe()</code> that returns a reference to the finite
+       element being used on this cell. The result is of course identical to
+       what a call to <code
+       class="member">iterator-&gt;get_dof_handler().get_fe()</code> would have
+       yielded.
+       <br>
+       (WB 2003/01/10)
+       </p>
+
   <li> <p> 
        New: Checked in new <code class="class">GridGenerator</code> 
        member function <code class="member">half_hyper_ball</code>,

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.