]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New class DoFTools
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 7 Jul 1999 18:36:48 +0000 (18:36 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 7 Jul 1999 18:36:48 +0000 (18:36 +0000)
component_to_base now in FiniteElementBase

git-svn-id: https://svn.dealii.org/trunk@1552 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_tools.h [new file with mode: 0644]
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/source/dofs/dof_tools.cc [new file with mode: 0644]
deal.II/deal.II/source/fe/fe.cc

diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h
new file mode 100644 (file)
index 0000000..92fea62
--- /dev/null
@@ -0,0 +1,42 @@
+// $Id$
+// Copyright Guido Kanschat, 1999
+
+#ifndef __deal_dof_tools_H
+#define __deal_dof_tools_H
+
+#include <bvector.h>
+
+/**
+ * Operations on DoF-numbers.
+ * This is a collection of functions manipulating the numbers of
+ * degrees of freedom. The documentation of the member functions will
+ * provide more information.
+ *
+ * All member functions are static, so there is no need to create an
+ * object of class #DoFTools#.
+ * @author Guido Kanschat, 1999
+ */
+class DoFTools
+{
+  public:
+                                    /**
+                                     * Extract DoFs of components.
+                                     * The bit vector #select#
+                                     * defines, which components of an
+                                     * #FESystem# are to be extracted
+                                     * from the DoFHandler #dof#. The
+                                     * numbers of these dofs are then
+                                     * entered consecutively into
+                                     * #selected_dofs#.
+                                     *
+                                     * A prior ordering of dof-values
+                                     * is not destroyed by this process.
+                                     */
+    template<int dim>
+    static void extract_dofs(const DoFHandler<dim>& dof,
+                            const bit_vector& select,
+                            bit_vector& selected_dofs);
+};
+
+
+#endif
index bccd93472161e88b5e27e5e96b92ae54fa73a81d..54bb895b1be965395f81f6f0a4ec967397d669f3 100644 (file)
@@ -274,6 +274,20 @@ class FiniteElementBase : public Subscriptor,
                                      * component.
                                      */
     pair<unsigned int,unsigned int> face_system_to_component_index (unsigned int index) const; 
+                                    /**
+                                     * The base element establishing a
+                                     * component.
+                                     *
+                                     * This table converts a
+                                     * component number to the
+                                     * #base_element# number. While
+                                     * component information contains
+                                     * multiplicity of base elements,
+                                     * the result allows access to
+                                     * shape functions of the base
+                                     * element.
+                                     */
+    unsigned int component_to_base(unsigned int index) const;
     
 
                                     /**
@@ -470,6 +484,20 @@ class FiniteElementBase : public Subscriptor,
                                      * Map between component and linear dofs on a face.
                                      */
     vector< vector<unsigned int> > face_component_to_system_table;
+                                    /**
+                                     * The base element establishing a
+                                     * component.
+                                     *
+                                     * This table converts a
+                                     * component number to the
+                                     * #base_element# number. While
+                                     * component information contains
+                                     * multiplicity of base elements,
+                                     * the result allows access to
+                                     * shape functions of the base
+                                     * element.
+                                     */
+    vector<unsigned int> component_to_base_table;
 };
 
 
@@ -1475,6 +1503,17 @@ FiniteElementBase<dim>::face_system_to_component_index (unsigned int index) cons
   return face_system_to_component_table[index];
 }
 
+template <int dim>  
+inline
+unsigned int
+FiniteElementBase<dim>::component_to_base (unsigned int index) const
+{
+  if (n_components == 1)
+    return 0;
+  Assert(index < component_to_base_table.size(), ExcInvalidIndex(index));
+  return component_to_base_table[index];
+}
+
 /*----------------------------   fe.h     ---------------------------*/
 /* end of #ifndef __fe_H */
 #endif
index fff981229cc6f2b39fa647adb95ac479e4eb7de8..d447286a52d75b0cf5a92cca550d38e1dd279229 100644 (file)
@@ -405,20 +405,6 @@ class FESystem : public FiniteElement<dim>
                                      */
     vector< ElementPair > base_elements;
 
