]> https://gitweb.dealii.org/ - dealii.git/commitdiff
New hp::DoFHandler::set/get_active_fe_indices functions.
authorRalf Hartmann <Ralf.Hartmann@dlr.de>
Fri, 12 Jun 2009 06:07:59 +0000 (06:07 +0000)
committerRalf Hartmann <Ralf.Hartmann@dlr.de>
Fri, 12 Jun 2009 06:07:59 +0000 (06:07 +0000)
git-svn-id: https://svn.dealii.org/trunk@18929 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/hp/dof_handler.h
deal.II/deal.II/source/hp/dof_handler.cc
deal.II/doc/news/changes.h

index 7636d273797e5719b84e953ead6e5bb3d6b1ce4d..4a34f2a4b380d96e2fb62d7635decdf67c5c72d6 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2005, 2006, 2007, 2008 by the deal.II authors
+//    Copyright (C) 2005, 2006, 2007, 2008, 2009 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -172,8 +172,8 @@ namespace hp
                                         * Go through the triangulation and
                                         * distribute the degrees of freedoms
                                         * needed for the given finite element
-                                        * according to the given distribution
-                                        * method.
+                                        * according to the current distribution
+                                        * of active fe indices.
                                         *
                                         * A pointer of the transferred
                                         * finite element is
@@ -189,6 +189,23 @@ namespace hp
                                         */
       virtual void distribute_dofs (const hp::FECollection<dim,spacedim> &fe);
 
+                                       /**
+                                        * Go through the triangulation and set
+                                        * the active FE indices of all active
+                                        * cells to the values given in @p
+                                        * active_fe_indices.
+                                        */
+      void set_active_fe_indices (const std::vector<unsigned int>& active_fe_indices);
+
+                                       /**
+                                        * Go through the triangulation and
+                                        * store the active FE indices of all
+                                        * active cells to the vector @p
+                                        * active_fe_indices. This vector is
+                                        * resized, if necessary.
+                                        */
+      void get_active_fe_indices (std::vector<unsigned int>& active_fe_indices) const;
+
                                        /**
                                         * Clear all data of this object and
                                         * especially delete the lock this object
@@ -1084,7 +1101,7 @@ namespace hp
                                         * Store the number of dofs
                                         * created last time.
                                         */
-      unsigned int              used_dofs;
+      unsigned int used_dofs;
 
                                        /**
                                         * Array to store the indices
index 3af05721cb34bcb42d1013515838226a4e23c997..5d312bcfc7cf48ce069c4efb7e57c371a2cf154f 100644 (file)
@@ -3397,6 +3397,44 @@ namespace hp
        }
   }
   
+
+  
+  template <int dim, int spacedim>
+  void DoFHandler<dim,spacedim>::set_active_fe_indices (const std::vector<unsigned int>& active_fe_indices)
+  {
+    Assert(active_fe_indices.size()==get_tria().n_active_cells(),
+          ExcDimensionMismatch(active_fe_indices.size(), get_tria().n_active_cells()));
+
+    create_active_fe_table ();
+                                    // we could set the values directly, since
+                                    // they are stored as protected data of
+                                    // this object, but for simplicity we use
+                                    // the cell-wise access. this way we also
+                                    // have to pass some debug-mode tests which
+                                    // we would have to duplicate ourselves
+                                    // otherwise
+    active_cell_iterator cell=begin_active(),
+                        endc=end();
+    for (unsigned int i=0; cell!=endc; ++cell, ++i)
+      cell->set_active_fe_index(active_fe_indices[i]);
+  }
+
+  
+
+  template <int dim, int spacedim>
+  void DoFHandler<dim,spacedim>::get_active_fe_indices (std::vector<unsigned int>& active_fe_indices) const
+  {
+    active_fe_indices.resize(get_tria().n_active_cells());
+
+                                    // we could try to extract the values directly, since
+                                    // they are stored as protected data of
+                                    // this object, but for simplicity we use
+                                    // the cell-wise access.
+    active_cell_iterator cell=begin_active(),
+                        endc=end();
+    for (unsigned int i=0; cell!=endc; ++cell, ++i)
+      active_fe_indices[i]=cell->active_fe_index();
+  }
   
 
   
index 9c867c80939fe65c07b7c409deaa9e4e470500b0..e9ad4a1c0657d0ec453d4798e1b5c46db569b922 100644 (file)
@@ -135,6 +135,20 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
+  <li>
+  <p>
+  New: The new hp::DoFHandler::set_active_fe_indices function allows
+  to distribute all active FE indices at once based on a given
+  vector. This might be useful if this information is stored
+  somewhere and has to be reconstructed or else if two DoFHandler
+  objects with the same FE index distribution should be created.
+  There is now also a corresponding
+  hp::DoFHandler::get_active_fe_indices function.  
+  <br>
+  (Tobias Leicht, RH 2009/06/12)
+  </p>
+  </li>
+
   <li>
   <p>
   Fix: The projection of quadrature points to subfaces in

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.