From: guido Date: Fri, 12 Feb 1999 21:33:03 +0000 (+0000) Subject: Approaching component tables X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ef28a2f3e1fe125c8fc3d2eb831800310180890f;p=dealii-svn.git Approaching component tables git-svn-id: https://svn.dealii.org/trunk@795 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index e4ae48e198..0579249e64 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -67,7 +67,7 @@ struct FiniteElementData<1> { /** * Number of components. */ - const unsigned int n_components; + const unsigned int n_components; /** * Default constructor. Constructs an element @@ -226,7 +226,7 @@ struct FiniteElementData<2> { * @author Wolfgang Bangerth, 1998 */ template -struct FiniteElementBase : +class FiniteElementBase : public Subscriptor, public FiniteElementData @@ -293,6 +293,10 @@ struct FiniteElementBase : /** * Compute component and index from system index. + * + * Return value contains first + * component and second index in + * component. */ pair system_to_component_index (unsigned index) const; @@ -471,6 +475,11 @@ struct FiniteElementBase : * space dimension. */ dFMatrix interface_constraints; + + /** + * Map between linear dofs and component dofs. + */ + vector< pair > system_to_component_table; }; @@ -1164,8 +1173,22 @@ class FiniteElement : public FiniteElementBase { DeclException0 (ExcJacobiDeterminantHasWrongSign); }; +template +inline unsigned +FiniteElementBase::component_to_system_index (unsigned component, + unsigned component_index) const +{ + Assert(false, ExcNotImplemented()); + return component_index; +} - +template +inline pair +FiniteElementBase::system_to_component_index (unsigned index) const +{ + Assert(index < system_to_component_table.size(), ExcInvalidIndex(index)); + return system_to_component_table[index]; +} /*---------------------------- fe.h ---------------------------*/ /* end of #ifndef __fe_H */ diff --git a/deal.II/deal.II/include/fe/fe_system.h b/deal.II/deal.II/include/fe/fe_system.h index 01651fdfb4..df5fa8cbce 100644 --- a/deal.II/deal.II/include/fe/fe_system.h +++ b/deal.II/deal.II/include/fe/fe_system.h @@ -1,12 +1,13 @@ /*---------------------------- fe_lib.system.h ---------------------------*/ /* $Id$ */ -#ifndef __fe_system_H -#define __fe_system_H +#ifndef __deal_fe_system_H +#define __deal_fe_system_H /*---------------------------- fe_lib.system.h ---------------------------*/ #include - +#include +#include /** * This class provides an interface to group several equal elements together @@ -67,11 +68,10 @@ template class FESystem : public FiniteElement { - /** - * Copy constructor prohibited. - */ - FESystem(const FESystem&); - + /** + * Copy constructor prohibited. + */ + FESystem(const FESystem&); public: @@ -303,55 +303,63 @@ class FESystem : public FiniteElement const unsigned int face_no, const unsigned int subface_no, const vector > &unit_points, - vector > &normal_vectors) const; + vector > &normal_vectors) const; - /** - * Number of subelements of this object. + /** + *Number of different composing + * elements of this object. + * * Since these objects can have - * subobjects themselves, this may be - * smaller than the total number of finite - * elements composed into this structure. - * This is definitely not what you'd - * usally intend, so don't do it! + * multiplicity and subobjects + * themselves, this may be + * smaller than the total number + * of finite elements composed + * into this structure. */ - const unsigned int n_sub_elements; - - /** - * Access to the single valued element. - * - * If you assemble your system - * matrix, you usually will not - * want to have an FEValues object - * with a lot of equal entries. Ok, - * so initialize your FEValues with - * the #base_element# yuo get by - * this function. - * - */ - const FiniteElement& get_base_element() const; - - /** - * Calculate the actual position. - * - * For a given #component# - * (e.g. u,v,w in the example - * above) of the #base# function of - * the #base_element#, return the - * actual index in the local - * degrees of freedom vector of the - * system. - * - */ - unsigned index(unsigned component, unsigned base) const; - + unsigned n_component_elements() const; + + /** + * How often is a composing element used. + * + */ + unsigned element_multiplicity(unsigned index) const; + + /** + * Access to a composing element. + * + * If you assemble your system + * matrix, you usually will not + * want to have an FEValues object + * with a lot of equal entries. Ok, + * so initialize your FEValues with + * the #base_element# you get by + * this function. In a mixed + * discretization, you can choose + * the different base element types + * by index. + * + */ + const FiniteElement& base_element(unsigned index) const; + private: /** - * Pointer to an object of the underlying - * finite element class. This object is - * created by the constructor. + * Pairs of multiplicity and element type. + */ + typedef pair *, unsigned > ElementPair; + + /** + * Pointer to underlying finite + * element classes. + * + * This object contains a pointer + * to each contributing element + * of a mixed discretization and + * its multiplicity. It is + * created by the constructor and + * constant afterwards. */ - const FiniteElement *const base_element; + vector< ElementPair > base_elements; /** @@ -365,12 +373,14 @@ class FESystem : public FiniteElement * transformation from unit to real * cell. */ - static FiniteElementData multiply_dof_numbers (const FiniteElementData &fe_data, - const unsigned int N); + static FiniteElementData + multiply_dof_numbers (const FiniteElementData &fe_data, + const unsigned int N); /** * This function is simply singled out of - * the constructor; it sets up the + * the constructor. It sets up the + * index table for the system as well as * #restriction# and #prolongation# * matrices. Since the operation of this * function can be done without explicit @@ -380,7 +390,7 @@ class FESystem : public FiniteElement * the general template definition in * the #.h# file. */ - void initialize_matrices (); + void initialize(); }; @@ -389,35 +399,45 @@ class FESystem : public FiniteElement /* ------------------------- template functions ------------------------- */ -template -inline const FiniteElement& -FESystem::get_base_element() const +template +inline unsigned +FESystem::n_component_elements() const { - return *base_element; + return base_elements.size(); } -template + +template inline unsigned -FESystem::index(unsigned component, unsigned base) const +FESystem::element_multiplicity(unsigned index) const { - return n_sub_elements * base + component; + return base_elements[index].second; } + template -template +inline const FiniteElement& +FESystem::base_element(unsigned index) const +{ + return *base_elements[index].first; +} + + +template +template FESystem::FESystem (const FE &fe, const unsigned int n_elements) : FiniteElement (multiply_dof_numbers(fe, n_elements)), - n_sub_elements (n_elements), - base_element (new FE()) + base_elements(1) { - base_element->subscribe (); - initialize_matrices (); + base_elements[0] = ElementPair(&fe, n_elements); + base_elements[0].first -> subscribe (); + initialize (); }; -/*---------------------------- fe_lib.system.h ---------------------------*/ +/*---------------------------- fe_lib.system.h ---------------------------*/ /* end of #ifndef __fe_system_H */ #endif -/*---------------------------- fe_lib.system.h ---------------------------*/ +/*---------------------------- fe_lib.system.h ---------------------------*/ diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index 2f2597f18b..c0c7901767 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -127,7 +127,6 @@ bool FiniteElementData<2>::operator== (const FiniteElementData<2> &f) const { /*------------------------------- FiniteElementBase ----------------------*/ - #if deal_II_dimension == 1 template <> @@ -142,6 +141,9 @@ FiniteElementBase<1>::FiniteElementBase (const FiniteElementData<1> &fe_data) : }; interface_constraints.reinit (1,1); interface_constraints(0,0)=1.; + + for (unsigned int j=0 ; j(0,j); }; #endif @@ -161,6 +163,8 @@ FiniteElementBase<2>::FiniteElementBase (const FiniteElementData<2> &fe_data) : }; interface_constraints.reinit (dofs_per_vertex+2*dofs_per_line, 2*dofs_per_vertex+dofs_per_line); + for (unsigned int j=0 ; j(0,j); }; #endif @@ -216,28 +220,6 @@ template FiniteElement::FiniteElement (const FiniteElementData &fe_data) : FiniteElementBase (fe_data) {}; - - -template -unsigned -FiniteElementBase::component_to_system_index (unsigned component, - unsigned component_index) const -{ - Assert(false, ExcInvalidIndex(component)); - return component_index; -} - - -template -pair -FiniteElementBase::system_to_component_index (unsigned index) const -{ - Assert(false, ExcInvalidIndex(index)); - return pair(0,0); -} - - - #if deal_II_dimension == 1 // declare this function to be explicitely specialized before first use diff --git a/deal.II/deal.II/source/fe/fe_system.cc b/deal.II/deal.II/source/fe/fe_system.cc index 1a43388b42..db56eae77a 100644 --- a/deal.II/deal.II/source/fe/fe_system.cc +++ b/deal.II/deal.II/source/fe/fe_system.cc @@ -11,15 +11,56 @@ template FESystem::~FESystem () { - base_element->unsubscribe (); - delete base_element; + for (unsigned i=0;i unsubscribe (); + delete base_elements[i].first; + } }; template -void FESystem::initialize_matrices () +void FESystem::initialize () { + // Initialize index table + // Multi-component base elements have to be thought of. + // 1. Vertices + unsigned total_index = 0; + for(unsigned comp = 0; comp < n_component_elements() ; ++comp) + { + for (unsigned local_index = 0 ; + local_index < base_element(comp).dofs_per_vertex ; + ++local_index) + { + system_to_component_table[total_index++] = pair + (comp, 0); + } + } + // 2. Lines + for(unsigned comp = 0; comp < n_component_elements() ; ++comp) + { + for (unsigned local_index = 0 ; + local_index < base_element(comp).dofs_per_vertex ; + ++local_index) + { + system_to_component_table[total_index++] = pair + (comp, 0); + } + } + // 3. Quads + for(unsigned comp = 0; comp < n_component_elements() ; ++comp) + { + for (unsigned local_index = 0 ; + local_index < base_element(comp).dofs_per_vertex ; + ++local_index) + { + system_to_component_table[total_index++] = pair + (comp, 0); + } + } + + // distribute the matrices of the base // finite element to the matrices of // this object @@ -67,7 +108,8 @@ FESystem<1>::multiply_dof_numbers (const FiniteElementData<1> &fe_data, template <> FiniteElementData<2> FESystem<2>::multiply_dof_numbers (const FiniteElementData<2> &fe_data, - const unsigned int N) { + const unsigned int N) +{ return FiniteElementData<2> (fe_data.dofs_per_vertex * N, fe_data.dofs_per_line * N, fe_data.dofs_per_quad * N, @@ -79,12 +121,14 @@ FESystem<2>::multiply_dof_numbers (const FiniteElementData<2> &fe_data, +/* template double FESystem::shape_value (const unsigned int i, - const Point &p) const { + const Point &p) const +{ Assert((ishape_value (i / n_sub_elements, p); +// return base_element->shape_value (i / n_sub_elements, p); }; @@ -95,7 +139,7 @@ FESystem::shape_grad (const unsigned int i, const Point &p) const { Assert((ishape_grad (i / n_sub_elements, p); +// return base_element->shape_grad (i / n_sub_elements, p); }; @@ -106,7 +150,7 @@ FESystem::shape_grad_grad (const unsigned int i, const Point &p) const { Assert((ishape_grad_grad (i / n_sub_elements, p); +// return base_element->shape_grad_grad (i / n_sub_elements, p); }; @@ -245,7 +289,7 @@ void FESystem::get_normal_vectors (const DoFHandler::cell_iterator &ce base_element->get_normal_vectors (cell, face_no, subface_no, unit_points, normal_vectors); }; - +*/