]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add unit_support_point() and face_unit_support_point() for cases where not all suppor...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 25 Nov 2002 17:30:05 +0000 (17:30 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 25 Nov 2002 17:30:05 +0000 (17:30 +0000)
git-svn-id: https://svn.dealii.org/trunk@6782 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_base.h
deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/source/fe/fe.cc
deal.II/deal.II/source/fe/fe_system.cc
deal.II/deal.II/source/numerics/vectors.cc
deal.II/doc/news/2002/c-3-4.html

index 3f4ca77227d802ac233ac1f2491c038a15ce708c..6d5f4db352c43932e1073878d2239e09c837d42b 100644 (file)
@@ -900,8 +900,57 @@ class FiniteElementBase : public Subscriptor,
                                      * then a call to the
                                      * @p{get_unit_support_points}
                                      * yields a non-empty array.
+                                     *
+                                     * The result may be false if an
+                                     * element is not defined by
+                                     * interpolating shape functions,
+                                     * for example by P-elements on
+                                     * quadrilaterals. It will
+                                     * usually only be true if the
+                                     * element constructs its shape
+                                     * functions by the requirement
+                                     * that they be one at a certain
+                                     * point and zero at all the
+                                     * points associated with the
+                                     * other shape functions.
+                                     *
+                                     * In composed elements (i.e. for
+                                     * the @ref{FESystem} class, the
+                                     * result will be true if all all
+                                     * the base elements have defined
+                                     * support points.
                                      */
     bool has_support_points () const;
+
+                                     /**
+                                      * Return the position of the
+                                      * support point of the
+                                      * @p{index}th shape function. If
+                                      * it does not exist, raise an
+                                      * exception.
+                                      *
+                                      * The default implementation
+                                      * simply returns the respective
+                                      * element from the array you get
+                                      * from
+                                      * @ref{get_unit_support_points},
+                                      * but derived elements may
+                                      * overload this function. In
+                                      * particular, note that the
+                                      * @ref{FESystem} class overloads
+                                      * it so that it can return the
+                                      * support points of individual
+                                      * base elements, of not all the
+                                      * base elements define support
+                                      * points. In this way, you can
+                                      * still ask for certain support
+                                      * points, even if
+                                      * @p{get_unit_support_points}
+                                      * only returns an empty array.
+                                      */
+    virtual
+    Point<dim>
+    unit_support_point (const unsigned int index) const;
     
                                     /**
                                      * Return the support points of
@@ -969,9 +1018,24 @@ class FiniteElementBase : public Subscriptor,
                                      * is true, then a call to the
                                      * @p{get_unit_support_points}
                                      * yields a non-empty array.
+                                     *
+                                     * For more information, see the
+                                     * documentation for the
+                                     * @ref{has_support_points}
+                                     * function.
                                      */
     bool has_face_support_points () const;
 
+                                     /**
+                                      * The function corresponding to
+                                      * the @p{unit_support_point}
+                                      * function, but for faces. See
+                                      * there for more information.
+                                      */
+    virtual
+    Point<dim-1>
+    unit_face_support_point (const unsigned int index) const;
+    
                                     /**
                                      * Return in which of the vector
                                      * components of this finite
index da53e2df548d96a7b1a6574ff1ac059a79e2e9a5..06a9bcb701a1207d4ef05b840cc0a18611b1d766 100644 (file)
@@ -312,7 +312,8 @@ class FESystem : public FiniteElement<dim>
     base_element (const unsigned int index) const;
     
                                     /**
-                                     * Check for non-zero values on a face.
+                                     * Check for non-zero values on a
+                                     * face.
                                      *
                                      * This function returns
                                      * @p{true}, if the shape
@@ -327,6 +328,24 @@ class FESystem : public FiniteElement<dim>
     virtual bool has_support_on_face (const unsigned int shape_index,
                                      const unsigned int face_index) const;
 
+                                     /**
+                                      * Implementation of the
+                                      * respective function in the
+                                      * base class.
+                                      */
+    virtual
+    Point<dim>
+    unit_support_point (const unsigned int index) const;
+    
+                                     /**
+                                      * Implementation of the
+                                      * respective function in the
+                                      * base class.
+                                      */
+    virtual
+    Point<dim-1>
+    unit_face_support_point (const unsigned int index) const;
+    
                                     /**
                                      * Determine an estimate for the
                                      * memory consumption (in bytes)
index 3029738e9f4e57c2cafe9facc222281c4f06f985..366f43d06daa505f7fa1c68765e5172e7b4018c2 100644 (file)
@@ -343,6 +343,19 @@ FiniteElementBase<dim>::has_support_points () const
 
 
 
+template <int dim>
+Point<dim>
+FiniteElementBase<dim>::unit_support_point (const unsigned index) const
+{
+  Assert (index < this->dofs_per_cell,
+          ExcIndexRange (index, 0, this->dofs_per_cell));
+  Assert (unit_support_points.size() == this->dofs_per_cell,
+          ExcFEHasNoSupportPoints ());
+  return unit_support_points[index];
+};
+
+
+
 template <int dim>
 const std::vector<Point<dim-1> > &
 FiniteElementBase<dim>::get_unit_face_support_points () const
@@ -368,6 +381,19 @@ FiniteElementBase<dim>::has_face_support_points () const
 
 
 
+template <int dim>
+Point<dim-1>
+FiniteElementBase<dim>::unit_face_support_point (const unsigned index) const
+{
+  Assert (index < this->dofs_per_face,
+          ExcIndexRange (index, 0, this->dofs_per_face));
+  Assert (unit_face_support_points.size() == this->dofs_per_face,
+          ExcFEHasNoSupportPoints ());
+  return unit_face_support_points[index];
+};
+
+
+
 template <int dim>
 unsigned int
 FiniteElementBase<dim>::memory_consumption () const
index 0757665a42aebd76c2e7440f7346e2ad9046b5d5..b2272e83405263c35562f3e52f80ef06568f9765 100644 (file)
@@ -2229,6 +2229,54 @@ FESystem<dim>::has_support_on_face (const unsigned int shape_index,
 
 
 
+template <int dim>
+Point<dim>
+FESystem<dim>::unit_support_point (const unsigned index) const
+{
+  Assert (index < this->dofs_per_cell,
+          ExcIndexRange (index, 0, this->dofs_per_cell));
+  Assert ((unit_support_points.size() == this->dofs_per_cell) ||
+          (unit_support_points.size() == 0),
+          typename FiniteElementBase<dim>::ExcFEHasNoSupportPoints ());
+
+                                   // let's see whether we have the
+                                   // information pre-computed
+  if (this->unit_support_points.size() != 0)
+    return this->unit_support_points[index];
+  else
+                                     // no. ask the base element
+                                     // whether it would like to
+                                     // provide this information
+    return (base_element(system_to_base_index(index).first.first)
+            .unit_support_point(system_to_base_index(index).second));
+};
+
+
+
+template <int dim>
+Point<dim-1>
+FESystem<dim>::unit_face_support_point (const unsigned index) const
+{
+  Assert (index < this->dofs_per_face,
+          ExcIndexRange (index, 0, this->dofs_per_face));
+  Assert ((unit_face_support_points.size() == this->dofs_per_face) ||
+          (unit_face_support_points.size() == 0),
+          typename FiniteElementBase<dim>::ExcFEHasNoSupportPoints ());
+
+                                   // let's see whether we have the
+                                   // information pre-computed
+  if (this->unit_face_support_points.size() != 0)
+    return this->unit_face_support_points[index];
+  else
+                                     // no. ask the base element
+                                     // whether it would like to
+                                     // provide this information
+    return (base_element(face_system_to_base_index(index).first.first)
+            .unit_face_support_point(face_system_to_base_index(index).second));
+};
+
+
+
 
 template <int dim>
 unsigned int
index d4c4ae8a0d598f26c7e70f4b0ad6dcbfc56b6b5c..16acb5cf98a9c3c2f720ca5a4dcfe0cb13052afd 100644 (file)
@@ -764,11 +764,12 @@ VectorTools::interpolate_boundary_values (const Mapping<1>              &mapping
 
 template <int dim>
 void
-VectorTools::interpolate_boundary_values (const Mapping<dim>            &mapping,
-                                         const DoFHandler<dim>         &dof,
-                                         const typename FunctionMap<dim>::type &function_map,
-                                         std::map<unsigned int,double> &boundary_values,
-                                         const std::vector<bool>       &component_mask_)
+VectorTools::
+interpolate_boundary_values (const Mapping<dim>            &mapping,
+                             const DoFHandler<dim>         &dof,
+                             const typename FunctionMap<dim>::type &function_map,
+                             std::map<unsigned int,double> &boundary_values,
+                             const std::vector<bool>       &component_mask_)
 {
                                   // if for whatever reason we were
                                   // passed an empty map, return
@@ -816,14 +817,36 @@ VectorTools::interpolate_boundary_values (const Mapping<dim>            &mapping
                                   // support points. this wil be used
                                   // to obtain the quadrature points
                                   // on the real cell's face
-  const std::vector<Point<dim-1> >
-    unit_support_points = fe.get_unit_face_support_points();
+  std::vector<Point<dim-1> >
+    unit_support_points = fe.get_unit_face_support_points();
   
                                   // check whether there are support
-                                  // points on the face, if not, then
-                                  // this FE does not allow to be
-                                  // used in this function
-  Assert (unit_support_points.size() != 0, ExcNonInterpolatingFE());
+                                  // points on the face. if not, then
+                                  // we should try a more clever
+                                  // way. the idea is that a finite
+                                  // element may not offer support
+                                  // points for all its shape
+                                  // functions, but maybe only
+                                  // some. if it offers support
+                                  // points for the components we are
+                                  // interested in in this function,
+                                  // then that's fine. if not, the
+                                  // function we call in the finite
+                                  // element will raise an
+                                  // exception. the support points
+                                  // for the other shape functions
+                                  // are left uninitialized (well,
+                                  // initialized by the default
+                                  // constructor), since we don't
+                                  // need them anyway.
+  if (unit_support_points.size() == 0)
+    {
+      unit_support_points.resize (fe.dofs_per_face);
+      for (unsigned int i=0; i<fe.dofs_per_face; ++i)
+        if (component_mask[fe.face_system_to_component_index(i).first]
+            == true)
+          unit_support_points[i] = fe.unit_face_support_point(i);
+    };
 
   Quadrature<dim-1> aux_quad (unit_support_points);
   FEFaceValues<dim> fe_values (mapping, fe, aux_quad, update_q_points);
index 6b5393de423290a700f484b8387ecb9698ee7bec..dc5f074c5e6562d9f4264079a3810d66177f23da 100644 (file)
@@ -467,6 +467,19 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p> 
+       New: For finite element classes, the functions
+       <code class="member">unit_support_point</code> and
+       <code class="member">unit_face_support_point</code> return the position
+       of an individual support point. This is necessary when you want to get
+       information about the support points of certain components in a composed
+       finite element, where not all components provide support points, and the
+       composed element thus does not fill the array the 
+       <code class="member">get_unit_support_points</code> function returns.
+       <br>
+       (WB 2002/11/23)
+       </p>
+
   <li> <p> 
        Fixed: Vectors could not be given as input and output vectors to the
        <code class="class">SolutionTransfer</code> class at the same time, but

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.