]> https://gitweb.dealii.org/ - dealii.git/commitdiff
new class BlockInfo
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Sun, 13 Sep 2009 17:37:30 +0000 (17:37 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Sun, 13 Sep 2009 17:37:30 +0000 (17:37 +0000)
git-svn-id: https://svn.dealii.org/trunk@19449 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/block_info.h [new file with mode: 0644]
deal.II/deal.II/source/dofs/block_info.cc [new file with mode: 0644]
deal.II/doc/news/changes.h

diff --git a/deal.II/deal.II/include/dofs/block_info.h b/deal.II/deal.II/include/dofs/block_info.h
new file mode 100644 (file)
index 0000000..93e1a3c
--- /dev/null
@@ -0,0 +1,189 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2009 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+#ifndef __deal2__block_info_h
+#define __deal2__block_info_h
+
+#include <base/subscriptor.h>
+#include <lac/block_indices.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+// Forward declarations
+
+template <int dim, int spacedim> class DoFHandler;
+template <int dim, int spacedim> class MGDoFHandler;
+namespace hp
+{
+  template <int dim, int spacedim> class DoFHandler;
+}
+
+
+/**
+ * @brief A small class collecting the different BlockIndices involved in
+ * global, multilevel and local computations.
+ *
+ * In addition to grouping all the block indices into a single
+ * instance, this class also provides functions to fill them
+ * automatically.
+ *
+ * @todo Extend the functions local() and renumber() to the concept to
+ * hpDoFHandler.
+ *
+ * @author Guido Kanschat, 2009
+ */
+class BlockInfo : public Subscriptor
+{
+  public:
+                                    /**
+                                     * @brief Fill the object with values
+                                     * describing block structure
+                                     * of the DoFHandler.
+                                     *
+                                     * This function will also clear
+                                     * the local() indices.
+                                     */
+    template <int dim, int spacedim>
+    void initialize(const DoFHandler<dim, spacedim>&);
+    
+                                    /**
+                                     * @brief Fill the object with values
+                                     * describing level block
+                                     * structure of the MGDoFHandler.
+                                     *
+                                     * This function will also clear
+                                     * the local() indices.
+                                     */
+    template <int dim, int spacedim>
+    void initialize(const MGDoFHandler<dim, spacedim>&);
+    
+                                    /**
+                                     * @brief Initialize block structure
+                                     * on cells and compute
+                                     * renumbering between cell
+                                     * dofs and block cell dofs.
+                                     */
+    template <int dim, int spacedim>
+    void initialize_local(const DoFHandler<dim, spacedim>&);
+    
+                                    /**
+                                     * Access the BlockIndices
+                                     * structure of the global
+                                     * system.
+                                     */
+    const BlockIndices& global() const;
+    
+                                    /**
+                                     * Access BlockIndices for the
+                                     * local system on a cell.
+                                     */
+    const BlockIndices& local() const;
+    
+                                    /**
+                                     * Access the BlockIndices
+                                     * structure of a level in the
+                                     * multilevel hierarchy.
+                                     */
+    const BlockIndices& level(unsigned int level) const;
+
+                                    /**
+                                     * Return the index after local
+                                     * renumbering.
+                                     *
+                                     * The input of this function is
+                                     * an index between zero and the
+                                     * number of dofs per cell,
+                                     * numbered in local block
+                                     * ordering, that is first all
+                                     * indices of the first system
+                                     * block, then all of the second
+                                     * block and so forth. The
+                                     * function then outpus the index
+                                     * in the standard local
+                                     * numbering of DoFAccessor.
+                                     */
+    unsigned int renumber (unsigned int i) const;
+  private:
+                                    /**
+                                     * @brief The block structure
+                                     * of the global system.
+                                     */
+    BlockIndices bi_global;
+                                    /**
+                                     * @brief The multilevel block structure.
+                                     */
+    std::vector<BlockIndices> levels;
+    
+                                    /**
+                                     * @brief The block structure
+                                     * of the cell systems.
+                                     */
+    BlockIndices bi_local;
+    
+                                    /**
+                                     * The base element associated
+                                     * with each block.
+                                     */
+    std::vector<unsigned int> base_element;
+    
+                                    /**
+                                     * A vector containing the
+                                     * renumbering from the
+                                     * standard order of degrees of
+                                     * freedom on a cell to a
+                                     * component wise
+                                     * ordering. Filled by
+                                     * initialize().
+                                     */
+    std::vector<unsigned int> local_renumbering;
+};
+
+
+
+//----------------------------------------------------------------------//
+
+inline
+const BlockIndices&
+BlockInfo::global() const
+{
+  return bi_global;
+}
+
+
+inline
+const BlockIndices&
+BlockInfo::local() const
+{
+  return bi_local;
+}
+
+
+inline
+const BlockIndices&
+BlockInfo::level(unsigned int l) const
+{
+  AssertIndexRange(l, levels.size());
+  return levels[l];
+}
+
+
+inline
+unsigned int BlockInfo::renumber(unsigned int i) const
+{
+  AssertIndexRange(i, local_renumbering.size());
+  return local_renumbering[i];
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
diff --git a/deal.II/deal.II/source/dofs/block_info.cc b/deal.II/deal.II/source/dofs/block_info.cc
new file mode 100644 (file)
index 0000000..f6b27d9
--- /dev/null
@@ -0,0 +1,77 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 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
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+
+
+#include <dofs/block_info.h>
+#include <dofs/dof_handler.h>
+#include <dofs/dof_tools.h>
+#include <multigrid/mg_dof_handler.h>
+#include <multigrid/mg_tools.h>
+#include <fe/fe.h>
+#include <fe/fe_tools.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+
+template <int dim, int spacedim>
+void
+BlockInfo::initialize(const DoFHandler<dim, spacedim>& dof)
+{
+  const FiniteElement<dim, spacedim>& fe = dof.get_fe();
+  std::vector<unsigned int> sizes(fe.n_blocks());
+  DoFTools::count_dofs_per_block(dof, sizes);
+  bi_global.reinit(sizes);
+}
+
+
+template <int dim, int spacedim>
+void
+BlockInfo::initialize_local(const DoFHandler<dim, spacedim>& dof)
+{
+  const FiniteElement<dim, spacedim>& fe = dof.get_fe();
+  std::vector<unsigned int> sizes(fe.n_blocks());
+
+  base_element.resize(fe.n_blocks());
+  
+  for (unsigned int i=0;i<base_element.size();++i)
+    base_element[i] = fe.block_to_base_index(i).first;
+  
+  local_renumbering.resize(fe.n_dofs_per_cell());
+  FETools::compute_block_renumbering(fe,
+                                    local_renumbering,
+                                    sizes, false);
+  bi_local.reinit(sizes);
+}
+
+
+template <int dim, int spacedim>
+void
+BlockInfo::initialize(const MGDoFHandler<dim, spacedim>& dof)
+{
+  initialize((DoFHandler<dim, spacedim>&) dof);
+  
+  std::vector<std::vector<unsigned int> > sizes;
+  MGTools::count_dofs_per_block(dof, sizes);
+  
+  levels.resize(sizes.size());
+  for (unsigned int i=0;i<sizes.size();++i)
+    levels[i].reinit(sizes[i]);
+}
+
+
+template void BlockInfo::initialize(const DoFHandler<deal_II_dimension,deal_II_dimension>&);
+template void BlockInfo::initialize(const MGDoFHandler<deal_II_dimension,deal_II_dimension>&);
+template void BlockInfo::initialize_local(const DoFHandler<deal_II_dimension,deal_II_dimension>&);
+
+
+DEAL_II_NAMESPACE_CLOSE
index 280beefc6843bee5a0e6bdce005c18c1643feab2..cdf57d8249f323c571f27bd4c1fd81a4efb19008 100644 (file)
@@ -420,6 +420,19 @@ inconvenience this causes.
 
 <ol>
 
+  <li><p>New: The class BlockInfo stores and computes all BlockIndices objects related to
+  DoFHandler and MGDoFHandler based on FESystem.
+  <br>
+  (GK 2009/09/13)
+  </p>
+  </li>  
+
+  <li><p>Improved: MGDoFHandler now also has a typedef for MGDoFHandler::Container
+  <br>
+  (GK 2009/09/13)
+  </p>
+  </li>  
+
   <li><p> FETools::compute_block_renumbering() can nor return block sizes instead of
   start indices.
   <br>

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.