]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add VectorFunctionFromScalarFunctionObject.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 22 Nov 2011 15:14:35 +0000 (15:14 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 22 Nov 2011 15:14:35 +0000 (15:14 +0000)
git-svn-id: https://svn.dealii.org/trunk@24758 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

index 5a8a6d76467da7e1d780498d590d027e5bf4b835..2f578603f69eb8319ea812dc0bdfa42789abb690 100644 (file)
@@ -82,6 +82,7 @@ This is in accordance to the other matrix classes.
 <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.
+The VectorFunctionFromScalarFunctionObject class does a similar thing.
 <br>
 (Wolfgang Bangerth, 2011/11/15)
 
index 5073ee55ceb6a992bca518ef6239b3b16d35ef03..e356655847cdb666e0abc40a0b09ea38422ac4a7 100644 (file)
@@ -607,7 +607,7 @@ class ConstantFunction : public ZeroFunction<dim>
  * This is a constant vector-valued function, in which one or more
  * components of the vector have a constant value and all other
  * components are zero.  It is especially useful as a weight function
- * for <tt>VectorTools::integrate_difference</tt>, where it allows to
+ * for VectorTools::integrate_difference, where it allows to
  * integrate only one or a few vector components, rather than the
  * entire vector-valued solution. See the step-20
  * tutorial program for a detailed explanation.
@@ -813,6 +813,84 @@ private:
 };
 
 
+/**
+ * This class is similar to the ScalarFunctionFromFunctionObject class in that it
+ * allows for the easy conversion of a function object to something that
+ * satisfies the interface of the Function base class. The difference is
+ * that here, the given function object is still a scalar function
+ * (i.e. it has a single value at each space point) but that the Function
+ * object generated is vector valued. The number of vector components
+ * is specified in the constructor, where one also selectes a single
+ * one of these vector components that should be filled by the passed
+ * object. The result is a vector Function object that returns zero in
+ * each component except the single selected one where it returns the
+ * value returned by the given as the first argument to the constructor.
+ * 
+ * @note In the above discussion, note the difference between the
+ * (scalar) "function object" (i.e., a C++ object <code>x</code> that can
+ * be called as in <code>x(p)</code>) and the capitalized (vector valued)
+ * "Function object" (i.e., an object of a class that is derived from
+ * the Function base class).
+ * 
+ * @author Wolfgang Bangerth, 2011
+ */
+template <int dim>
+class VectorFunctionFromScalarFunctionObject : 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.
+   * 
+   * @param function_object The scalar function that will form one component
+   *     of the resulting Function object.
+   * @param n_components The total number of vector components of the
+   *     resulting Function object.
+   * @param selected_component The single component that should be
+   *     filled by the first argument.
+   **/
+  VectorFunctionFromScalarFunctionObject (const std_cxx1x::function<double (const Point<dim> &)> &function_object,
+                                          const unsigned int n_components,
+                                          const unsigned int selected_component);
+
+  /**
+   * 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;
+
+  /**
+   * Return all components of a
+   * vector-valued function at a
+   * given point.
+   *
+   * <tt>values</tt> shall have the right
+   * size beforehand,
+   * i.e. #n_components.
+   */
+  virtual void vector_value (const Point<dim>   &p,
+                             Vector<double>     &values) 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;
+  
+  /**
+   * The vector component whose value is to be filled by the
+   * given scalar function.
+   */
+  const unsigned int selected_component;
+};
+
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index f558d533e6c33ef69bedbb8fc6f6bd45a7b62078..47296a638fad839e0c7d0196c498e28fb2efd416 100644 (file)
@@ -66,9 +66,8 @@ double Function<dim>::value (const Point<dim> &,
 
 
 template <int dim>
-void Function<dim>::vector_value (
-  const Point<dim>& p,
-  Vector<double>& v) const
+void Function<dim>::vector_value (const Point<dim>& p,
+                                  Vector<double>& v) const
 {
   AssertDimension(v.size(), this->n_components);
   for (unsigned int i=0;i<this->n_components;++i)
@@ -538,6 +537,50 @@ ScalarFunctionFromFunctionObject<dim>::value (const Point<dim> &p,
 
 
 
+template <int dim>
+VectorFunctionFromScalarFunctionObject<dim>::
+VectorFunctionFromScalarFunctionObject (const std_cxx1x::function<double (const Point<dim> &)> &function_object,
+                                        const unsigned int n_components,
+                                        const unsigned int selected_component)
+:
+Function<dim>(n_components),
+function_object (function_object),
+selected_component (selected_component)
+{}
+
+
+
+template <int dim>
+double
+VectorFunctionFromScalarFunctionObject<dim>::value (const Point<dim> &p,
+                                                    const unsigned int component) const
+{
+  Assert (component < this->n_components,
+          ExcIndexRange (component, 0, this->n_components));
+  
+  if (component == selected_component)
+    return function_object (p);
+  else
+    return 0;
+}
+
+
+
+template <int dim>
+void
+VectorFunctionFromScalarFunctionObject<dim>::
+vector_value (const Point<dim>   &p,
+              Vector<double>     &values) const
+{
+  AssertDimension(values.size(), this->n_components);
+
+  // set everything to zero, and then the right component to its correct value
+  values = 0;
+  values(selected_component) = function_object (p);
+}
+
+
+
 // explicit instantiations
 
 template class Function<1>;
@@ -545,18 +588,21 @@ template class ZeroFunction<1>;
 template class ConstantFunction<1>;
 template class ComponentSelectFunction<1>;
 template class ScalarFunctionFromFunctionObject<1>;
+template class VectorFunctionFromScalarFunctionObject<1>;
 
 template class Function<2>;
 template class ZeroFunction<2>;
 template class ConstantFunction<2>;
 template class ComponentSelectFunction<2>;
 template class ScalarFunctionFromFunctionObject<2>;
+template class VectorFunctionFromScalarFunctionObject<2>;
 
 template class Function<3>;
 template class ZeroFunction<3>;
 template class ConstantFunction<3>;
 template class ComponentSelectFunction<3>;
 template class ScalarFunctionFromFunctionObject<3>;
+template class VectorFunctionFromScalarFunctionObject<3>;
 
 
 DEAL_II_NAMESPACE_CLOSE

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.