]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Replace FECollection::get_fe() and n_finite_elements() by operator[] and size()
authorRalf Hartmann <Ralf.Hartmann@dlr.de>
Fri, 10 Feb 2006 12:07:23 +0000 (12:07 +0000)
committerRalf Hartmann <Ralf.Hartmann@dlr.de>
Fri, 10 Feb 2006 12:07:23 +0000 (12:07 +0000)
git-svn-id: https://svn.dealii.org/trunk@12291 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_accessor.templates.h
deal.II/deal.II/include/fe/fe_collection.h
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/dofs/hp_dof_handler.cc
deal.II/deal.II/source/fe/fe_values.cc

index e75d761c05774b7fdbb81233097e2183df01aaf5..762bf9663ed61bbacde502158972f447e6e57080 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -1203,21 +1203,21 @@ inline
 const FiniteElement<1> &
 DoFObjectAccessor<1,1,hp::DoFHandler>::get_fe () const
 {
-  return dof_handler->finite_elements->get_fe(active_fe_index ());
+  return (*dof_handler->finite_elements)[active_fe_index ()];
 }
 template <>
 inline
 const FiniteElement<2> &
 DoFObjectAccessor<1,2,hp::DoFHandler>::get_fe () const
 {
-  return dof_handler->finite_elements->get_fe(active_fe_index ());
+  return (*dof_handler->finite_elements)[active_fe_index ()];
 }
 template <>
 inline
 const FiniteElement<3> &
 DoFObjectAccessor<1,3,hp::DoFHandler>::get_fe () const
 {
-  return dof_handler->finite_elements->get_fe(active_fe_index ());
+  return (*dof_handler->finite_elements)[active_fe_index ()];
 }
 
 
@@ -1405,7 +1405,7 @@ inline
 const FiniteElement<2> &
 DoFObjectAccessor<2,2,hp::DoFHandler>::get_fe () const
 {
-  return dof_handler->finite_elements->get_fe(active_fe_index ());
+  return (*dof_handler->finite_elements)[active_fe_index ()];
 }
 
 
@@ -1415,7 +1415,7 @@ inline
 const FiniteElement<3> &
 DoFObjectAccessor<2,3,hp::DoFHandler>::get_fe () const
 {
-  return dof_handler->finite_elements->get_fe(active_fe_index ());
+  return (*dof_handler->finite_elements)[active_fe_index ()];
 }
 
 
@@ -1539,7 +1539,7 @@ inline
 const FiniteElement<3> &
 DoFObjectAccessor<3,3,hp::DoFHandler>::get_fe () const
 {
-  return dof_handler->finite_elements->get_fe(active_fe_index ());
+  return (*dof_handler->finite_elements)[active_fe_index ()];
 }
 
 
index 69d26c6ba5e99758518ec065d18dcb2e1a2318bf..0f2fe622aa39d7c5eb2ec880f5999b0ff8671722 100644 (file)
@@ -72,13 +72,13 @@ namespace hp
                                         * collection.
                                         */
       const FiniteElement<dim> &
-      get_fe (const unsigned int index) const;
+      operator[] (const unsigned int index) const;
 
                                        /**
                                         * Return the number of finite element
                                         * objects stored in this collection.
                                         */
