]> https://gitweb.dealii.org/ - dealii.git/commitdiff
change inheritance structure of ZeroFunction and ConstantFunction
authorDenis Davydov <davydden@gmail.com>
Sat, 19 Aug 2017 14:22:32 +0000 (16:22 +0200)
committerDenis Davydov <davydden@gmail.com>
Sat, 19 Aug 2017 14:22:32 +0000 (16:22 +0200)
include/deal.II/base/function.h
include/deal.II/base/function.templates.h

index 8602e357a7c49983699b2baba4179e1ab090ce88..5b6fab9e41499be2a25d527d59f89254ebf6cc4c 100644 (file)
@@ -350,82 +350,15 @@ public:
 namespace Functions
 {
 
-
-  /**
-   * Provide a function which always returns zero. Obviously, also the derivatives
-   * 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>
-  {
-  public:
-    /**
-     * Constructor. The number of components is preset to one.
-     */
-    ZeroFunction (const unsigned int n_components = 1);
-
-    /**
-     * Virtual destructor; absolutely necessary in this case.
-     *
-     */
-    virtual ~ZeroFunction ();
-
-    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>            &values,
-                             const unsigned int              component = 0) const;
-
-    virtual void vector_value_list (const std::vector<Point<dim> > &points,
-                                    std::vector<Vector<Number> >   &values) const;
-
-    virtual Tensor<1,dim, Number> gradient (const Point<dim> &p,
-                                            const unsigned int component = 0) const;
-
-    virtual void vector_gradient (const Point<dim>            &p,
-                                  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 derivatives 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>
+  class ConstantFunction : public Function<dim, Number>
   {
   public:
     /**
@@ -471,6 +404,19 @@ namespace Functions
     virtual void vector_value_list (const std::vector<Point<dim> > &points,
                                     std::vector<Vector<Number> >   &return_values) const;
 
+    virtual Tensor<1,dim, Number> gradient (const Point<dim> &p,
+                                            const unsigned int component = 0) const;
+
+    virtual void vector_gradient (const Point<dim>            &p,
+                                  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;
+
     std::size_t memory_consumption () const;
 
   protected:
@@ -480,6 +426,37 @@ namespace Functions
     std::vector<Number> function_value_vector;
   };
 
+
+
+  /**
+   * Provide a function which always returns zero. Obviously, also the derivatives
+   * 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 ConstantFunction<dim, Number>
+  {
+  public:
+    /**
+     * Constructor. The number of components is preset to one.
+     */
+    ZeroFunction (const unsigned int n_components = 1);
+
+    /**
+     * Destructor.
+     *
+     */
+    virtual ~ZeroFunction ();
+
+  };
+
 }
 
 /**
index 14e87c92ca2963282c27336a78130d7c25e32823..cfac4bd9212fa94394fac6f7a9e660709214ce34 100644 (file)
@@ -303,7 +303,7 @@ namespace Functions
   template <int dim, typename Number>
   ZeroFunction<dim, Number>::ZeroFunction (const unsigned int n_components)
     :
-    Function<dim, Number> (n_components)
+    ConstantFunction<dim, Number> (Number(), n_components)
   {}
 
 
@@ -311,108 +311,6 @@ namespace Functions
   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 ();
-      };
-  }
-
 }
 
 //---------------------------------------------------------------------------
@@ -424,7 +322,7 @@ namespace Functions
   ConstantFunction<dim, Number>::ConstantFunction (const Number value,
                                                    const unsigned int n_components)
     :
-    ZeroFunction<dim, Number> (n_components),
+    Function<dim, Number> (n_components),
     function_value_vector (n_components, value)
   {}
 
@@ -432,7 +330,7 @@ namespace Functions
   ConstantFunction<dim, Number>::
   ConstantFunction (const std::vector<Number> &values)
     :
-    ZeroFunction<dim, Number> (values.size()),
+    Function<dim, Number> (values.size()),
     function_value_vector (values)
   {}
 
@@ -441,7 +339,7 @@ namespace Functions
   ConstantFunction<dim, Number>::
   ConstantFunction (const Vector<Number> &values)
     :
-    ZeroFunction<dim, Number> (values.size()),
+    Function<dim, Number> (values.size()),
     function_value_vector (values.size())
   {
     Assert (values.size() == function_value_vector.size(),
@@ -454,7 +352,7 @@ namespace Functions
   ConstantFunction<dim, Number>::
   ConstantFunction (const Number *begin_ptr, const unsigned int n_components)
     :
-    ZeroFunction<dim, Number> (n_components),
+    Function<dim, Number> (n_components),
     function_value_vector (n_components)
   {
     Assert (begin_ptr != nullptr, ExcMessage ("Null pointer encountered!"));
@@ -539,6 +437,60 @@ namespace Functions
     return (sizeof(*this) + this->n_components*sizeof(Number));
   }
 
+
+
+  template <int dim, typename Number>
+  Tensor<1,dim,Number> ConstantFunction<dim, Number>::gradient (const Point<dim> &,
+      const unsigned int) const
+  {
+    return Tensor<1,dim,Number>();
+  }
+
+
+  template <int dim, typename Number>
+  void ConstantFunction<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 ConstantFunction<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 ConstantFunction<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 ();
+      };
+  }
+
+
 }
 
 //---------------------------------------------------------------------------

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.