#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
};
+/**
+ * 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
}
+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
--- /dev/null
+//-----------------------------------------------------------------------------
+// $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> ();
+}
+
+