]> https://gitweb.dealii.org/ - dealii.git/commitdiff
move ConstantFunction and ZeroFunction to Functions:: namespace
authorDenis Davydov <davydden@gmail.com>
Fri, 18 Aug 2017 13:44:35 +0000 (15:44 +0200)
committerDenis Davydov <davydden@gmail.com>
Sat, 19 Aug 2017 05:40:54 +0000 (07:40 +0200)
include/deal.II/base/function.h
include/deal.II/base/function.templates.h
source/base/function.inst.in

index a377407d4cfc0822beb23d4063e8e937d11b266d..2b618b8dd41307a16eb4210e6c3efcb79b8a61fd 100644 (file)
@@ -347,135 +347,158 @@ public:
 };
 
 
-
-/**
- * Provide a function which always returns zero. Obviously, also the derivates
- * of this function are zero. Also, it returns zero on all components in case
- * the function is not a scalar one, which can be obtained by passing the
- * constructor the appropriate number of components.
- *
- * This function is of use when you want to implement homogeneous boundary
- * conditions, or zero initial conditions.
- *
- * @ingroup functions
- * @author Wolfgang Bangerth, 1998, 1999
- */
-template <int dim, typename Number=double>
-class ZeroFunction : public Function<dim, Number>
+namespace Functions
 {
-public:
-  /**
-   * Constructor. The number of components is preset to one.
-   */
-  ZeroFunction (const unsigned int n_components = 1);
+
 
   /**
-   * Virtual destructor; absolutely necessary in this case.
+   * Provide a function which always returns zero. Obviously, also the derivates
+   * of this function are zero. Also, it returns zero on all components in case
+   * the function is not a scalar one, which can be obtained by passing the
+   * constructor the appropriate number of components.
    *
+   * This function is of use when you want to implement homogeneous boundary
+   * conditions, or zero initial conditions.
+   *
+   * @ingroup functions
+   * @author Wolfgang Bangerth, 1998, 1999
    */
-  virtual ~ZeroFunction ();
+  template <int dim, typename Number=double>
+  class ZeroFunction : public Function<dim, Number>
+  {
+  public:
+    /**
+     * Constructor. The number of components is preset to one.
+     */
+    ZeroFunction (const unsigned int n_components = 1);
 
-  virtual Number value (const Point<dim>   &p,
-                        const unsigned int  component) const;
+    /**
+     * Virtual destructor; absolutely necessary in this case.
+     *
+     */
+    virtual ~ZeroFunction ();
 
-  virtual void vector_value (const Point<dim> &p,
-                             Vector<Number>   &return_value) const;
+    virtual Number value (const Point<dim>   &p,
+                          const unsigned int  component = 0) const;
 
-  virtual void value_list (const std::vector<Point<dim> > &points,
-                           std::vector<Number>            &values,
-                           const unsigned int              component = 0) const;
+    virtual void vector_value (const Point<dim> &p,
+                               Vector<Number>   &return_value) const;
 
-  virtual void vector_value_list (const std::vector<Point<dim> > &points,
-                                  std::vector<Vector<Number> >   &values) const;
+    virtual void value_list (const std::vector<Point<dim> > &points,
+                             std::vector<Number>            &values,
+                             const unsigned int              component = 0) const;
 
-  virtual Tensor<1,dim, Number> gradient (const Point<dim> &p,
-                                          const unsigned int component = 0) const;
+    virtual void vector_value_list (const std::vector<Point<dim> > &points,
+                                    std::vector<Vector<Number> >   &values) const;
 
-  virtual void vector_gradient (const Point<dim>            &p,
-                                std::vector<Tensor<1,dim, Number> > &gradients) const;
+    virtual Tensor<1,dim, Number> gradient (const Point<dim> &p,
+                                            const unsigned int component = 0) const;
 
-  virtual void gradient_list (const std::vector<Point<dim> > &points,
-                              std::vector<Tensor<1,dim, Number> >    &gradients,
-                              const unsigned int              component = 0) const;
+    virtual void vector_gradient (const Point<dim>            &p,
+                                  std::vector<Tensor<1,dim, Number> > &gradients) const;
 
-  virtual void vector_gradient_list (const std::vector<Point<dim> >            &points,
-                                     std::vector<std::vector<Tensor<1,dim, Number> > > &gradients) const;
-};
+    virtual void gradient_list (const std::vector<Point<dim> > &points,
+                                std::vector<Tensor<1,dim, Number> >    &gradients,
+                                const unsigned int              component = 0) const;
+
+    virtual void vector_gradient_list (const std::vector<Point<dim> >            &points,
+                                       std::vector<std::vector<Tensor<1,dim, Number> > > &gradients) const;
+  };
 
 
 
