]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implementation of all FE::shape_value, shape_grad and shape_grad_grad functions.
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 8 May 2001 10:02:13 +0000 (10:02 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 8 May 2001 10:02:13 +0000 (10:02 +0000)
git-svn-id: https://svn.dealii.org/trunk@4560 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_base.h
deal.II/deal.II/include/fe/fe_dgq.h
deal.II/deal.II/include/fe/fe_q.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_dgq.cc
deal.II/deal.II/source/fe/fe_q.cc
deal.II/deal.II/source/fe/fe_system.cc

index 9bad222e319a508485c570577036db66053a97a5..3b1361b4322ff15e88f379a1d04c3887866dcf22 100644 (file)
@@ -379,6 +379,64 @@ class FiniteElementBase : public Subscriptor,
     FiniteElementBase (const FiniteElementData<dim> &fe_data,
                       const std::vector<bool> &restriction_is_additive_flags);
 
+                                    /**
+                                     * Return the value of the
+                                     * @p{i}th shape function at the
+                                     * point @p{p}.  @p{p} is a point
+                                     * on the reference element.
+                                     *
+                                     * An
+                                     * @p{ExcUnitShapeValuesDoNotExist}
+                                     * is thrown if the shape values
+                                     * of the @p{FiniteElement} under
+                                     * consideration depend on the
+                                     * shape of the cell in real
+                                     * space.
+                                     */
+    virtual double shape_value (const unsigned int i,
+                               const Point<dim> &p) const;
+    
+                                    /**
+                                     * Return the gradient of the
+                                     * @p{i}th shape function at the
+                                     * point @p{p}. @p{p} is a point
+                                     * on the reference element, and
+                                     * likewise the gradient is the
+                                     * gradient on the unit cell with
+                                     * respect to unit cell
+                                     * coordinates.
+                                     *
+                                     * An
+                                     * @p{ExcUnitShapeValuesDoNotExist}
+                                     * is thrown if the shape values
+                                     * of the @p{FiniteElement} under
+                                     * consideration depend on the
+                                     * shape of the cell in real
+                                     * space.
+                                     */
+    virtual Tensor<1,dim> shape_grad (const unsigned int  i,
+                                     const Point<dim>   &p) const;
+
+                                    /**
+                                     * Return the tensor of second
+                                     * derivatives of the @p{i}th
+                                     * shape function at point @p{p}
+                                     * on the unit cell. The
+                                     * derivatives are derivatives on
+                                     * the unit cell with respect to
+                                     * unit cell coordinates.
+                                     *
+                                     * An
+                                     * @p{ExcUnitShapeValuesDoNotExist}
+                                     * is thrown if the shape values
+                                     * of the @p{FiniteElement} under
+                                     * consideration depend on the
+                                     * shape of the cell in real
+                                     * space.
+                                     */
+    virtual Tensor<2,dim> shape_grad_grad (const unsigned int  i,
+                                          const Point<dim> &p) const;
+
                                     /**
                                      * Projection from a fine grid
                                      * space onto a coarse grid
@@ -663,6 +721,11 @@ class FiniteElementBase : public Subscriptor,
                                      */
     unsigned int memory_consumption () const;
 
+                                    /**
+                                     * Exception
+                                     */
+    DeclException0 (ExcUnitShapeValuesDoNotExist);
+
                                     /**
                                      * Exception
                                      */
index ac2b3307da8ee084c747b76c236691c280c3b004..86e8a5425f882e1720cc2045422779a592636d79 100644 (file)
@@ -43,7 +43,41 @@ class FE_DGQ : public FiniteElement<dim>
                                      * Destructor.
                                      */
     ~FE_DGQ ();
+    
+                                    /**
+                                     * Return the value of the
+                                     * @p{i}th shape function at the
+                                     * point @p{p}.  @p{p} is a point
+                                     * on the reference element.
+                                     */
+    virtual double shape_value (const unsigned int i,
+                               const Point<dim> &p) const;
+    
+                                    /**
+                                     * Return the gradient of the
+                                     * @p{i}th shape function at the
+                                     * point @p{p}. @p{p} is a point
+                                     * on the reference element, and
+                                     * likewise the gradient is the
+                                     * gradient on the unit cell with
+                                     * respect to unit cell
+                                     * coordinates.
+                                     */
+    virtual Tensor<1,dim> shape_grad (const unsigned int  i,
+                                     const Point<dim>   &p) const;
 
+                                    /**
+                                     * Return the tensor of second
+                                     * derivatives of the @p{i}th
+                                     * shape function at point @p{p}
+                                     * on the unit cell. The
+                                     * derivatives are derivatives on
+                                     * the unit cell with respect to
+                                     * unit cell coordinates.
+                                     */
+    virtual Tensor<2,dim> shape_grad_grad (const unsigned int  i,
+                                          const Point<dim> &p) const;
+    
                                     /**
                                      * Return the polynomial degree
                                      * of this finite element,
index 3b7bb883c446e0e73f4c66ab58eda1a6bfb9cf62..67b096c783e54d49a70260b6149e26ff2149330d 100644 (file)
@@ -248,11 +248,54 @@ class FE_Q : public FiniteElement<dim>
                                      * Destructor.
                                      */
     ~FE_Q ();
+    
+                                    /**
+                                     * Return the value of the
+                                     * @p{i}th shape function at the
+                                     * point @p{p}.  @p{p} is a point
+                                     * on the reference element.
+                                     */
+    virtual double shape_value (const unsigned int i,
+                               const Point<dim> &p) const;
+    
+                                    /**
+                                     * Return the gradient of the
+                                     * @p{i}th shape function at the
+                                     * point @p{p}. @p{p} is a point
+                                     * on the reference element, and
+                                     * likewise the gradient is the
+                                     * gradient on the unit cell with
+                                     * respect to unit cell
+                                     * coordinates.
+                                     */
+    virtual Tensor<1,dim> shape_grad (const unsigned int  i,
+                                     const Point<dim>   &p) const;
+
+                                    /**
+                                     * Return the tensor of second
+                                     * derivatives of the @p{i}th
+                                     * shape function at point @p{p}
+                                     * on the unit cell. The
+                                     * derivatives are derivatives on
+                                     * the unit cell with respect to
+                                     * unit cell coordinates.
+                                     */
+    virtual Tensor<2,dim> shape_grad_grad (const unsigned int  i,
+                                          const Point<dim> &p) const;
 
 //TODO:[RH] make get_renumber private or move it some other place    
                                     /**
                                      * Read-only access to the
                                      * renumber vector.
+                                     *
+                                     * This function shouldn't be
+                                     * used, i.e. the users code
+                                     * shouldn't rely on the actual
+                                     * numbering of the degrees of
+                                     * dreedom on each cell. This,
+                                     * because this internal
+                                     * numbering might change in
+                                     * future.
                                      */
     const std::vector<unsigned int> & get_renumber() const;
 
@@ -514,7 +557,15 @@ class FE_Q : public FiniteElement<dim>
                                      * shape function numbering.
                                      */
     std::vector<unsigned int> renumber;
-    
+
+                                    /**
+                                     * Inverse renumber
+                                     * vector. i.e. mapping from
+                                     * shape function numbering to
+                                     * lexicographic numbering.
+                                     */
+    std::vector<unsigned int> renumber_inverse;
+             
                                     /**
                                      * Mapping from lexicographic to
                                      * shape function numbering on first face.
index dbeb5bca16217353ad8f0a49ed68e835750946a8..2252f569912b5d804b20f38e7bc2098bf9484c59 100644 (file)
@@ -122,6 +122,66 @@ class FESystem : public FiniteElement<dim>
                                      */
     virtual ~FESystem ();
 
+                                    /**
+                                     * Return the value of the
+                                     * @p{i}th shape function at the
+                                     * point @p{p}.  @p{p} is a point
+                                     * on the reference element.
+                                     *
+                                     * An
+                                     * @p{ExcUnitShapeValuesDoNotExist}
+                                     * is thrown if the shape values
+                                     * of the @p{FiniteElement}
+                                     * (corresponding to the @p{i}th
+                                     * shape function) depend on the
+                                     * shape of the cell in real
+                                     * space.
+                                     */
+    virtual double shape_value (const unsigned int i,
+                               const Point<dim> &p) const;
+    
+                                    /**
+                                     * Return the gradient of the
+                                     * @p{i}th shape function at the
+                                     * point @p{p}. @p{p} is a point
+                                     * on the reference element, and
+                                     * likewise the gradient is the
+                                     * gradient on the unit cell with
+                                     * respect to unit cell
+                                     * coordinates.
+                                     *
+                                     * An
+                                     * @p{ExcUnitShapeValuesDoNotExist}
+                                     * is thrown if the shape values
+                                     * of the @p{FiniteElement}
+                                     * (corresponding to the @p{i}th
+                                     * shape function) depend on the
+                                     * shape of the cell in real
+                                     * space.
+                                     */
+    virtual Tensor<1,dim> shape_grad (const unsigned int  i,
+                                     const Point<dim>   &p) const;
+
+                                    /**
+                                     * Return the tensor of second
+                                     * derivatives of the @p{i}th
+                                     * shape function at point @p{p}
+                                     * on the unit cell. The
+                                     * derivatives are derivatives on
+                                     * the unit cell with respect to
+                                     * unit cell coordinates.
+                                     *
+                                     * An
+                                     * @p{ExcUnitShapeValuesDoNotExist}
+                                     * is thrown if the shape values
+                                     * of the @p{FiniteElement}
+                                     * (corresponding to the @p{i}th
+                                     * shape function) depend on the
+                                     * shape of the cell in real
+                                     * space.
+                                     */
+    virtual Tensor<2,dim> shape_grad_grad (const unsigned int  i,
+                                          const Point<dim> &p) const;
 
                                     /** 
                                      * Number of different base
index 5a2378127c932eca785d1674fc2a4828d003a86d..c3877ab1f28797760a7c4e164afc7a02a7860ac1 100644 (file)
@@ -167,6 +167,37 @@ FiniteElementBase<dim>::FiniteElementBase (const FiniteElementData<dim> &fe_data
 };
 
 
+template <int dim>
+double
+FiniteElementBase<dim>::shape_value (const unsigned int,
+                                    const Point<dim> &) const
+{
+  AssertThrow(false, ExcUnitShapeValuesDoNotExist());
+  return 0.;
+}
+
+
+
+template <int dim>
+Tensor<1,dim>
+FiniteElementBase<dim>::shape_grad (const unsigned int,
+                                   const Point<dim> &) const
+{
+  AssertThrow(false, ExcUnitShapeValuesDoNotExist());
+  return Tensor<1,dim> ();
+}
+
+
+
+template <int dim>
+Tensor<2,dim>
+FiniteElementBase<dim>::shape_grad_grad (const unsigned int,
+                                        const Point<dim> &) const
+{
+  AssertThrow(false, ExcUnitShapeValuesDoNotExist());
+  return Tensor<2,dim> ();
+}
+
 
 template <int dim>
 const FullMatrix<double> &
index e4fe8910e7cae205a39aa5d4fc3300cf518fd136..87a428ec69cae47afe1db693f968cb8ae1baae6c 100644 (file)
@@ -214,6 +214,34 @@ FE_DGQ<dim>::clone() const
 
 
 
+template <int dim>
+double
+FE_DGQ<dim>::shape_value (const unsigned int i,
+                         const Point<dim> &p) const
+{
+  return poly->compute_value(i, p);
+}
+
+
+
+template <int dim>
+Tensor<1,dim>
+FE_DGQ<dim>::shape_grad (const unsigned int i,
+                        const Point<dim> &p) const
+{
+  return poly->compute_grad(i, p);
+}
+
+
+
+template <int dim>
+Tensor<2,dim>
+FE_DGQ<dim>::shape_grad_grad (const unsigned int i,
+                             const Point<dim> &p) const
+{
+  return poly->compute_grad_grad(i, p);
+}
+
 
 //----------------------------------------------------------------------
 // Auxiliary functions
index 4b0828fa2eff4df9d3f86bc46748cf28622a2f86..f21d11d117f443b5d448b5dd3a1b462cb2e73cd2 100644 (file)
@@ -34,6 +34,7 @@ FE_Q<dim>::FE_Q (const unsigned int degree)
                                    std::vector<bool> (1,false)),
                degree(degree),
                renumber(dofs_per_cell, 0),
+               renumber_inverse(dofs_per_cell, 0),
                face_renumber(dofs_per_face, 0),
                poly(0)
 {
@@ -52,6 +53,11 @@ FE_Q<dim>::FE_Q (const unsigned int degree)
                                   // face function is empty
   build_renumbering (*this, degree, renumber);
   build_face_renumbering (degree, face_renumber);
+  
+                                  // build inverse of renumbering
+                                  // vector
+  for (unsigned int i=0; i<dofs_per_cell; ++i)
+    renumber_inverse[renumber[i]]=i;
 
                                   // copy constraint matrices if they
                                   // are defined. otherwise set them
@@ -596,6 +602,35 @@ FE_Q<dim>::clone() const
 
 
 
+template <int dim>
+double
+FE_Q<dim>::shape_value (const unsigned int i,
+                       const Point<dim> &p) const
+{
+  return poly->compute_value(renumber_inverse[i], p);
+}
+
+
+
+template <int dim>
+Tensor<1,dim>
+FE_Q<dim>::shape_grad (const unsigned int i,
+                      const Point<dim> &p) const
+{
+  return poly->compute_grad(renumber_inverse[i], p);
+}
+
+
+
+template <int dim>
+Tensor<2,dim>
+FE_Q<dim>::shape_grad_grad (const unsigned int i,
+                           const Point<dim> &p) const
+{
+  return poly->compute_grad_grad(renumber_inverse[i], p);
+}
+
+
 //----------------------------------------------------------------------
 // Auxiliary functions
 //----------------------------------------------------------------------
index 58b134bd04116ffa58e415dd0679fb4d970082e9..d9323270a3876454b09125c722fa60405936d8e1 100644 (file)
@@ -232,6 +232,52 @@ FESystem<dim>::clone() const
 
 
 
+template <int dim>
+double
+FESystem<dim>::shape_value (const unsigned int i,
+                           const Point<dim> &p) const
+{
+  Assert((i<dofs_per_cell), ExcIndexRange(i, 0, dofs_per_cell));
+
+  pair<unsigned,unsigned> comp = system_to_component_index(i);
+  
+  return base_element(component_to_base_table[comp.first])
+    .shape_value(comp.second, p);
+}
+
+
+
+template <int dim>
+Tensor<1,dim>
+FESystem<dim>::shape_grad (const unsigned int i,
+                          const Point<dim> &p) const
+{
+  Assert((i<dofs_per_cell), ExcIndexRange(i, 0, dofs_per_cell));
+
+  pair<unsigned,unsigned> comp = system_to_component_index(i);
+  
+  return base_element(component_to_base_table[comp.first])
+    .shape_grad(comp.second, p);
+}
+
+
+
+template <int dim>
+Tensor<2,dim>
+FESystem<dim>::shape_grad_grad (const unsigned int i,
+                               const Point<dim> &p) const
+{
+  Assert((i<dofs_per_cell), ExcIndexRange(i, 0, dofs_per_cell));
+
+
+  pair<unsigned,unsigned> comp = system_to_component_index(i);
+  
+  return base_element(component_to_base_table[comp.first])
+    .shape_grad_grad(comp.second, p);
+}
+
+
+
 template <int dim>
 UpdateFlags
 FESystem<dim>::update_once (const UpdateFlags flags) const

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.