]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Provide ScalarFunctionFromFunctionObject.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 16 Nov 2011 14:33:54 +0000 (14:33 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 16 Nov 2011 14:33:54 +0000 (14:33 +0000)
git-svn-id: https://svn.dealii.org/trunk@24750 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/base/function.h
deal.II/source/base/function.cc
tests/base/functions_05.cc [new file with mode: 0644]
tests/base/functions_05/cmp/generic [new file with mode: 0644]

index 8747872fc0617afc7377dc287772e9c13ab80333..09a6e01a97300d5d8b87d05a8946f33773068948 100644 (file)
@@ -78,6 +78,12 @@ enabled due to a missing include file in file
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: The class ScalarFunctionFromFunctionObject provides a quick way to
+convert an arbitrary function into a function object that can be passed
+to part of the library that require the Function@<dim@> interface.
+<br>
+(Wolfgang Bangerth, 2011/10/30)
+
 <li> New: FETools::get_fe_from_name() can now return objects of type FE_Nothing.
 <br>
 (Scott Miller, Jonathan Pitt, 2011/11/10)
index c9daa892be8e5c2c57854ef389f25c79366d9d82..5073ee55ceb6a992bca518ef6239b3b16d35ef03 100644 (file)
@@ -19,6 +19,8 @@
 #include <deal.II/base/subscriptor.h>
 #include <deal.II/base/tensor.h>
 #include <deal.II/base/point.h>
+#include <deal.II/base/std_cxx1x/function.h>
+
 #include <vector>
 
 DEAL_II_NAMESPACE_OPEN
@@ -705,6 +707,112 @@ class ComponentSelectFunction : public ConstantFunction<dim>
 };
 
 
+/**
+ * This class provides a way to convert a scalar function of the kind
+ * @code
+ *   double foo (const Point<dim> &);
+ * @endcode
+ * into an object of type Function@<dim@>. 
+ * Since the argument returns a scalar, the result is clearly a
+ * Function object for which <code>function.n_components==1</code>.
+ * The class works by storing a pointer to the given function and
+ * every time <code>function.value(p,component)</code> is called,
+ * calls <code>foo(p)</code> and returns the corresponding value. It
+ * also makes sure that <code>component</code> is in fact zero, as needs
+ * be for scalar functions.
+ * 
+ * The class provides an easy way to turn a simple global function into
+ * something that has the required Function@<dim@> interface for operations
+ * like VectorTools::interpolate_boundary_values() etc., and thereby
+ * allows for simpler experimenting without having to write all the
+ * boiler plate code of declaring a class that is derived from Function
+ * and implementing the Function::value() function.
+ * 
+ * The class gains additional expressvive power because the argument it
+ * takes does not have to be a pointer to an actual function. Rather, it is
+ * a function object, i.e., it can also be the result of call to 
+ * std::bind (or boost::bind) or some other object that can be called with
+ * a single argument. For example, if you need a Function object that
+ * returns the norm of a point, you could write it like so:
+ * @code
+ *   template <int dim> 
+ *   class Norm : public Function<dim> {
+ *     public:
+ *       virtual double value (const Point<dim> &p,
+ *                             const unsigned int component) const {
+ *         Assert (component == 0, ExcMessage ("This object is scalar!"));
+ *         return p.norm();
+ *       }
+ *    };
+ * 
+ *    Norm<2> my_norm_object;
+ * @endcode
+ * and then pass the <code>my_norm_object</code> around, or you could write it
+ * like so:
+ * @code
+ *   ScalarFunctionFromFunctionObject<dim> my_norm_object (&Point<dim>::norm);
+ * @endcode
+ * 
+ * Similarly, to generate an object that computes the distance to a point
+ * <code>q</code>, we could do this:
+ * @code
+ *   template <int dim> 
+ *   class DistanceTo : public Function<dim> {
+ *     public:
+ *       DistanceTo (const Point<dim> &q) : q(q) {}
+ *       virtual double value (const Point<dim> &p,
+ *                             const unsigned int component) const {
+ *         Assert (component == 0, ExcMessage ("This object is scalar!"));
+ *         return q.distance(p);
+ *       }
+ *     private:
+ *       const Point<dim> q;
+ *    };
+ * 
+ *    Point<2> q (2,3);
+ *    DistanceTo<2> my_distance_object;
+ * @endcode
+ * or we could write it like so:
+ * @code
+ *    ScalarFunctionFromFunctionObject<dim>
+ *      my_distance_object (std_cxx1x::bind (&Point<dim>::distance,
+ *                                          q,
+ *                                          std_cxx1x::_1));
+ * @endcode
+ * The savings in work to write this are apparent.
+ * 
+ * @author Wolfgang Bangerth, 2011
+ */
+template <int dim>
+class ScalarFunctionFromFunctionObject : public Function<dim>
+{
+public:
+  /**
+   * Given a function object that takes a Point and returns a double
+   * value, convert this into an object that matches the Function<dim>
+   * interface.
+   */
+  ScalarFunctionFromFunctionObject (const std_cxx1x::function<double (const Point<dim> &)> &function_object);
+
+  /**
+   * Return the value of the
+   * function at the given
+   * point. Returns the value the
+   * function given to the constructor
+   * produces for this point.
+   */
+  virtual double value (const Point<dim>   &p,
+                       const unsigned int  component = 0) const;
+
+private:
+  /**
+   * The function object which we call when this class's value() or
+   * value_list() functions are called.
+   **/
+  const std_cxx1x::function<double (const Point<dim> &)> function_object;
+};
+
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index 2ac740564749766d245c5ee109cf350ef0d793ad..f558d533e6c33ef69bedbb8fc6f6bd45a7b62078 100644 (file)
@@ -516,19 +516,47 @@ ComponentSelectFunction<dim>::memory_consumption () const
 }
 
 