-                                    /**
-                                     * The base element establishing a
-                                     * component.
-                                     *
-                                     * This table converts a
-                                     * component number to the
-                                     * #base_element# number. While
-                                     * component information contains
-                                     * multiplicity of base elements,
-                                     * the result allows access to
-                                     * shape functions of the base
-                                     * element.
-                                     */
-    vector<unsigned int> component_to_base_table;
     
                                     /**
                                      * Helper function used in the constructor:
@@ -534,8 +520,7 @@ template <int dim>
 template <class FE>
 FESystem<dim>::FESystem (const FE &fe, const unsigned int n_elements) :
                FiniteElement<dim> (multiply_dof_numbers(fe, n_elements)),
-               base_elements(1),
-               component_to_base_table(n_components)
+               base_elements(1)
 {
   base_elements[0] = ElementPair(new FE, n_elements);
   base_elements[0].first -> subscribe ();
@@ -550,8 +535,7 @@ FESystem<dim>::FESystem (const FE1 &fe1, const unsigned int n1,
                         const FE2 &fe2, const unsigned int n2)
                :
                FiniteElement<dim> (multiply_dof_numbers(fe1, n1, fe2, n2)),
-               base_elements(2),
-               component_to_base_table(n_components)
+               base_elements(2)
 {
   Assert(fe1.n_transform_functions == fe2.n_transform_functions,
         ExcElementTransformNotEqual());
@@ -574,8 +558,7 @@ FESystem<dim>::FESystem (const FE1 &fe1, const unsigned int n1,
                FiniteElement<dim> (multiply_dof_numbers(fe1, n1,
                                                         fe2, n2,
                                                         fe3, n3)),
-               base_elements(3),
-               component_to_base_table(n_components)
+               base_elements(3)
 {
   Assert(fe1.n_transform_functions == fe2.n_transform_functions,
         ExcElementTransformNotEqual());
diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc
new file mode 100644 (file)
index 0000000..c8e7dc8
--- /dev/null
@@ -0,0 +1,41 @@
+// $Id$
+// Copyright Guido Kanschat, 1999
+
+#include <grid/dof.h>
+#include <grid/dof_accessor.h>
+#include <grid/tria_iterator.h>
+#include <fe/fe.h>
+#include <fe/fe_system.h>
+#include <basic/dof_tools.h>
+
+template<int dim>
+void
+DoFTools::extract_dofs(const DoFHandler<dim>& dof,
+                      const bit_vector& local_select,
+                      bit_vector& selected_dofs)
+{
+  const FiniteElement<dim>& fe = dof.get_fe();
+  Assert(local_select.size() == fe.n_components,
+        ExcDimensionMismatch(local_select.size(), fe.n_components));
+  Assert(selected_dofs.size() == dof.n_dofs(),
+        ExcDimensionMismatch(selected_dofs.size(), dof.n_dofs()));
+
+  vector<int> indices(fe.total_dofs);
+  
+  DoFHandler<dim>::active_cell_iterator c;
+  for (c = dof.begin_active() ; c != dof.end() ; ++ c)
+    {
+      c->get_dof_indices(indices);
+      for (unsigned int i=0;i<fe.total_dofs;++i)
+       {
+         pair<unsigned int, unsigned int> comp
+           = fe.system_to_component_index(i);
+         selected_dofs[indices[i]] = local_select[comp.first];
+       }
+    }
+}
+
+template void DoFTools::extract_dofs(const DoFHandler<deal_II_dimension>& dof,
+                                    const bit_vector& local_select,
+                                    bit_vector& selected_dofs);
+//template DoFTools<deal_II_dimension>;
index 728b1089235b47a6c85db0f0a4fed36b6f90d3cd..48a6c38eebcad2ea77d3480cddfd74d875de67d4 100644 (file)
@@ -148,7 +148,8 @@ FiniteElementBase<dim>::FiniteElementBase (const FiniteElementData<dim> &fe_data
                system_to_component_table(total_dofs),
                face_system_to_component_table(dofs_per_face),
                component_to_system_table(n_components, vector<unsigned>(total_dofs)),
-               face_component_to_system_table(n_components, vector<unsigned>(dofs_per_face))
+               face_component_to_system_table(n_components, vector<unsigned>(dofs_per_face)),
+               component_to_base_table(n_components)
 {
   for (unsigned int i=0; i<GeometryInfo<dim>::children_per_cell; ++i) 
     {

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.