]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
check new update_flags for higher derivatives and andd new functions and
authorjanssen <janssen@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 26 May 2012 22:26:25 +0000 (22:26 +0000)
committerjanssen <janssen@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 26 May 2012 22:26:25 +0000 (22:26 +0000)
variables for higher derivatives

git-svn-id: https://svn.dealii.org/branches/branch_higher_derivatives@25548 0785d39b-7218-0410-832d-ea1e28bc413d

12 files changed:
deal.II/include/deal.II/base/polynomial_space.h
deal.II/include/deal.II/base/tensor_product_polynomials.h
deal.II/include/deal.II/fe/fe.h
deal.II/include/deal.II/fe/fe_poly.h
deal.II/include/deal.II/fe/fe_poly.templates.h
deal.II/include/deal.II/fe/fe_poly_tensor.h
deal.II/include/deal.II/fe/fe_update_flags.h
deal.II/source/base/polynomial_space.cc
deal.II/source/base/tensor_product_polynomials.cc
deal.II/source/fe/fe_poly_tensor.cc
tests/fe/update_flags_higher_derivatives.cc [new file with mode: 0644]
tests/fe/update_flags_higher_derivatives/cmp/generic [new file with mode: 0644]

index 1caa9e997a6892437832f73dfdf9aa35b31f23a5..a67996ca297ff26bb1dd77bb204c0168d92bff25 100644 (file)
@@ -21,6 +21,7 @@
 #include <deal.II/base/smartpointer.h>
 
 #include <vector>
+#include <boost/any.hpp>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -186,6 +187,10 @@ class PolynomialSpace
     Tensor<2,dim> compute_hessian (const unsigned int i,
                                      const Point<dim> &p) const;
 
+    boost::any compute_nth_derivative (const unsigned int i,
+                                     const Point<dim> &p,
+                                     const unsigned int nth_derivative) const;
+
                                      /**
                                       * Return the number of
                                       * polynomials spanning the space
index 8108852dc9627dcd4526ed4f43a76f6aadc2d095..73c5f21eedf596f981b1de30dd249d3afa9994b6 100644 (file)
@@ -20,6 +20,7 @@
 #include <deal.II/base/polynomial.h>
 
 #include <vector>
+#include <boost/any.hpp>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -218,6 +219,10 @@ class TensorProductPolynomials
     Tensor<2,dim> compute_hessian (const unsigned int i,
                                      const Point<dim> &p) const;
 
+    boost::any compute_nth_derivative (const unsigned int i,
+                                        const Point<dim> &,
+                                        const unsigned int nth_derivative) const;
+
                                      /**
                                       * Returns the number of tensor
                                       * product polynomials. For <i>n</i>
@@ -450,6 +455,9 @@ class AnisotropicPolynomials
     Tensor<2,dim> compute_hessian (const unsigned int i,
                                      const Point<dim> &p) const;
 
+    template<unsigned int nth_derivative>
+    Tensor<nth_derivative, dim> compute_nth_derivatives (const unsigned int i,
+                                                         const Point<dim> &) const;
                                      /**
                                       * Returns the number of tensor
                                       * product polynomials. It is the
index cd32d2eca10d07fa78216f7e5a1d06d1f3f5ffea..b2eb49ffce65e4f495e96c8b0b4c99d44811b859 100644 (file)
@@ -628,11 +628,12 @@ class FiniteElement : public Subscriptor,
                                                      const unsigned int component) const;
 
     /**
-     *
-     * @param i
+     * For computing the
+     * @param nth_derivative derivatives of the
+     * @param i th shape_function at point
      * @param p
-     * @param nth_derivative
-     * @return ...say that the any-type can only store Tensor<nth_derivative,dim>
+     * @param nth_derivative. The return value
+     * @return can only store Tensor<nth_derivative,dim>.
      */
     virtual
     boost::any
@@ -640,6 +641,15 @@ class FiniteElement : public Subscriptor,
                           const Point<dim> &p,
                           const unsigned int nth_derivative) const;
 
