From: Wolfgang Bangerth Date: Mon, 27 Apr 1998 09:27:36 +0000 (+0000) Subject: Small changes and doc updates. X-Git-Tag: v8.0.0~23059 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0fe3a9e25bafee20d961cd5f2559d79dfddc9c50;p=dealii.git Small changes and doc updates. git-svn-id: https://svn.dealii.org/trunk@195 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_handler.h b/deal.II/deal.II/include/dofs/dof_handler.h index 741b41aa9b..08fcc2d66b 100644 --- a/deal.II/deal.II/include/dofs/dof_handler.h +++ b/deal.II/deal.II/include/dofs/dof_handler.h @@ -5,7 +5,6 @@ /*---------------------------- dof.h ---------------------------*/ #include -#include #include @@ -27,6 +26,8 @@ template class TriaActiveIterator; template class Triangulation; +template class FiniteElementBase; + class dVector; class dSMatrix; class dSMatrixStruct; @@ -468,7 +469,7 @@ class DoFHandler : public DoFDimensionInfo { * don't need them after calling this * function, or if so store them. */ - void distribute_dofs (const FiniteElement &); + void distribute_dofs (const FiniteElementBase &); /** * Renumber the degrees of freedom according @@ -999,7 +1000,7 @@ class DoFHandler : public DoFDimensionInfo { * Return a constant reference to the * selected finite element object. */ - const FiniteElement & get_selected_fe () const; + const FiniteElementBase & get_selected_fe () const; /** * Return a constant reference to the @@ -1065,18 +1066,24 @@ class DoFHandler : public DoFDimensionInfo { * Reserve enough space in the * #levels[]# objects to store the * numbers of the degrees of freedom - * needed for the given element. + * needed for the given element. The + * given element is that one which + * was selected when calling + * #distribute_dofs# the last time. */ - void reserve_space (const FiniteElement &); + void reserve_space (); /** * Distribute dofs on the given cell, * with new dofs starting with index * #next_free_dof#. Return the next - * unused index number. + * unused index number. The finite + * element used is the one given to + * #distribute_dofs#, which is copied + * to #selected_fe#. + * */ int distribute_dofs_on_cell (active_cell_iterator &cell, - const FiniteElement &fe, unsigned int next_free_dof); /** @@ -1121,15 +1128,13 @@ class DoFHandler : public DoFDimensionInfo { /** * Store a copy of the finite element given * latest for the distribution of dofs. In - * fact, since the FE class itself has not - * much functionality, this object only + * fact, since the FE base class itself has + * not much functionality, this object only * stores numbers such as the number of - * dofs per line etc. Calling any of the - * more specific functions will result in - * an error about calling a pure virtual - * function. + * dofs per line, but also restriction + * and prolongation matrices, etc. */ - FiniteElement *selected_fe; + const FiniteElementBase *selected_fe; /** * Store the number of dofs created last diff --git a/deal.II/deal.II/include/grid/tria_boundary.h b/deal.II/deal.II/include/grid/tria_boundary.h index 4df6ee48e2..ae9f05541d 100644 --- a/deal.II/deal.II/include/grid/tria_boundary.h +++ b/deal.II/deal.II/include/grid/tria_boundary.h @@ -56,6 +56,11 @@ class Boundary { Specialisation of \Ref{Boundary}, which places the new point right into the middle of the given points. The middle is defined as the arithmetic mean of the points. + + This class does not really describe a boundary in the usual sense. By + placing new points in teh middle of old ones, it rather assumes that the + boundary of the domain is given by the polygon/polyhedron defined by the + boundary of the initial coarse triangulation. */ template class StraightBoundary : public Boundary { diff --git a/deal.II/deal.II/include/numerics/assembler.h b/deal.II/deal.II/include/numerics/assembler.h index 6ca23f3edd..20cca022bb 100644 --- a/deal.II/deal.II/include/numerics/assembler.h +++ b/deal.II/deal.II/include/numerics/assembler.h @@ -109,7 +109,7 @@ struct AssemblerData { dVector &rhs_vector, const Quadrature &quadrature, const FiniteElement &fe, - const UpdateFields &update_flags, + const UpdateFlags &update_flags, const Boundary &boundary); /** @@ -155,7 +155,7 @@ struct AssemblerData { * FEValues object need to be reinitialized * on each cell. */ - const UpdateFields update_flags; + const UpdateFlags update_flags; /** * Store a pointer to the object describing diff --git a/deal.II/deal.II/include/numerics/base.h b/deal.II/deal.II/include/numerics/base.h index b4cf45b036..03d682b35a 100644 --- a/deal.II/deal.II/include/numerics/base.h +++ b/deal.II/deal.II/include/numerics/base.h @@ -8,6 +8,7 @@ #include #include #include +#include #include @@ -25,7 +26,6 @@ template class StraightBoundary; - /** Denote which norm/integral is to be computed. The following possibilities are implemented: @@ -83,8 +83,7 @@ enum NormType { since for higher order finite elements, we may be tempted to use curved faces of cells for better approximation of the boundary. In this case, the transformation from the unit cell to the real cell requires knowledge of - the exact boundary of the domain. By default it is assumed that the - boundary be approximated by straight segments. + the exact boundary of the domain. {\bf Solving} @@ -106,7 +105,7 @@ enum NormType { Usually, all other boundary conditions, such as inhomogeneous Neumann values or mixed boundary conditions are handled in the weak formulation. No attempt - is made to include this into the process of assemblage therefore. + is made to include these into the process of assemblage therefore. The inclusion into the assemblage process is as follows: when the matrix and vectors are set up, a list of nodes subject to dirichlet bc is made and @@ -165,6 +164,14 @@ enum NormType { This can be done by overloading the virtual function #make_boundary_value_list# which must return a list of boundary dofs and their corresponding values. + + You should be aware that the boundary function may be evaluated at nodes + on the interior of faces. These, however, need not be on the true + boundary, but rather are on the approximation of the boundary represented + by teh mapping of the unit cell to the real cell. Since this mapping will + in most cases not be the exact one at the face, the boundary function is + evaluated at points which are not on the boundary and you should make + sure that the returned values are reasonable in some sense anyway. {\bf Computing errors} @@ -301,7 +308,7 @@ class ProblemBase { virtual void assemble (const Equation &equation, const Quadrature &q, const FiniteElement &fe, - const UpdateFields &update_flags, + const UpdateFlags &update_flags, const DirichletBC &dirichlet_bc = DirichletBC(), const Boundary &boundary = StraightBoundary());