+  /**
+   * Provide a function which always returns the constant values handed to the
+   * constructor.
+   *
+   * Obviously, the derivates of this function are zero, which is why we derive
+   * this class from <tt>ZeroFunction</tt>: we then only have to overload the
+   * value functions, not all the derivatives. In some way, it would be more
+   * obvious to do the derivation in the opposite direction, i.e. let
+   * <tt>ZeroFunction</tt> be a more specialized version of
+   * <tt>ConstantFunction</tt>; however, this would be less efficient, since we
+   * could not make use of the fact that the function value of the
+   * <tt>ZeroFunction</tt> is known at compile time and need not be looked up
+   * somewhere in memory.
+   *
+   * @ingroup functions
+   * @author Wolfgang Bangerth, 1998, 1999, Lei Qiao, 2015
+   */
+  template <int dim, typename Number=double>
+  class ConstantFunction : public ZeroFunction<dim, Number>
+  {
+  public:
+    /**
+     * Constructor; set values of all components to the provided one. The
+     * default number of components is one.
+     */
+    ConstantFunction (const Number       value,
+                      const unsigned int n_components = 1);
+
+    /**
+     * Constructor; takes an <tt>std::vector<Number></tt> object as an argument.
+     * The number of components is determined by <tt>values.size()</tt>.
+     */
+    ConstantFunction (const std::vector<Number> &values);
+
+    /**
+     * Constructor; takes an <tt>Vector<Number></tt> object as an argument. The
+     * number of components is determined by <tt>values.size()</tt>.
+     */
+    ConstantFunction (const Vector<Number> &values);
+
+    /**
+     * Constructor; uses whatever stores in [begin_ptr, begin_ptr+n_components)
+     * to initialize a new object.
+     */
+    ConstantFunction (const Number *begin_ptr, const unsigned int n_components);
+
+    /**
+     * Virtual destructor; absolutely necessary in this case.
+     */
+    virtual ~ConstantFunction ();
+
+    virtual Number value (const Point<dim>   &p,
+                          const unsigned int  component = 0) const;
+
+    virtual void vector_value (const Point<dim> &p,
+                               Vector<Number>   &return_value) const;
+
+    virtual void value_list (const std::vector<Point<dim> > &points,
+                             std::vector<Number>            &return_values,
+                             const unsigned int              component = 0) const;
+
+    virtual void vector_value_list (const std::vector<Point<dim> > &points,
+                                    std::vector<Vector<Number> >   &return_values) const;
+
+    std::size_t memory_consumption () const;
+
+  protected:
+    /**
+     * Store the constant function value vector.
+     */
+    std::vector<Number> function_value_vector;
+  };
+
+}
+
 /**
  * Provide a function which always returns the constant values handed to the
  * constructor.
  *
- * Obviously, the derivates of this function are zero, which is why we derive
- * this class from <tt>ZeroFunction</tt>: we then only have to overload the
- * value functions, not all the derivatives. In some way, it would be more
- * obvious to do the derivation in the opposite direction, i.e. let
- * <tt>ZeroFunction</tt> be a more specialized version of
- * <tt>ConstantFunction</tt>; however, this would be less efficient, since we
- * could not make use of the fact that the function value of the
- * <tt>ZeroFunction</tt> is known at compile time and need not be looked up
- * somewhere in memory.
- *
- * @ingroup functions
- * @author Wolfgang Bangerth, 1998, 1999, Lei Qiao, 2015
+ * @deprecated use Functions::ConstantFunction instead.
  */
