From: guido Date: Wed, 7 Jul 1999 18:36:48 +0000 (+0000) Subject: New class DoFTools X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=36801c2fcda68a08d99c4464cd09649dfee2a33b;p=dealii-svn.git New class DoFTools component_to_base now in FiniteElementBase git-svn-id: https://svn.dealii.org/trunk@1552 0785d39b-7218-0410-832d-ea1e28bc413d --- 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 index 0000000000..92fea62c7d --- /dev/null +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -0,0 +1,42 @@ +// $Id$ +// Copyright Guido Kanschat, 1999 + +#ifndef __deal_dof_tools_H +#define __deal_dof_tools_H + +#include + +/** + * 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 + static void extract_dofs(const DoFHandler& dof, + const bit_vector& select, + bit_vector& selected_dofs); +}; + + +#endif diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index bccd934721..54bb895b1b 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -274,6 +274,20 @@ class FiniteElementBase : public Subscriptor, * component. */ pair 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 > 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 component_to_base_table; }; @@ -1475,6 +1503,17 @@ FiniteElementBase::face_system_to_component_index (unsigned int index) cons return face_system_to_component_table[index]; } +template +inline +unsigned int +FiniteElementBase::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 diff --git a/deal.II/deal.II/include/fe/fe_system.h b/deal.II/deal.II/include/fe/fe_system.h index fff981229c..d447286a52 100644 --- a/deal.II/deal.II/include/fe/fe_system.h +++ b/deal.II/deal.II/include/fe/fe_system.h @@ -405,20 +405,6 @@ class FESystem : public FiniteElement */ 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 component_to_base_table; /** * Helper function used in the constructor: @@ -534,8 +520,7 @@ template template FESystem::FESystem (const FE &fe, const unsigned int n_elements) : FiniteElement (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::FESystem (const FE1 &fe1, const unsigned int n1, const FE2 &fe2, const unsigned int n2) : FiniteElement (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::FESystem (const FE1 &fe1, const unsigned int n1, FiniteElement (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 index 0000000000..c8e7dc8217 --- /dev/null +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -0,0 +1,41 @@ +// $Id$ +// Copyright Guido Kanschat, 1999 + +#include +#include +#include +#include +#include +#include + +template +void +DoFTools::extract_dofs(const DoFHandler& dof, + const bit_vector& local_select, + bit_vector& selected_dofs) +{ + const FiniteElement& 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 indices(fe.total_dofs); + + DoFHandler::active_cell_iterator c; + for (c = dof.begin_active() ; c != dof.end() ; ++ c) + { + c->get_dof_indices(indices); + for (unsigned int i=0;i comp + = fe.system_to_component_index(i); + selected_dofs[indices[i]] = local_select[comp.first]; + } + } +} + +template void DoFTools::extract_dofs(const DoFHandler& dof, + const bit_vector& local_select, + bit_vector& selected_dofs); +//template DoFTools; diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index 728b108923..48a6c38eeb 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -148,7 +148,8 @@ FiniteElementBase::FiniteElementBase (const FiniteElementData &fe_data system_to_component_table(total_dofs), face_system_to_component_table(dofs_per_face), component_to_system_table(n_components, vector(total_dofs)), - face_component_to_system_table(n_components, vector(dofs_per_face)) + face_component_to_system_table(n_components, vector(dofs_per_face)), + component_to_base_table(n_components) { for (unsigned int i=0; i::children_per_cell; ++i) {