+    /**
+     * For computing the
+     * @param nth_derivative derivatives of the
+     * @param i th shape_function in the
+     * @param component component at point
+     * @param p
+     * @param nth_derivative. The return value
+     * @return can only store Tensor<nth_derivative,dim>.
+     */
     virtual
     boost::any
     shape_nth_derivative_component (const unsigned int i,
index ca0da2666ecb2247b9a86039b41d35f20ce6af4a..2f95eaa759fccf531012eaeb835a3f39ef2fb3a3 100644 (file)
@@ -39,6 +39,10 @@ DEAL_II_NAMESPACE_OPEN
  *
  * Tensor<2,dim> compute_hessian (const unsigned int i,
  *                                  const Point<dim> &p) const;
+ *
+ * boost::any compute_nth_derivative (const unsigned int i,
+ *                                  const Point<dim> &p,
+ *                                  const unsigned int nth_derivative) const;
  * @endcode
  * Example classes are TensorProductPolynomials, PolynomialSpace or
  * PolynomialsP.
@@ -194,6 +198,50 @@ class FE_Poly : public FiniteElement<dim,spacedim>
     virtual Tensor<2,dim> shape_hessian_component (const unsigned int i,
                                                      const Point<dim> &p,
                                                      const unsigned int component) const;
+
+                                     /**
+                                      * Return the tensor of the 
+                                      * @param nth_derivative th
+                                      * derivatives of the
+                                      * @param i th shape function at
+                                      * point 
+                                      * @param p on the unit
+                                      * cell. See the
+                                      * FiniteElement base class
+                                      * for more information about the
+                                      * semantics of this function.
+                                      */
+    virtual
+    boost::any
+    shape_nth_derivative (const unsigned int i,
+                          const Point<dim> &p,
+                          const unsigned int nth_derivative) const;
+
+                                     /**
+                                      * Return the tensor of the 
+                                      * @param nth_derivative th
+                                      * derivatives of the
+                                      * @param i th shape function at
+                                      * point 
+                                      * @param p on the unit
+                                      * cell. See the
+                                      * FiniteElement base class
+                                      * for more information about the
+                                      * semantics of this function.
+                                      *
+                                      * Since this element is scalar,
+                                      * the returned value is the same
+                                      * as if the function without the
+                                      * @param component suffix
+                                      * were called, provided that the
+                                      * specified component is zero.
+                                      */
+    virtual
+    boost::any
+    shape_nth_derivative_component (const unsigned int i,
+                                    const Point<dim> &p,
+                                    const unsigned int component,
+                                    const unsigned int nth_derivative) const;
   protected:
 
     virtual
@@ -368,6 +416,14 @@ class FE_Poly : public FiniteElement<dim,spacedim>
                                           * visiting an actual cell.
                                           */
         std::vector<std::vector<Tensor<1,dim> > > shape_gradients;
+        std::vector<std::vector<Tensor<2,dim> > > shape_hessians;
+        std::vector<std::vector<Tensor<3,dim> > > shape_3rd_derivatives;
+        std::vector<std::vector<Tensor<4,dim> > > shape_4th_derivatives;
+        std::vector<std::vector<Tensor<5,dim> > > shape_5th_derivatives;
+        std::vector<std::vector<Tensor<6,dim> > > shape_6th_derivatives;
+        std::vector<std::vector<Tensor<7,dim> > > shape_7th_derivatives;
+        std::vector<std::vector<Tensor<8,dim> > > shape_8th_derivatives;
+        std::vector<std::vector<Tensor<9,dim> > > shape_9th_derivatives;
     };
 
                                      /**
index 4fe053bdf99fffb4c19a3e1b5a48989934ddeb4c..53852119491f6aa1b1fe456e22e93a28d4ecc67c 100644 (file)
@@ -112,6 +112,31 @@ FE_Poly<POLY,dim,spacedim>::shape_hessian_component (const unsigned int i,
 
 
 
+template <class POLY, int dim, int spacedim>
+boost::any
+FE_Poly<POLY,dim,spacedim>::shape_nth_derivative (const unsigned int i,
+                                                  const Point<dim> &p,
+                                                  const unsigned int nth_derivative) const
+{
+  Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+  return poly_space.compute_nth_derivative(i, p, nth_derivative);
+}
+
+
+
+template <class POLY, int dim, int spacedim>
+boost::any
+FE_Poly<POLY,dim,spacedim>::shape_nth_derivative_component (const unsigned int i,
+                                                            const Point<dim> &p,
+                                                            const unsigned int component,
+                                                            const unsigned int nth_derivative) const
+{
+  Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+  Assert (component == 0, ExcIndexRange (component, 0, 1));
+  return poly_space.compute_nth_derivative(i, p, nth_derivative);
+}
+
+
 //---------------------------------------------------------------------------
 // Auxiliary functions
 //---------------------------------------------------------------------------
@@ -143,6 +168,20 @@ FE_Poly<POLY,dim,spacedim>::update_each (const UpdateFlags flags) const
     out |= update_gradients | update_covariant_transformation;
   if (flags & update_hessians)
     out |= update_hessians | update_covariant_transformation;
+  if (flags & update_3rd_derivatives)
+    out |= update_3rd_derivatives | update_covariant_transformation;
+  if (flags & update_4th_derivatives)
+    out |= update_4th_derivatives | update_covariant_transformation;
+  if (flags & update_5th_derivatives)
+    out |= update_5th_derivatives | update_covariant_transformation;
+  if (flags & update_6th_derivatives)
+    out |= update_6th_derivatives | update_covariant_transformation;
+  if (flags & update_7th_derivatives)
+    out |= update_7th_derivatives | update_covariant_transformation;
+  if (flags & update_8th_derivatives)
+    out |= update_8th_derivatives | update_covariant_transformation;
+  if (flags & update_9th_derivatives)
+    out |= update_9th_derivatives | update_covariant_transformation;
   if (flags & update_cell_normal_vectors)
     out |= update_cell_normal_vectors | update_JxW_values;
 
@@ -180,30 +219,86 @@ FE_Poly<POLY,dim,spacedim>::get_data (const UpdateFlags      update_flags,
   std::vector<double> values(0);
   std::vector<Tensor<1,dim> > grads(0);
   std::vector<Tensor<2,dim> > hessians(0);
+  std::vector<Tensor<3,dim> > third_derivatives(0);
+  std::vector<Tensor<4,dim> > fourth_derivatives(0);
+  std::vector<Tensor<5,dim> > fifth_derivatives(0);
+  std::vector<Tensor<6,dim> > sixth_derivatives(0);
+  std::vector<Tensor<7,dim> > seventh_derivatives(0);
+  std::vector<Tensor<8,dim> > eighth_derivatives(0);
+  std::vector<Tensor<9,dim> > ninth_derivatives(0);
 
                                    // initialize fields only if really
                                    // necessary. otherwise, don't
                                    // allocate memory
   if (flags & update_values)
-    {
-      values.resize (this->dofs_per_cell);
-      data->shape_values.resize (this->dofs_per_cell,
-                                 std::vector<double> (n_q_points));
-    }
+  {
+    values.resize (this->dofs_per_cell);
+    data->shape_values.resize (this->dofs_per_cell,
+        std::vector<double> (n_q_points));
+  }
 
   if (flags & update_gradients)
-    {
-      grads.resize (this->dofs_per_cell);
-      data->shape_gradients.resize (this->dofs_per_cell,
-                                    std::vector<Tensor<1,dim> > (n_q_points));
-    }
+  {
+    grads.resize (this->dofs_per_cell);
+    data->shape_gradients.resize (this->dofs_per_cell,
+        std::vector<Tensor<1,dim> > (n_q_points));
+  }
 
-                                   // if second derivatives through
-                                   // finite differencing is required,
-                                   // then initialize some objects for
-                                   // that
   if (flags & update_hessians)
-    data->initialize_2nd (this, mapping, quadrature);
+  {
+    hessians.resize (this->dofs_per_cell);
+    data->shape_hessians.resize (this->dofs_per_cell,
+        std::vector<Tensor<2,dim> >(n_q_points));
+  }
+
+  if (flags & update_3rd_derivatives)
+  {
+    third_derivatives.resize (this->dofs_per_cell);
+    data->shape_3rd_derivatives.resize (this->dofs_per_cell,
+        std::vector<Tensor<3,dim> >(n_q_points));
+  }
+
+  if (flags & update_4th_derivatives)
+  {
+    fourth_derivatives.resize (this->dofs_per_cell);
+    data->shape_4th_derivatives.resize (this->dofs_per_cell,
+        std::vector<Tensor<4,dim> >(n_q_points));
+  }
+
+  if (flags & update_5th_derivatives)
+  {
+    fifth_derivatives.resize (this->dofs_per_cell);
+    data->shape_5th_derivatives.resize (this->dofs_per_cell,
+        std::vector<Tensor<5,dim> >(n_q_points));
+  }
+
+  if (flags & update_6th_derivatives)
+  {
+    sixth_derivatives.resize (this->dofs_per_cell);
+    data->shape_6th_derivatives.resize (this->dofs_per_cell,
+        std::vector<Tensor<6,dim> >(n_q_points));
+  }
+
+  if (flags & update_7th_derivatives)
+  {
+    seventh_derivatives.resize (this->dofs_per_cell);
+    data->shape_7th_derivatives.resize (this->dofs_per_cell,
+        std::vector<Tensor<7,dim> >(n_q_points));
+  }
+
+  if (flags & update_8th_derivatives)
+  {
+    eighth_derivatives.resize (this->dofs_per_cell);
+    data->shape_8th_derivatives.resize (this->dofs_per_cell,
+        std::vector<Tensor<8,dim> >(n_q_points));
+  }
+
+  if (flags & update_9th_derivatives)
+  {
+    ninth_derivatives.resize (this->dofs_per_cell);
+    data->shape_9th_derivatives.resize (this->dofs_per_cell,
+        std::vector<Tensor<9,dim> >(n_q_points));
+  }
 
                                    // next already fill those fields
                                    // of which we have information by
@@ -212,6 +307,7 @@ FE_Poly<POLY,dim,spacedim>::get_data (const UpdateFlags      update_flags,
                                    // unit cell, and need to be
                                    // transformed when visiting an
                                    // actual cell
+  //TODO[BJ] add higher derivatives here
   if (flags & (update_values | update_gradients))
     for (unsigned int i=0; i<n_q_points; ++i)
       {
index 065b4f1cef812976d9c85953d7eb0259d7ecadde..14a24361a8523c3708bd9f8e08cc4be709e4fa27 100644 (file)
@@ -163,6 +163,20 @@ class FE_PolyTensor : public FiniteElement<dim,spacedim>
     virtual Tensor<2,dim> shape_hessian_component (const unsigned int i,
                                                      const Point<dim> &p,
                                                      const unsigned int component) const;
+    
+                                     /**
+                                      * Since these elements are
+                                      * vector valued, an exception is
+                                      * thrown.
+                                      */
+    virtual boost::any shape_nth_derivative (const unsigned int  i,
+                                             const Point<dim> &p,
+                                             const unsigned int nth_derivative) const;
+
+    virtual boost::any shape_nth_derivative_component (const unsigned int i,
+                                                const Point<dim> &p,
+                                                const unsigned int component,
+                                                const unsigned int nth_derivative) const;
 
                                      /**
                                       * Given <tt>flags</tt>,
@@ -332,6 +346,16 @@ class FE_PolyTensor : public FiniteElement<dim,spacedim>
                                       * shape_hessian_component().
                                       */
     mutable std::vector<Tensor<3,dim> > cached_hessians;