+template <int dim>
+ScalarFunctionFromFunctionObject<dim>::
+ScalarFunctionFromFunctionObject (const std_cxx1x::function<double (const Point<dim> &)> &function_object)
+:
+Function<dim>(1),
+function_object (function_object)
+{}
+
+
+
+template <int dim>
+double
+ScalarFunctionFromFunctionObject<dim>::value (const Point<dim> &p,
+                                             const unsigned int component) const
+{
+  Assert (component == 0,
+         ExcMessage ("This object represents only scalar functions"));
+  return function_object (p);
+}
+
+
+
 // explicit instantiations
 
 template class Function<1>;
 template class ZeroFunction<1>;
 template class ConstantFunction<1>;
 template class ComponentSelectFunction<1>;
+template class ScalarFunctionFromFunctionObject<1>;
+
 template class Function<2>;
 template class ZeroFunction<2>;
 template class ConstantFunction<2>;
 template class ComponentSelectFunction<2>;
+template class ScalarFunctionFromFunctionObject<2>;
+
 template class Function<3>;
 template class ZeroFunction<3>;
 template class ConstantFunction<3>;
 template class ComponentSelectFunction<3>;
+template class ScalarFunctionFromFunctionObject<3>;
+
 
 DEAL_II_NAMESPACE_CLOSE
diff --git a/tests/base/functions_05.cc b/tests/base/functions_05.cc
new file mode 100644 (file)
index 0000000..43a8673
--- /dev/null
@@ -0,0 +1,83 @@
+//-----------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2007, 2008, 2010, 2011 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
+//    further information on this license.
+//
+//-----------------------------------------------------------------------------
+
+// Test ScalarFunctionFromFunctionObject
+
+#include "../tests.h"
+#include <deal.II/base/function.h>
+
+
+template <int dim>
+void check1 ()
+{
+  ScalarFunctionFromFunctionObject<dim>
+    object (&Point<dim>::norm);
+
+  for (unsigned int i=0; i<10; ++i)
+    {
+      Point<dim> p;
+      for (unsigned int d=0; d<dim; ++d)
+       p[d] = i+d;
+
+      Assert (object.value(p) == p.norm(),
+             ExcInternalError());
+    }
+
+  deallog << "OK" << std::endl;
+}
+
+
+template <int dim>
+void check2 ()
+{
+  Point<dim> q;
+  for (unsigned int d=0; d<dim; ++d)
+    q[d] = d;
+
+  ScalarFunctionFromFunctionObject<dim>
+    object (std_cxx1x::bind (&Point<dim>::distance,
+                            q,
+                            std_cxx1x::_1));
+
+  for (unsigned int i=0; i<10; ++i)
+    {
+      Point<dim> p;
+      for (unsigned int d=0; d<dim; ++d)
+       p[d] = i+d;
+
+      Assert (object.value(p) == q.distance (p),
+             ExcInternalError());
+    }
+
+  deallog << "OK" << std::endl;
+}
+
+
+int main()
+{
+  std::string logname = JobIdentifier::base_name(__FILE__) + std::string("/output");
+  std::ofstream logfile(logname.c_str());
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  check1<1> ();
+  check1<2> ();
+  check1<3> ();
+
+  check2<1> ();
+  check2<2> ();
+  check2<3> ();
+}
+
+
diff --git a/tests/base/functions_05/cmp/generic b/tests/base/functions_05/cmp/generic
new file mode 100644 (file)
index 0000000..0aa61ff
--- /dev/null
@@ -0,0 +1,7 @@
+
+DEAL::OK
+DEAL::OK
+DEAL::OK
+DEAL::OK
+DEAL::OK
+DEAL::OK

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.