-template <int dim, typename Number=double>
-class ConstantFunction : public ZeroFunction<dim, Number>
-{
-public:
-  /**
-   * Constructor; set values of all components to the provided one. The
-   * default number of components is one.
-   */
-  ConstantFunction (const Number       value,
-                    const unsigned int n_components = 1);
-
-  /**
-   * Constructor; takes an <tt>std::vector<Number></tt> object as an argument.
-   * The number of components is determined by <tt>values.size()</tt>.
-   */
-  ConstantFunction (const std::vector<Number> &values);
-
-  /**
-   * Constructor; takes an <tt>Vector<Number></tt> object as an argument. The
-   * number of components is determined by <tt>values.size()</tt>.
-   */
-  ConstantFunction (const Vector<Number> &values);
-
-  /**
-   * Constructor; uses whatever stores in [begin_ptr, begin_ptr+n_components)
-   * to initialize a new object.
-   */
-  ConstantFunction (const Number *begin_ptr, const unsigned int n_components);
-
-  /**
-   * Virtual destructor; absolutely necessary in this case.
-   */
-  virtual ~ConstantFunction ();
+template<int dim, typename Number=double>
+using ConstantFunction DEAL_II_DEPRECATED = Functions::ConstantFunction<dim,Number>;
 
-  virtual Number value (const Point<dim>   &p,
-                        const unsigned int  component) const;
-
-  virtual void vector_value (const Point<dim> &p,
-                             Vector<Number>   &return_value) const;
-
-  virtual void value_list (const std::vector<Point<dim> > &points,
-                           std::vector<Number>            &return_values,
-                           const unsigned int              component = 0) const;
-
-  virtual void vector_value_list (const std::vector<Point<dim> > &points,
-                                  std::vector<Vector<Number> >   &return_values) const;
-
-  std::size_t memory_consumption () const;
+/**
+ * Provide a function which always returns zero.
+ *
+ * @deprecated use Functions::ZeroFunction instead.
+ */
+template<int dim, typename Number=double>
+using ZeroFunction DEAL_II_DEPRECATED = Functions::ZeroFunction<dim,Number>;
 
-protected:
-  /**
-   * Store the constant function value vector.
-   */
-  std::vector<Number> function_value_vector;
-};
 
 
 /**
index bd6deb9ab10db6426a79134e8a1f4fb38408f81d..14e87c92ca2963282c27336a78130d7c25e32823 100644 (file)
@@ -297,240 +297,248 @@ Function<dim, Number>::memory_consumption () const
 
 //---------------------------------------------------------------------------
 
-
-
-template <int dim, typename Number>
-ZeroFunction<dim, Number>::ZeroFunction (const unsigned int n_components)
-  :
-  Function<dim, Number> (n_components)
-{}
-
-
-template <int dim, typename Number>
-ZeroFunction<dim, Number>::~ZeroFunction ()
-{}
-
-
-template <int dim, typename Number>
-Number ZeroFunction<dim, Number>::value (const Point<dim> &,
-                                         const unsigned int) const
-{
-  return 0.;
-}
-
-
-template <int dim, typename Number>
-void ZeroFunction<dim, Number>::vector_value (const Point<dim> &,
-                                              Vector<Number>   &return_value) const
+namespace Functions
 {
-  Assert (return_value.size() == this->n_components,
-          ExcDimensionMismatch (return_value.size(), this->n_components));
-
-  std::fill (return_value.begin(), return_value.end(), 0.0);
-}
-
-
-template <int dim, typename Number>
-void ZeroFunction<dim, Number>::value_list (
-  const std::vector<Point<dim> > &points,
-  std::vector<Number>            &values,
-  const unsigned int              /*component*/) const
-{
-  (void)points;
-  Assert (values.size() == points.size(),
-          ExcDimensionMismatch(values.size(), points.size()));
 
-  std::fill (values.begin(), values.end(), 0.);
-}
+  template <int dim, typename Number>
+  ZeroFunction<dim, Number>::ZeroFunction (const unsigned int n_components)
+    :
+    Function<dim, Number> (n_components)
+  {}
+
+
+  template <int dim, typename Number>
+  ZeroFunction<dim, Number>::~ZeroFunction ()
+  {}
+
+
+  template <int dim, typename Number>
+  Number ZeroFunction<dim, Number>::value (const Point<dim> &,
+                                           const unsigned int) const
+  {
+    return 0.;
+  }
+
+
+  template <int dim, typename Number>
+  void ZeroFunction<dim, Number>::vector_value (const Point<dim> &,
+                                                Vector<Number>   &return_value) const
+  {
+    Assert (return_value.size() == this->n_components,
+            ExcDimensionMismatch (return_value.size(), this->n_components));
+
+    std::fill (return_value.begin(), return_value.end(), 0.0);
+  }
+
+
+  template <int dim, typename Number>
+  void ZeroFunction<dim, Number>::value_list (
+    const std::vector<Point<dim> > &points,
+    std::vector<Number>            &values,
+    const unsigned int              /*component*/) const
+  {
+    (void)points;
+    Assert (values.size() == points.size(),
+            ExcDimensionMismatch(values.size(), points.size()));
+
+    std::fill (values.begin(), values.end(), 0.);
+  }
+
+
+  template <int dim, typename Number>
+  void ZeroFunction<dim, Number>::vector_value_list (
+    const std::vector<Point<dim> > &points,
+    std::vector<Vector<Number> >   &values) const
+  {
+    Assert (values.size() == points.size(),
+            ExcDimensionMismatch(values.size(), points.size()));
+
+    for (unsigned int i=0; i<points.size(); ++i)
+      {
+        Assert (values[i].size() == this->n_components,
+                ExcDimensionMismatch(values[i].size(), this->n_components));
+        std::fill (values[i].begin(), values[i].end(), 0.);
+      };
+  }
+
+
+  template <int dim, typename Number>
+  Tensor<1,dim,Number> ZeroFunction<dim, Number>::gradient (const Point<dim> &,
+                                                            const unsigned int) const
+  {
+    return Tensor<1,dim,Number>();
+  }
+
+
+  template <int dim, typename Number>
+  void ZeroFunction<dim, Number>::vector_gradient (
+    const Point<dim> &,
+    std::vector<Tensor<1,dim,Number> > &gradients) const
+  {
+    Assert (gradients.size() == this->n_components,
+            ExcDimensionMismatch(gradients.size(), this->n_components));
+
+    for (unsigned int c=0; c<this->n_components; ++c)
+      gradients[c].clear ();
+  }
+
+
+  template <int dim, typename Number>
+  void ZeroFunction<dim, Number>::gradient_list (
+    const std::vector<Point<dim> >     &points,
+    std::vector<Tensor<1,dim,Number> > &gradients,
+    const unsigned int                  /*component*/) const
+  {
+    Assert (gradients.size() == points.size(),
+            ExcDimensionMismatch(gradients.size(), points.size()));
+
+    for (unsigned int i=0; i<points.size(); ++i)
+      gradients[i].clear ();
+  }
+
+
+  template <int dim, typename Number>
+  void ZeroFunction<dim, Number>::vector_gradient_list (
+    const std::vector<Point<dim> >                   &points,
+    std::vector<std::vector<Tensor<1,dim,Number> > > &gradients) const
+  {
+    Assert (gradients.size() == points.size(),
+            ExcDimensionMismatch(gradients.size(), points.size()));
+    for (unsigned int i=0; i<points.size(); ++i)
+      {
+        Assert (gradients[i].size() == this->n_components,
+                ExcDimensionMismatch(gradients[i].size(), this->n_components));
+        for (unsigned int c=0; c<this->n_components; ++c)
+          gradients[i][c].clear ();
+      };
+  }
 