+
+    mutable std::vector<std::vector<boost::any> > cached_nth_derivatives;
+
+    mutable std::vector<Tensor<4,dim> > cached_3rd_derivatives;
+    mutable std::vector<Tensor<5,dim> > cached_4th_derivatives;
+    mutable std::vector<Tensor<6,dim> > cached_5th_derivatives;
+    mutable std::vector<Tensor<7,dim> > cached_6th_derivatives;
+    mutable std::vector<Tensor<8,dim> > cached_7th_derivatives;
+    mutable std::vector<Tensor<9,dim> > cached_8th_derivatives;
+    mutable std::vector<Tensor<10,dim> > cached_9th_derivatives;
 };
 
 DEAL_II_NAMESPACE_CLOSE
index f20c662cf23e0bdbc060985c190797d0cdc6adb8..cc7bcb26b3adc63f10df7f6ef508fb485f33c71d 100644 (file)
@@ -322,7 +322,37 @@ enum UpdateFlags
       update_piola = update_volume_elements | update_contravariant_transformation
 };
 
-UpdateFlags update_nth_derivatives (const unsigned int nth_derivative);
+inline
+UpdateFlags update_nth_derivatives (const unsigned int nth_derivative)
+{
+  switch (nth_derivative)
+  {
+    case 0:
+      return update_values;
+    case 1:
+      return update_gradients;
+    case 2:
+      return update_hessians;
+    case 3:
+      return update_3rd_derivatives;
+    case 4:
+      return update_4th_derivatives;
+    case 5:
+      return update_5th_derivatives;
+    case 6:
+      return update_6th_derivatives;
+    case 7:
+      return update_7th_derivatives;
+    case 8:
+      return update_8th_derivatives;
+    case 9:
+      return update_9th_derivatives;
+    default:
+      Assert (nth_derivative<10, ExcNotImplemented());
+  } 
+  return UpdateFlags();
+}
+
 
 
 /**
@@ -338,6 +368,13 @@ STREAM& operator << (STREAM& s, UpdateFlags u)
   if (u & update_values)                       s << "values|";
   if (u & update_gradients)                    s << "gradients|";
   if (u & update_hessians)                     s << "hessians|";
+  if (u & update_3rd_derivatives)              s << "3rd_derivatives|";
+  if (u & update_4th_derivatives)              s << "4th_derivatives|";
+  if (u & update_5th_derivatives)              s << "5th_derivatives|";
+  if (u & update_6th_derivatives)              s << "6th_derivatives|";
+  if (u & update_7th_derivatives)              s << "7th_derivatives|";
+  if (u & update_8th_derivatives)              s << "8th_derivatives|";
+  if (u & update_9th_derivatives)              s << "9th_derivatives|";
   if (u & update_quadrature_points)            s << "quadrature_points|";
   if (u & update_JxW_values)                   s << "JxW_values|";
   if (u & update_normal_vectors)               s << "normal_vectors|";
@@ -426,6 +463,111 @@ operator &= (UpdateFlags &f1, UpdateFlags f2)
   return f1;
 }
 
+inline
+UpdateFlags update_up_to_nth_derivatives (const unsigned int nth_derivative)
+{
+  UpdateFlags return_flags;
+  return_flags = update_values;
+  if(nth_derivative<1)
+    return return_flags;
+
+  return_flags |= update_gradients;
+  if(nth_derivative<2)
+    return return_flags;
+
+  return_flags |= update_hessians;
+  if(nth_derivative<3)
+    return return_flags;
+
+  return_flags |= update_3rd_derivatives;
+  if(nth_derivative<4)
+    return return_flags;
+
+  return_flags |= update_4th_derivatives;
+  if(nth_derivative<5)
+    return return_flags;
+
+  return_flags |= update_5th_derivatives;
+  if(nth_derivative<6)
+    return return_flags;
+
+  return_flags |= update_6th_derivatives;
+  if(nth_derivative<7)
+    return return_flags;
+
+  return_flags |= update_7th_derivatives;
+  if(nth_derivative<8)
+    return return_flags;
+
+  return_flags |= update_8th_derivatives;
+  if(nth_derivative<9)
+    return return_flags;
+
+  return_flags |= update_9th_derivatives;
+  if(nth_derivative<10)
+    return return_flags;
+
+  Assert (nth_derivative<10, ExcNotImplemented());
+  return UpdateFlags();
+}
+
+
+
+inline
+UpdateFlags update_derivatives (const unsigned int nth_derivative, const unsigned int mth_derivative)
+{
+  UpdateFlags return_flags = update_nth_derivatives(nth_derivative);
+  if(mth_derivative<1)
+    return return_flags;
+
+  if(nth_derivative<1)
+    return_flags |= update_gradients;
+  if(mth_derivative<2)
+    return return_flags;
+
+  if(nth_derivative<2)
+  return_flags |= update_hessians;
+  if(mth_derivative<3)
+    return return_flags;
+
+  if(nth_derivative<3)
+  return_flags |= update_3rd_derivatives;
+  if(mth_derivative<4)
+    return return_flags;
+
+  if(nth_derivative<4)
+  return_flags |= update_4th_derivatives;
+  if(mth_derivative<5)
+    return return_flags;
+
+  if(nth_derivative<5)
+  return_flags |= update_5th_derivatives;
+  if(mth_derivative<6)
+    return return_flags;
+
+  if(nth_derivative<6)
+  return_flags |= update_6th_derivatives;
+  if(mth_derivative<7)
+    return return_flags;
+
+  if(nth_derivative<7)
+  return_flags |= update_7th_derivatives;
+  if(mth_derivative<8)
+    return return_flags;
+
+  if(nth_derivative<8)
+  return_flags |= update_8th_derivatives;
+  if(mth_derivative<9)
+    return return_flags;
+
+  if(nth_derivative<9)
+  return_flags |= update_9th_derivatives;
+  if(mth_derivative<10)
+    return return_flags;
+
+  Assert (mth_derivative<10, ExcNotImplemented());
+  return UpdateFlags();
+}
 
 
 /**
index 1199671b9a6372b3ef7800c8afe4130169147486..df7dda0307f62f2545b8fa23a236f3294042ece6 100644 (file)
@@ -203,6 +203,46 @@ PolynomialSpace<dim>::compute_hessian (const unsigned int i,
 
 
 
+template <int dim>
+boost::any
+PolynomialSpace<dim>::compute_nth_derivative (const unsigned int i,
+                                         const Point<dim>  &p,
+                                         const unsigned int nth_derivative) const
+{
+  /*
+  unsigned int ix[dim];
+  compute_index(i,ix);
+
+  Tensor<2,dim> result;
+  for (unsigned int d=0; d<dim; ++d)
+    for (unsigned int d1=0; d1<dim; ++d1)
+      result[d][d1] = 1.;
+
+                                   // get value, first and second
+                                   // derivatives
+  std::vector<double> v(3);
+  for (unsigned int d=0; d<dim; ++d)
+    {
+      polynomials[ix[d]].value(p(d), v);
+      result[d][d] *= v[2];
+      for (unsigned int d1=0; d1<dim; ++d1)
+        {
+          if (d1 != d)
+            {
+              result[d][d1] *= v[1];
+              result[d1][d] *= v[1];
+              for (unsigned int d2=0; d2<dim; ++d2)
+                if (d2 != d)
+                  result[d1][d2] *= v[0];
+            }
+        }
+    }
+  return result;
+  */
+  return boost::any();
+}
+
+
 
 template <int dim>
 void