-      unsigned int n_finite_elements () const;
+      unsigned int size () const;
 
                                        /**
                                         * Return the number of vector
@@ -157,7 +157,7 @@ namespace hp
   template <int dim>
   inline
   unsigned int
-  FECollection<dim>::n_finite_elements () const 
+  FECollection<dim>::size () const 
   {
     return finite_elements.size();
   }
@@ -176,7 +176,7 @@ namespace hp
   template <int dim>
   inline
   const FiniteElement<dim> &
-  FECollection<dim>::get_fe (const unsigned int index) const 
+  FECollection<dim>::operator[] (const unsigned int index) const 
   {
     Assert (index < finite_elements.size(),
             ExcIndexRange (index, 0, finite_elements.size()));
index f54525037b7e8afa46d99a7a60754a590d7f88c0..871417e4d9f3f498c241ff20f0c22769fec1c480 100644 (file)
@@ -3706,12 +3706,12 @@ DoFTools::convert_couplings_to_blocks (
   const Table<2, Coupling>& table,
   std::vector<Table<2,Coupling> >& tables_by_block)
 {
-  const hp::FECollection<dim>& coll = dof_handler.get_fe();
-  tables_by_block.resize(coll.n_finite_elements());
+  const hp::FECollection<dim>& fe_collection = dof_handler.get_fe();
+  tables_by_block.resize(fe_collection.size());
 
-  for (unsigned int f=0;f<coll.n_finite_elements();++f)
+  for (unsigned int f=0;f<fe_collection.size();++f)
     {
-      const FiniteElement<dim>& fe = coll.get_fe(f);
+      const FiniteElement<dim>& fe = fe_collection[f];
       
       const unsigned int nb = fe.n_blocks();
       tables_by_block[f].reinit(nb, nb);
index d4d0f23fe5f2675d43bfe90b5f20011fa6024870..52eacb8bf4ae963e3833f76851e03c62b61dad05 100644 (file)
@@ -1934,7 +1934,7 @@ namespace hp
                                      // their size
     clear_space ();
 
-    vertex_dofs.resize(tria->vertices.size()*finite_elements->get_fe(0).dofs_per_vertex,
+    vertex_dofs.resize(tria->vertices.size()*(*finite_elements)[0].dofs_per_vertex,
                        invalid_dof_index);
 
     for (unsigned int i=0; i<tria->n_levels(); ++i) 
@@ -1948,8 +1948,8 @@ namespace hp
         for (unsigned int j = 0; j < tria->levels[i]->lines.lines.size(); ++j)
           {
             levels.back()->dof_line_index_offset[j] = dofs_for_lines;
-            dofs_for_lines += finite_elements->
-                              get_fe (levels.back ()->active_fe_indices[j]).dofs_per_line;
+            dofs_for_lines += (*finite_elements)[
+             levels.back ()->active_fe_indices[j]].dofs_per_line;
           }
 
         levels.back()->line_dofs = std::vector<unsigned int>(dofs_for_lines,
@@ -1984,7 +1984,7 @@ namespace hp
                                      // their size
     clear_space ();
   
-    vertex_dofs.resize(tria->vertices.size()*finite_elements->get_fe(0).dofs_per_vertex,
+    vertex_dofs.resize(tria->vertices.size()*(*finite_elements)[0].dofs_per_vertex,
                        invalid_dof_index);
 
 //TODO[?] Does not work for continuous FEs. Problem is at faces, which might have a different
@@ -2065,18 +2065,18 @@ namespace hp
              active_fe2 = active_fe1;
 
             if (active_fe1 == active_fe2)
-             dofs_for_lines += finite_elements->get_fe (active_fe1).dofs_per_line + 2;
+             dofs_for_lines += (*finite_elements)[active_fe1].dofs_per_line + 2;
             else
-             dofs_for_lines += finite_elements->get_fe (active_fe1).dofs_per_line +
-                                finite_elements->get_fe (active_fe2).dofs_per_line + 4;
+             dofs_for_lines += (*finite_elements)[active_fe1].dofs_per_line +
+                                (*finite_elements)[active_fe2].dofs_per_line + 4;
           }
 
         unsigned int dofs_for_quads = 0;
         for (unsigned int j = 0; j < tria->levels[i]->quads.quads.size (); ++j)
           {
             levels.back()->dof_quad_index_offset[j] = dofs_for_quads;
-            dofs_for_quads += finite_elements->
-                              get_fe (levels.back ()->active_fe_indices[j]).dofs_per_quad;
+            dofs_for_quads += (*finite_elements)[
+             levels.back ()->active_fe_indices[j]].dofs_per_quad;
           }
 
         levels.back()->line_dofs = std::vector<unsigned int> (dofs_for_lines,
@@ -2109,8 +2109,8 @@ namespace hp
                 static_cast<unsigned int>(-1);
             else
               {
-                unsigned int start_next_fe = dofs_for_lines + 
-                                             finite_elements->get_fe (active_fe1).dofs_per_line + 2;
+                const unsigned int start_next_fe =
+                 dofs_for_lines + (*finite_elements)[active_fe1].dofs_per_line + 2;
 
                 levels.back()->line_dofs[dofs_for_lines + 1] = start_next_fe;
 
@@ -2148,7 +2148,7 @@ namespace hp
                                      // their size
     clear_space ();
   
-    vertex_dofs.resize(tria->vertices.size()*finite_elements->get_fe(0).dofs_per_vertex,
+    vertex_dofs.resize(tria->vertices.size()*(*finite_elements)[0].dofs_per_vertex,
                        invalid_dof_index);
 
 //TODO[?] Does not work for continuous FEs. Problem is at faces, which might have a different
@@ -2189,8 +2189,8 @@ namespace hp
         for (unsigned int j=0; j<tria->levels[i]->hexes.hexes.size(); ++j)
           {
             levels.back()->dof_hex_index_offset[j] = dofs_for_hexes;
-            dofs_for_hexes += finite_elements->
-                              get_fe (levels.back ()->active_fe_indices[j]).dofs_per_hex;
+            dofs_for_hexes += (*finite_elements)[
+             levels.back ()->active_fe_indices[j]].dofs_per_hex;
           }
 
         levels.back()->line_dofs = std::vector<unsigned int> (dofs_for_lines,
index 476ca7b90d3c539bf6738ecaab3833083b842f7a..3f1732986ed8195443bb41cd951a056d3680f717 100644 (file)
@@ -1424,7 +1424,8 @@ void FEFaceValues<dim>::reinit (const typename hp::DoFHandler<dim>::cell_iterato
                                   // used by the DoFHandler used by
                                   // this cell, are the same
   Assert (static_cast<const FiniteElementData<dim>&>(*this->fe) ==
-         static_cast<const FiniteElementData<dim>&>(cell->get_dof_handler().get_fe().get_fe (cell->active_fe_index ())),
+         static_cast<const FiniteElementData<dim>&>(
+           cell->get_dof_handler().get_fe()[cell->active_fe_index ()]),
          typename FEValuesBase<dim>::ExcFEDontMatch());
 
   Assert (face_no < GeometryInfo<dim>::faces_per_cell,
@@ -1639,7 +1640,8 @@ void FESubfaceValues<dim>::reinit (const typename hp::DoFHandler<dim>::cell_iter
                                   // used by the hp::DoFHandler used by
                                   // this cell, are the same
   Assert (static_cast<const FiniteElementData<dim>&>(*this->fe) ==
-         static_cast<const FiniteElementData<dim>&>(cell->get_dof_handler().get_fe().get_fe (cell->active_fe_index ())),
+         static_cast<const FiniteElementData<dim>&>(
+           cell->get_dof_handler().get_fe()[cell->active_fe_index ()]),
          typename FEValuesBase<dim>::ExcFEDontMatch());
   Assert (face_no < GeometryInfo<dim>::faces_per_cell,
          ExcIndexRange (face_no, 0, GeometryInfo<dim>::faces_per_cell));

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.