-
-template <int dim, typename Number>
-void ZeroFunction<dim, Number>::vector_value_list (
-  const std::vector<Point<dim> > &points,
-  std::vector<Vector<Number> >   &values) const
-{
-  Assert (values.size() == points.size(),
-          ExcDimensionMismatch(values.size(), points.size()));
-
-  for (unsigned int i=0; i<points.size(); ++i)
-    {
-      Assert (values[i].size() == this->n_components,
-              ExcDimensionMismatch(values[i].size(), this->n_components));
-      std::fill (values[i].begin(), values[i].end(), 0.);
-    };
-}
-
-
-template <int dim, typename Number>
-Tensor<1,dim,Number> ZeroFunction<dim, Number>::gradient (const Point<dim> &,
-                                                          const unsigned int) const
-{
-  return Tensor<1,dim,Number>();
-}
-
-
-template <int dim, typename Number>
-void ZeroFunction<dim, Number>::vector_gradient (
-  const Point<dim> &,
-  std::vector<Tensor<1,dim,Number> > &gradients) const
-{
-  Assert (gradients.size() == this->n_components,
-          ExcDimensionMismatch(gradients.size(), this->n_components));
-
-  for (unsigned int c=0; c<this->n_components; ++c)
-    gradients[c].clear ();
-}
-
-
-template <int dim, typename Number>
-void ZeroFunction<dim, Number>::gradient_list (
-  const std::vector<Point<dim> >     &points,
-  std::vector<Tensor<1,dim,Number> > &gradients,
-  const unsigned int                  /*component*/) const
-{
-  Assert (gradients.size() == points.size(),
-          ExcDimensionMismatch(gradients.size(), points.size()));
-
-  for (unsigned int i=0; i<points.size(); ++i)
-    gradients[i].clear ();
-}
-
-
-template <int dim, typename Number>
-void ZeroFunction<dim, Number>::vector_gradient_list (
-  const std::vector<Point<dim> >                   &points,
-  std::vector<std::vector<Tensor<1,dim,Number> > > &gradients) const
-{
-  Assert (gradients.size() == points.size(),
-          ExcDimensionMismatch(gradients.size(), points.size()));
-  for (unsigned int i=0; i<points.size(); ++i)
-    {
-      Assert (gradients[i].size() == this->n_components,
-              ExcDimensionMismatch(gradients[i].size(), this->n_components));
-      for (unsigned int c=0; c<this->n_components; ++c)
-        gradients[i][c].clear ();
-    };
 }
 
 //---------------------------------------------------------------------------
 