index 62d0c3d72372e4f19a84961acb78cbb29f8e722e..7b052f0f6730a9db278ba60bd101beb1c815894e 100644 (file)
@@ -206,6 +206,52 @@ TensorProductPolynomials<dim>::compute_hessian (const unsigned int i,
   return hessian;
 }
 
+template <int dim>
+boost::any
+TensorProductPolynomials<dim>::compute_nth_derivative (const unsigned int i,
+                                                  const Point<dim> &p,
+                                                  const unsigned int nth_derivative) const
+{
+  /*
+  unsigned int indices[dim];
+  compute_index (i, indices);
+
+  double v [dim][3];
+  {
+    std::vector<double> tmp (3);
+    for (unsigned int d=0; d<dim; ++d)
+      {
+        polynomials[indices[d]].value (p(d), tmp);
+        v[d][0] = tmp[0];
+        v[d][1] = tmp[1];
+        v[d][2] = tmp[2];
+      }
+  }
+
+  Tensor<2,dim> hessian;
+  for (unsigned int d1=0; d1<dim; ++d1)
+    for (unsigned int d2=0; d2<dim; ++d2)
+      {
+        hessian[d1][d2] = 1.;
+        for (unsigned int x=0; x<dim; ++x)
+          {
+            unsigned int derivative=0;
+            if (d1==x || d2==x)
+              {
+                if (d1==d2)
+                  derivative=2;
+                else
+                  derivative=1;
+              }
+            hessian[d1][d2] *= v[x][derivative];
+          }
+      }
+
+  return hessian;
+      */
+  return boost::any();
+}
+
 
 
 
index f19f4d788cdbe99f77134fdfe6bbb57615f3f82c..63de560b0d85febeeafeb118e2fa8e86c54610fa 100644 (file)
@@ -236,6 +236,45 @@ FE_PolyTensor<POLY,dim,spacedim>::shape_hessian_component (const unsigned int i,
 }
 
 
+template <class POLY, int dim, int spacedim>
+boost::any
+FE_PolyTensor<POLY,dim,spacedim>::shape_nth_derivative (const unsigned int, const Point<dim> &, const unsigned int) const
+{
+  typedef    FiniteElement<dim,spacedim> FEE;
+  Assert(false, typename FEE::ExcFENotPrimitive());
+  return Tensor<2,dim>();
+}
+
+
+
+template <class POLY, int dim, int spacedim>
+boost::any
+FE_PolyTensor<POLY,dim,spacedim>::shape_nth_derivative_component (const unsigned int i,
+                                                    const Point<dim> &p,
+                                                    const unsigned int component,
+                                                    const unsigned int nth_derivative) const
+{
+  Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+  Assert (component < dim, ExcIndexRange (component, 0, dim));
+
+  if (cached_point != p || cached_nth_derivatives[nth_derivative].size() == 0)
+    {
+      cached_point = p;
+      cached_nth_derivatives[nth_derivative].resize(poly_space.n());
+      //poly_space.compute(p, cached_values, cached_grads, cached_hessians);
+    }
+
+  boost::any s;
+  /*
+  if (inverse_node_matrix.n_cols() == 0)
+    return cached_nth_derivatives[nth_derivative][i][component];
+  else
+    for (unsigned int j=0;j<inverse_node_matrix.n_cols();++j)
+      s += inverse_node_matrix(i,j) * cached_nth_derivatives[nth_derivative][j][component];
+      */
+
+  return s;
+}
 
 //---------------------------------------------------------------------------
 // Data field initialization