-template <int dim, typename Number>
-ConstantFunction<dim, Number>::ConstantFunction (const Number value,
-                                                 const unsigned int n_components)
-  :
-  ZeroFunction<dim, Number> (n_components),
-  function_value_vector (n_components, value)
-{}
-
-template <int dim, typename Number>
-ConstantFunction<dim, Number>::
-ConstantFunction (const std::vector<Number> &values)
-  :
-  ZeroFunction<dim, Number> (values.size()),
-  function_value_vector (values)
-{}
-
-
-template <int dim, typename Number>
-ConstantFunction<dim, Number>::
-ConstantFunction (const Vector<Number> &values)
-  :
-  ZeroFunction<dim, Number> (values.size()),
-  function_value_vector (values.size())
-{
-  Assert (values.size() == function_value_vector.size(),
-          ExcDimensionMismatch (values.size(), function_value_vector.size()));
-  std::copy (values.begin(),values.end(),function_value_vector.begin());
-}
-
-
-template <int dim, typename Number>
-ConstantFunction<dim, Number>::
-ConstantFunction (const Number *begin_ptr, const unsigned int n_components)
-  :
-  ZeroFunction<dim, Number> (n_components),
-  function_value_vector (n_components)
-{
-  Assert (begin_ptr != nullptr, ExcMessage ("Null pointer encountered!"));
-  std::copy (begin_ptr, begin_ptr+n_components, function_value_vector.begin());
-}
-
-
-
-template <int dim, typename Number>
-ConstantFunction<dim, Number>::~ConstantFunction ()
-{
-  function_value_vector.clear();
-}
-
-
-template <int dim, typename Number>
-Number ConstantFunction<dim, Number>::value (const Point<dim> &,
-                                             const unsigned int component) const
-{
-  Assert (component < this->n_components,
-          ExcIndexRange (component, 0, this->n_components));
-  return function_value_vector[component];
-}
+namespace Functions
+{
+
+  template <int dim, typename Number>
+  ConstantFunction<dim, Number>::ConstantFunction (const Number value,
+                                                   const unsigned int n_components)
+    :
+    ZeroFunction<dim, Number> (n_components),
+    function_value_vector (n_components, value)
+  {}
+
+  template <int dim, typename Number>
+  ConstantFunction<dim, Number>::
+  ConstantFunction (const std::vector<Number> &values)
+    :
+    ZeroFunction<dim, Number> (values.size()),
+    function_value_vector (values)
+  {}
+
+
+  template <int dim, typename Number>
+  ConstantFunction<dim, Number>::
+  ConstantFunction (const Vector<Number> &values)
+    :
+    ZeroFunction<dim, Number> (values.size()),
+    function_value_vector (values.size())
+  {
+    Assert (values.size() == function_value_vector.size(),
+            ExcDimensionMismatch (values.size(), function_value_vector.size()));
+    std::copy (values.begin(),values.end(),function_value_vector.begin());
+  }
+
+
+  template <int dim, typename Number>
+  ConstantFunction<dim, Number>::
+  ConstantFunction (const Number *begin_ptr, const unsigned int n_components)
+    :
+    ZeroFunction<dim, Number> (n_components),
+    function_value_vector (n_components)
+  {
+    Assert (begin_ptr != nullptr, ExcMessage ("Null pointer encountered!"));
+    std::copy (begin_ptr, begin_ptr+n_components, function_value_vector.begin());
+  }
+
+
+
+  template <int dim, typename Number>
+  ConstantFunction<dim, Number>::~ConstantFunction ()
+  {
+    function_value_vector.clear();
+  }
+
+
+  template <int dim, typename Number>
+  Number ConstantFunction<dim, Number>::value (const Point<dim> &,
+                                               const unsigned int component) const
+  {
+    Assert (component < this->n_components,
+            ExcIndexRange (component, 0, this->n_components));
+    return function_value_vector[component];
+  }
+
+
+
+  template <int dim, typename Number>
+  void ConstantFunction<dim, Number>::vector_value (const Point<dim> &,
+                                                    Vector<Number>   &return_value) const
+  {
+    Assert (return_value.size() == this->n_components,
+            ExcDimensionMismatch (return_value.size(), this->n_components));
+
+    std::copy (function_value_vector.begin(),function_value_vector.end(),
+               return_value.begin());
+  }
+
+
+
+  template <int dim, typename Number>
+  void ConstantFunction<dim, Number>::value_list (
+    const std::vector<Point<dim> > &points,
+    std::vector<Number>            &return_values,
+    const unsigned int              component) const
+  {
+    // To avoid warning of unused parameter
+    (void)points;
+    Assert (component < this->n_components,
+            ExcIndexRange (component, 0, this->n_components));
+    Assert (return_values.size() == points.size(),
+            ExcDimensionMismatch(return_values.size(), points.size()))
+
+    std::fill (return_values.begin(), return_values.end(), function_value_vector[component]);
+  }
+
+
+
+  template <int dim, typename Number>
+  void ConstantFunction<dim, Number>::vector_value_list (
+    const std::vector<Point<dim> > &points,
+    std::vector<Vector<Number> >   &return_values) const
+  {
+    Assert (return_values.size() == points.size(),
+            ExcDimensionMismatch(return_values.size(), points.size()));
+
+    for (unsigned int i=0; i<points.size(); ++i)
+      {
+        Assert (return_values[i].size() == this->n_components,
+                ExcDimensionMismatch(return_values[i].size(), this->n_components));
+        std::copy (function_value_vector.begin(),function_value_vector.end(),
+                   return_values[i].begin());
+      };
+  }
+
+
+
+  template <int dim, typename Number>
+  std::size_t
+  ConstantFunction<dim, Number>::memory_consumption () const
+  {
+    // Here we assume Number is a simple type.
+    return (sizeof(*this) + this->n_components*sizeof(Number));
+  }
 
-
-
-template <int dim, typename Number>
-void ConstantFunction<dim, Number>::vector_value (const Point<dim> &,
-                                                  Vector<Number>   &return_value) const
-{
-  Assert (return_value.size() == this->n_components,
-          ExcDimensionMismatch (return_value.size(), this->n_components));
-
-  std::copy (function_value_vector.begin(),function_value_vector.end(),
-             return_value.begin());
-}
-
-
-
-template <int dim, typename Number>
-void ConstantFunction<dim, Number>::value_list (
-  const std::vector<Point<dim> > &points,
-  std::vector<Number>            &return_values,
-  const unsigned int              component) const
-{
-  // To avoid warning of unused parameter
-  (void)points;
-  Assert (component < this->n_components,
-          ExcIndexRange (component, 0, this->n_components));
-  Assert (return_values.size() == points.size(),
-          ExcDimensionMismatch(return_values.size(), points.size()))
-
-  std::fill (return_values.begin(), return_values.end(), function_value_vector[component]);
-}
-
-
-
-template <int dim, typename Number>
-void ConstantFunction<dim, Number>::vector_value_list (
-  const std::vector<Point<dim> > &points,
-  std::vector<Vector<Number> >   &return_values) const
-{
-  Assert (return_values.size() == points.size(),
-          ExcDimensionMismatch(return_values.size(), points.size()));
-
-  for (unsigned int i=0; i<points.size(); ++i)
-    {
-      Assert (return_values[i].size() == this->n_components,
-              ExcDimensionMismatch(return_values[i].size(), this->n_components));
-      std::copy (function_value_vector.begin(),function_value_vector.end(),
-                 return_values[i].begin());
-    };
-}
-
-
-
-template <int dim, typename Number>
-std::size_t
-ConstantFunction<dim, Number>::memory_consumption () const
-{
-  // Here we assume Number is a simple type.
-  return (sizeof(*this) + this->n_components*sizeof(Number));
 }
 
 //---------------------------------------------------------------------------
index 3d83238438e614e0b6c5d66439897a094eb20bec..41bc637f27f32c6debb77d1c38143d64a0845690 100644 (file)
@@ -17,8 +17,8 @@
 for (S : REAL_SCALARS; dim : SPACE_DIMENSIONS)
 {
     template class Function<dim, S>;
-    template class ZeroFunction<dim, S>;
-    template class ConstantFunction<dim, S>;
+    template class Functions::ZeroFunction<dim, S>;
+    template class Functions::ConstantFunction<dim, S>;
     template class ComponentSelectFunction<dim, S>;
     template class ScalarFunctionFromFunctionObject<dim, S>;
     template class VectorFunctionFromScalarFunctionObject<dim, S>;
@@ -28,8 +28,8 @@ for (S : REAL_SCALARS; dim : SPACE_DIMENSIONS)
 for (S : COMPLEX_SCALARS; dim : SPACE_DIMENSIONS)
 {
     template class Function<dim, S>;
-    template class ZeroFunction<dim, S>;
-    template class ConstantFunction<dim, S>;
+    template class Functions::ZeroFunction<dim, S>;
+    template class Functions::ConstantFunction<dim, S>;
     template class ComponentSelectFunction<dim, S>;
     template class ScalarFunctionFromFunctionObject<dim, S>;
     template class VectorFunctionFromScalarFunctionObject<dim, S>;

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.