diff --git a/tests/fe/update_flags_higher_derivatives.cc b/tests/fe/update_flags_higher_derivatives.cc
new file mode 100644 (file)
index 0000000..71035cf
--- /dev/null
@@ -0,0 +1,60 @@
+//----------------------------  update_flags_higher_derivatives.cc  ---------------------------
+//    update_flags_higher_derivatives.cc,v 1.3 2003/06/09 21:55:00 wolf Exp
+//    Version: 
+//
+//    Copyright (C) 2003, 2005 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    fuqher information on this license.
+//
+//----------------------------  update_flags_higher_derivatives.cc  ---------------------------
+
+
+// check update_nth_derivative
+// check update_up_to_nth_derivative
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/fe/fe_update_flags.h>
+
+#include <fstream>
+#include <iomanip>
+
+#define PRECISION 5
+
+void test ()
+{
+  //Note that derivatives are implemented up to order 9 only
+  for (unsigned int d=0; d<10; ++d)
+    deallog.get_file_stream() << update_nth_derivatives(d) << std::endl; 
+
+  deallog << std::endl;
+  for (unsigned int d=0; d<10; ++d)
+    deallog.get_file_stream() << update_up_to_nth_derivatives(d) << std::endl; 
+
+  deallog << std::endl;
+  for (unsigned int d=0; d<10; ++d)
+    for (unsigned int c=d; c<10; ++c)
+      deallog.get_file_stream() << update_derivatives(d,c) << std::endl; 
+}
+
+
+int
+main()
+{
+  std::ofstream logfile ("update_flags_higher_derivatives/output");
+  deallog << std::setprecision(PRECISION);
+  deallog << std::fixed;  
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test ();
+
+  return 0;
+}
+
+
+
diff --git a/tests/fe/update_flags_higher_derivatives/cmp/generic b/tests/fe/update_flags_higher_derivatives/cmp/generic
new file mode 100644 (file)
index 0000000..4f6ddfa
--- /dev/null
@@ -0,0 +1,78 @@
+
+ UpdateFlags|values|
+ UpdateFlags|gradients|
+ UpdateFlags|hessians|
+ UpdateFlags|3rd_derivatives|
+ UpdateFlags|4th_derivatives|
+ UpdateFlags|5th_derivatives|
+ UpdateFlags|6th_derivatives|
+ UpdateFlags|7th_derivatives|
+ UpdateFlags|8th_derivatives|
+ UpdateFlags|9th_derivatives|
+DEAL::
+ UpdateFlags|values|
+ UpdateFlags|values|gradients|
+ UpdateFlags|values|gradients|hessians|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|9th_derivatives|
+DEAL::
+ UpdateFlags|values|
+ UpdateFlags|values|gradients|
+ UpdateFlags|values|gradients|hessians|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|
+ UpdateFlags|values|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|9th_derivatives|
+ UpdateFlags|gradients|
+ UpdateFlags|gradients|hessians|
+ UpdateFlags|gradients|hessians|3rd_derivatives|
+ UpdateFlags|gradients|hessians|3rd_derivatives|4th_derivatives|
+ UpdateFlags|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|
+ UpdateFlags|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|
+ UpdateFlags|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|
+ UpdateFlags|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|
+ UpdateFlags|gradients|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|9th_derivatives|
+ UpdateFlags|hessians|
+ UpdateFlags|hessians|3rd_derivatives|
+ UpdateFlags|hessians|3rd_derivatives|4th_derivatives|
+ UpdateFlags|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|
+ UpdateFlags|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|
+ UpdateFlags|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|
+ UpdateFlags|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|
+ UpdateFlags|hessians|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|9th_derivatives|
+ UpdateFlags|3rd_derivatives|
+ UpdateFlags|3rd_derivatives|4th_derivatives|
+ UpdateFlags|3rd_derivatives|4th_derivatives|5th_derivatives|
+ UpdateFlags|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|
+ UpdateFlags|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|
+ UpdateFlags|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|
+ UpdateFlags|3rd_derivatives|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|9th_derivatives|
+ UpdateFlags|4th_derivatives|
+ UpdateFlags|4th_derivatives|5th_derivatives|
+ UpdateFlags|4th_derivatives|5th_derivatives|6th_derivatives|
+ UpdateFlags|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|
+ UpdateFlags|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|
+ UpdateFlags|4th_derivatives|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|9th_derivatives|
+ UpdateFlags|5th_derivatives|
+ UpdateFlags|5th_derivatives|6th_derivatives|
+ UpdateFlags|5th_derivatives|6th_derivatives|7th_derivatives|
+ UpdateFlags|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|
+ UpdateFlags|5th_derivatives|6th_derivatives|7th_derivatives|8th_derivatives|9th_derivatives|
+ UpdateFlags|6th_derivatives|
+ UpdateFlags|6th_derivatives|7th_derivatives|
+ UpdateFlags|6th_derivatives|7th_derivatives|8th_derivatives|
+ UpdateFlags|6th_derivatives|7th_derivatives|8th_derivatives|9th_derivatives|
+ UpdateFlags|7th_derivatives|
+ UpdateFlags|7th_derivatives|8th_derivatives|
+ UpdateFlags|7th_derivatives|8th_derivatives|9th_derivatives|
+ UpdateFlags|8th_derivatives|
+ UpdateFlags|8th_derivatives|9th_derivatives|
+ UpdateFlags|9th_derivatives|

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.