--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2010 - 2019 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_function_restriction_h
+#define dealii_function_restriction_h
+
+#include <deal.II/base/function.h>
+#include <deal.II/base/point.h>
+#include <deal.II/base/smartpointer.h>
+#include <deal.II/base/tensor.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Functions
+{
+ /**
+ * This class takes a function in `dim + 1` dimensions and creates a new
+ * function in one dimension lower by restricting one of the coordinates to a
+ * given value. Mathematically this corresponds to taking a function
+ * $f = f(x, y, z)$,
+ * a fixed value, $Z$, and defining a new function (the restriction)
+ * $g = g(x, y) = f(x, y, Z)$.
+ * Using this class, this translates to
+ * @code
+ * Function<3> & function = ...
+ * double z = ...
+ * unsigned int restricted_direction = 2;
+ * CoordinateRestriction<2> restriction(function, restricted_direction, z);
+ * @endcode
+ *
+ * The `dim`-dimensional coordinates on the restriction are ordered starting
+ * from the restricted (`dim + 1`)-coordinate. In particular, this means that
+ * if the $y$-coordinate is locked to $Y$ in 3D, the coordinates are ordered
+ * as $(z, x)$ on the restriction:
+ * $g = g(z, x) = f(x, Y, z)$.
+ * This is the same convention as in BoundingBox::cross_section.
+ */
+ template <int dim>
+ class CoordinateRestriction : public Function<dim>
+ {
+ public:
+ /**
+ * Constructor, takes the (`dim + 1`)-coordinate direction and the value
+ * that the incoming function should be restricted to.
+ *
+ * A pointer to the incoming function is stored internally, so the function
+ * must have a longer lifetime than the created restriction.
+ */
+ CoordinateRestriction(const Function<dim + 1> &function,
+ const unsigned int direction,
+ const double coordinate_value);
+
+ double
+ value(const Point<dim> &point, const unsigned int component) const override;
+
+ Tensor<1, dim>
+ gradient(const Point<dim> & point,
+ const unsigned int component) const override;
+
+ SymmetricTensor<2, dim>
+ hessian(const Point<dim> & point,
+ const unsigned int component) const override;
+
+ private:
+ // The higher-dimensional function that has been restricted.
+ const SmartPointer<const Function<dim + 1>> function;
+
+ // The (`dim + 1`)-coordinate direction that has been restricted.
+ const unsigned int restricted_direction;
+
+ // Value of the restriced coordinate.
+ const double coordinate_value;
+ };
+
+
+
+ /**
+ * This class creates a 1-dimensional function from a `dim + 1` dimensional
+ * function by restricting `dim` of the coordinate values to a given point.
+ * Mathematically this corresponds to taking a function, $f = f(x, y, z)$, and
+ * a point $(Y, Z)$, and defining a new function $g = g(x) = f(x, Y, Z)$.
+ * Using this class, this translates to
+ * @code
+ * Function<3> & function = ...
+ * Point<2> point(y, z);
+ * unsigned int open_direction = 0;
+ * PointRestriction<2> restriction(function, open_direction, point);
+ * @endcode
+ *
+ * The coordinates of the point will be expanded in the higher-dimensional
+ * functions coordinates starting from the open-direction (and wrapping
+ * around). In particular, if we restrict to a point $(Z, X)$ and choose to
+ * keep the y-direction open, the restriction that is created is the function
+ * $g(y) = f(X, y, Z)$.
+ * This is consistent with the convention in BoundingBox::cross_section.
+ */
+ template <int dim>
+ class PointRestriction : public Function<1>
+ {
+ public:
+ /**
+ * Constructor, takes the point that the incoming function should be
+ * restricted to and which (`dim + 1`)-dimensional coordinate direction
+ * should be kept "open".
+ *
+ * A pointer to the incoming function is stored internally, so the function
+ * must have a longer lifetime than the created restriction.
+ */
+ PointRestriction(const Function<dim + 1> &function,
+ const unsigned int open_direction,
+ const Point<dim> & point);
+
+ double
+ value(const Point<1> &point, const unsigned int component) const override;
+
+ Tensor<1, 1>
+ gradient(const Point<1> & point,
+ const unsigned int component) const override;
+
+ SymmetricTensor<2, 1>
+ hessian(const Point<1> &point, const unsigned int component) const override;
+
+ private:
+ // The higher-dimensional function that has been restricted.
+ const SmartPointer<const Function<dim + 1>> function;
+
+ // The (`dim + 1`)-coordinate direction that is kept "open"
+ const unsigned int open_direction;
+
+ // The point that we have restricted the above function to.
+ const Point<dim> point;
+ };
+
+} // namespace Functions
+
+
+namespace internal
+{
+ /**
+ * Creates a (`dim + 1`)-dimensional point by copying over the coordinates of
+ * the incoming `dim`-dimensional point and setting the "missing"
+ * (`dim + 1`)-dimensional component to the incoming coordinate value.
+ *
+ * For example, given the input
+ * $\{(x, y), 2, z \}$ this function creates the point $(x, y, z)$.
+ *
+ * The coordinates of the `dim`-dimensional point are written to the
+ * coordinates of the (`dim + 1`)-dimensional point in the order of the
+ * convention given by the function coordinate_to_one_dim_higher. Thus, the
+ * order of coordinates on the lower-dimensional point are not preserved:
+ * $\{(z, x), 1, y \}$ creates the point $(x, y, z)$.
+ */
+ template <int dim>
+ Point<dim + 1>
+ create_higher_dim_point(const Point<dim> & point,
+ const unsigned int component_in_dim_plus_1,
+ const double coordinate_value);
+} // namespace internal
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif /* dealii_function_restriction_h */
function_lib.cc
function_lib_cutoff.cc
function_parser.cc
+ function_restriction.cc
function_spherical.cc
function_time.cc
function_tools.cc
data_out_base.inst.in
function.inst.in
function_level_set.inst.in
+ function_restriction.inst.in
function_time.inst.in
function_tools.inst.in
incremental_function.inst.in
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2010 - 2019 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/bounding_box.h>
+#include <deal.II/base/function_restriction.h>
+
+DEAL_II_NAMESPACE_OPEN
+namespace internal
+{
+ template <int dim>
+ Point<dim + 1>
+ create_higher_dim_point(const Point<dim> & point,
+ const unsigned int component_in_dim_plus_1,
+ const double coordinate_value)
+ {
+ AssertIndexRange(component_in_dim_plus_1, dim + 1);
+
+ Point<dim + 1> output;
+ output(component_in_dim_plus_1) = coordinate_value;
+ for (unsigned int d = 0; d < dim; ++d)
+ {
+ const unsigned int component_to_write_to =
+ dealii::internal::coordinate_to_one_dim_higher<dim>(
+ component_in_dim_plus_1, d);
+ output(component_to_write_to) = point(d);
+ }
+
+ return output;
+ }
+} // namespace internal
+
+
+
+namespace Functions
+{
+ template <int dim>
+ CoordinateRestriction<dim>::CoordinateRestriction(
+ const Function<dim + 1> &function,
+ const unsigned int direction,
+ const double coordinate_value)
+ : function(&function)
+ , restricted_direction(direction)
+ , coordinate_value(coordinate_value)
+ {
+ AssertIndexRange(restricted_direction, dim + 1);
+ }
+
+
+
+ template <int dim>
+ double
+ CoordinateRestriction<dim>::value(const Point<dim> & point,
+ const unsigned int component) const
+ {
+ const Point<dim + 1> full_point =
+ internal::create_higher_dim_point(point,
+ restricted_direction,
+ coordinate_value);
+
+ return function->value(full_point, component);
+ }
+
+
+
+ template <int dim>
+ Tensor<1, dim>
+ CoordinateRestriction<dim>::gradient(const Point<dim> & point,
+ const unsigned int component) const
+ {
+ const Point<dim + 1> full_point =
+ internal::create_higher_dim_point(point,
+ restricted_direction,
+ coordinate_value);
+
+ const Tensor<1, dim + 1> full_gradient =
+ function->gradient(full_point, component);
+
+ // CoordinateRestriction is constant in restricted direction. Through away
+ // the derivatives with respect to this direction and copy the other
+ // values.
+ Tensor<1, dim> grad;
+ for (unsigned int d = 0; d < dim; ++d)
+ {
+ const unsigned int index_to_write_from =
+ internal::coordinate_to_one_dim_higher<dim>(restricted_direction, d);
+ grad[d] = full_gradient[index_to_write_from];
+ }
+ return grad;
+ }
+
+
+
+ template <int dim>
+ SymmetricTensor<2, dim>
+ CoordinateRestriction<dim>::hessian(const Point<dim> & point,
+ const unsigned int component) const
+ {
+ const Point<dim + 1> full_point =
+ internal::create_higher_dim_point(point,
+ restricted_direction,
+ coordinate_value);
+
+ const Tensor<2, dim + 1> full_hessian =
+ function->hessian(full_point, component);
+
+ // CoordinateRestriction is constant in restricted direction. Through away
+ // the derivatives with respect to this direction and copy the other
+ // values.
+ SymmetricTensor<2, dim> hess;
+ for (unsigned int i = 0; i < dim; ++i)
+ {
+ const unsigned int i_to_write_from =
+ internal::coordinate_to_one_dim_higher<dim>(restricted_direction, i);
+ for (unsigned int j = 0; j < dim; ++j)
+ {
+ const unsigned int j_to_write_from =
+ internal::coordinate_to_one_dim_higher<dim>(restricted_direction,
+ j);
+ hess[i][j] = full_hessian[i_to_write_from][j_to_write_from];
+ }
+ }
+ return hess;
+ }
+
+
+
+ template <int dim>
+ PointRestriction<dim>::PointRestriction(const Function<dim + 1> &function,
+ const unsigned int open_direction,
+ const Point<dim> & point)
+ : function(&function)
+ , open_direction(open_direction)
+ , point(point)
+ {
+ AssertIndexRange(open_direction, dim + 1);
+ }
+
+
+
+ template <int dim>
+ double
+ PointRestriction<dim>::value(const Point<1> & point_1D,
+ const unsigned int component) const
+ {
+ const Point<dim + 1> full_point =
+ internal::create_higher_dim_point(point, open_direction, point_1D(0));
+ return function->value(full_point, component);
+ }
+
+
+
+ template <int dim>
+ Tensor<1, 1>
+ PointRestriction<dim>::gradient(const Point<1> & point_1D,
+ const unsigned int component) const
+ {
+ const Point<dim + 1> full_point =
+ internal::create_higher_dim_point(point, open_direction, point_1D(0));
+ const Tensor<1, dim + 1> full_gradient =
+ function->gradient(full_point, component);
+
+ // The PointRestrictions is constant in all but the open direction. Throw
+ // away the derivatives in all but this direction.
+ Tensor<1, 1> grad;
+ grad[0] = full_gradient[open_direction];
+ return grad;
+ }
+
+
+
+ template <int dim>
+ SymmetricTensor<2, 1>
+ PointRestriction<dim>::hessian(const Point<1> & point_1D,
+ const unsigned int component) const
+ {
+ const Point<dim + 1> full_point =
+ internal::create_higher_dim_point(point, open_direction, point_1D(0));
+ const Tensor<2, dim + 1> full_hessian =
+ function->hessian(full_point, component);
+
+ // The PointRestrictions is constant in all but the open direction. Throw
+ // away the derivatives in all but this direction.
+ SymmetricTensor<2, 1> hess;
+ hess[0][0] = full_hessian[open_direction][open_direction];
+ return hess;
+ }
+
+
+
+} // namespace Functions
+#include "function_restriction.inst"
+DEAL_II_NAMESPACE_CLOSE
--- /dev/null
+for (deal_II_dimension : DIMENSIONS)
+ {
+ template class Functions::PointRestriction<deal_II_dimension>;
+
+ template Point<deal_II_dimension> internal::create_higher_dim_point(
+ const Point<deal_II_dimension - 1> &, const unsigned int, const double);
+
+#if deal_II_dimension != 3
+ template class Functions::CoordinateRestriction<deal_II_dimension>;
+#endif
+ }
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2010 - 2019 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+/*
+ * Test that CoordinateRestriction returns the derivatives of its underlying
+ * function truncated correctly.
+ */
+
+
+#include <deal.II/base/function_restriction.h>
+
+#include "../tests.h"
+
+namespace
+{
+ using namespace dealii;
+
+ /*
+ * A simple test-function which is constant in space, but varies
+ * in each component.
+ */
+ template <int dim>
+ class TestFunction : public Function<dim>
+ {
+ public:
+ double
+ value(const Point<dim> &, const unsigned int) const override
+ {
+ return 1;
+ }
+
+ Tensor<1, dim>
+ gradient(const Point<dim> &, const unsigned int) const override
+ {
+ Tensor<1, dim> tensor;
+ for (unsigned int i = 0; i < dim; ++i)
+ tensor[i] = i + 1;
+ return tensor;
+ }
+
+ SymmetricTensor<2, dim>
+ hessian(const Point<dim> &, const unsigned int) const override
+ {
+ SymmetricTensor<2, dim> tensor;
+ for (unsigned int i = 0; i < dim; ++i)
+ for (unsigned int j = 0; j < dim; ++j)
+ tensor[i][j] = 1 + i + j;
+
+ return tensor;
+ }
+ };
+
+
+ /*
+ * Create a CoordinateRestriction that has the incoming direction of the
+ * above TestFunction restricted. Write value, gradient & hessian
+ * to deallog to check that they are truncated as expected.
+ */
+ template <int dim>
+ void
+ test_truncates_derivatives_correctly(const unsigned int restricted_direction)
+ {
+ deallog << "Restricted direction = " << restricted_direction << std::endl;
+
+ const TestFunction<dim> function;
+ const double coordinate_value = 0;
+
+ const Functions::CoordinateRestriction<dim - 1> restriction(
+ function, restricted_direction, coordinate_value);
+
+ const Point<dim - 1> point;
+ const unsigned int component = 0;
+
+ deallog << "value = " << restriction.value(point, component) << std::endl;
+ deallog << "gradient = " << restriction.gradient(point, component)
+ << std::endl;
+ deallog << "Hessian = " << restriction.hessian(point, component)
+ << std::endl;
+
+ deallog << std::endl;
+ }
+
+
+
+ template <int dim>
+ void
+ test_all_coordinate_directions()
+ {
+ deallog << "dim = " << dim << std::endl;
+ for (unsigned int direction = 0; direction < dim; direction++)
+ test_truncates_derivatives_correctly<dim>(direction);
+ }
+
+} // namespace
+
+int
+main()
+{
+ initlog();
+ test_all_coordinate_directions<2>();
+ test_all_coordinate_directions<3>();
+}
--- /dev/null
+
+DEAL::dim = 2
+DEAL::Restricted direction = 0
+DEAL::value = 1.00000
+DEAL::gradient = 2.00000
+DEAL::Hessian = 3.00000
+DEAL::
+DEAL::Restricted direction = 1
+DEAL::value = 1.00000
+DEAL::gradient = 1.00000
+DEAL::Hessian = 1.00000
+DEAL::
+DEAL::dim = 3
+DEAL::Restricted direction = 0
+DEAL::value = 1.00000
+DEAL::gradient = 2.00000 3.00000
+DEAL::Hessian = 3.00000 4.00000 4.00000 5.00000
+DEAL::
+DEAL::Restricted direction = 1
+DEAL::value = 1.00000
+DEAL::gradient = 3.00000 1.00000
+DEAL::Hessian = 5.00000 3.00000 3.00000 1.00000
+DEAL::
+DEAL::Restricted direction = 2
+DEAL::value = 1.00000
+DEAL::gradient = 1.00000 2.00000
+DEAL::Hessian = 1.00000 2.00000 2.00000 3.00000
+DEAL::
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2010 - 2019 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/base/function_restriction.h>
+
+#include "../tests.h"
+
+namespace
+{
+ using namespace dealii;
+
+ // Test the function create_higher_dim_point by making sure that we
+ // can construct the point (x, y, z),
+ const double x = 1, y = 2, z = 3;
+ // for each possible component. That is, we test that we get the following
+ // results:
+ // {point, component, coordinate} -> result
+ // {(y, z), 0, x} -> (x, y, z)
+ // {(z, x), 1, y} -> (x, y, z)
+ // {(y, z), 2, x} -> (x, y, z)
+ //
+ // The same thing is done for the 2D case.
+ void
+ test3D()
+ {
+ const int dim = 3;
+ deallog << "dim=" << dim << std::endl;
+ {
+ const unsigned int component = 0;
+ const Point<dim - 1> low_dim_point(y, z);
+ const Point<dim> point =
+ internal::create_higher_dim_point(low_dim_point, component, x);
+ deallog << "Point: " << point << std::endl;
+ }
+ {
+ const unsigned int component = 1;
+ // The order (z,x) here follows convention in
+ // coordinate_to_one_dim_higher.
+ const Point<dim - 1> low_dim_point(z, x);
+ const Point<dim> point =
+ internal::create_higher_dim_point(low_dim_point, component, y);
+ deallog << "Point: " << point << std::endl;
+ }
+ {
+ const unsigned int component = 2;
+ Point<dim - 1> low_dim_point(x, y);
+ Point<dim> point =
+ internal::create_higher_dim_point(low_dim_point, component, z);
+ deallog << "Point: " << point << std::endl;
+ }
+ }
+
+
+
+ void
+ test2D()
+ {
+ const int dim = 2;
+ deallog << "dim=" << dim << std::endl;
+ {
+ const unsigned int component = 0;
+ const Point<dim - 1> low_dim_point(y);
+ const Point<dim> point =
+ internal::create_higher_dim_point(low_dim_point, component, x);
+ deallog << "Point: " << point << std::endl;
+ }
+ {
+ const unsigned int component = 1;
+ const Point<dim - 1> low_dim_point(x);
+ const Point<dim> point =
+ internal::create_higher_dim_point(low_dim_point, component, y);
+ deallog << "Point: " << point << std::endl;
+ }
+ }
+
+
+
+ void
+ test1D()
+ {
+ const int dim = 1;
+ deallog << "dim=" << dim << std::endl;
+ const unsigned int component = 0;
+ const Point<dim - 1> low_dim_point;
+ const Point<dim> point =
+ internal::create_higher_dim_point(low_dim_point, component, x);
+ deallog << "Point: " << point << std::endl;
+ }
+} // namespace
+
+
+
+int
+main()
+{
+ initlog();
+ test1D();
+ test2D();
+ test3D();
+}
--- /dev/null
+
+DEAL::dim=1
+DEAL::Point: 1.00000
+DEAL::dim=2
+DEAL::Point: 1.00000 2.00000
+DEAL::Point: 1.00000 2.00000
+DEAL::dim=3
+DEAL::Point: 1.00000 2.00000 3.00000
+DEAL::Point: 1.00000 2.00000 3.00000
+DEAL::Point: 1.00000 2.00000 3.00000
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2010 - 2019 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+/*
+ * Test that PointRestriction returns the derivatives of its underlying function
+ * truncated correctly.
+ */
+
+
+#include <deal.II/base/function_restriction.h>
+
+#include "../tests.h"
+
+namespace
+{
+ using namespace dealii;
+
+ /*
+ * A simple test-function which is constant in space, but varies
+ * in each gradient component, and on the diagonal of the Hessian.
+ */
+ template <int dim>
+ class TestFunction : public Function<dim>
+ {
+ public:
+ double
+ value(const Point<dim> &, const unsigned int) const override
+ {
+ return 1;
+ }
+
+ Tensor<1, dim>
+ gradient(const Point<dim> &, const unsigned int) const override
+ {
+ Tensor<1, dim> tensor;
+ for (unsigned int i = 0; i < dim; ++i)
+ tensor[i] = i + 1;
+ return tensor;
+ }
+
+ SymmetricTensor<2, dim>
+ hessian(const Point<dim> &, const unsigned int) const override
+ {
+ SymmetricTensor<2, dim> tensor;
+ for (unsigned int i = 0; i < dim; ++i)
+ tensor[i][i] = i + 1;
+
+ return tensor;
+ }
+ };
+
+
+ /*
+ * Create a PointRestriction that has the incoming direction open and has
+ * the above TestFunction as underlying function. Check that the values that
+ *
+ * PointRestriction::values
+ * PointRestriction::gradient
+ * PointRestriction::hessian
+ *
+ * return are
+ *
+ * value
+ * gradient[open_direction]
+ * hessian[open_direction][open_direction]
+ *
+ * where value, gradient and hessian are the values from TestFunction.
+ */
+ template <int dim>
+ void
+ test_truncates_derivatives_correctly(const unsigned int open_direction)
+ {
+ deallog << "open direction = " << open_direction << std::endl;
+ const Point<dim - 1> locked_point;
+ const TestFunction<dim> function;
+
+ const Functions::PointRestriction<dim - 1> restriction(function,
+ open_direction,
+ locked_point);
+
+ const Point<1> point;
+ const unsigned int component = 0;
+
+ deallog << restriction.value(point, component) << std::endl;
+ deallog << restriction.gradient(point, component) << std::endl;
+ deallog << restriction.hessian(point, component) << std::endl;
+
+ deallog << std::endl;
+ }
+
+
+
+ template <int dim>
+ void
+ test_all_coordinate_directions()
+ {
+ deallog << "dim = " << dim << std::endl;
+
+ for (unsigned int direction = 0; direction < dim; direction++)
+ test_truncates_derivatives_correctly<dim>(direction);
+ }
+
+} // namespace
+
+int
+main()
+{
+ initlog();
+ test_all_coordinate_directions<2>();
+ test_all_coordinate_directions<3>();
+}
--- /dev/null
+
+DEAL::dim = 2
+DEAL::open direction = 0
+DEAL::1.00000
+DEAL::1.00000
+DEAL::1.00000
+DEAL::
+DEAL::open direction = 1
+DEAL::1.00000
+DEAL::2.00000
+DEAL::2.00000
+DEAL::
+DEAL::dim = 3
+DEAL::open direction = 0
+DEAL::1.00000
+DEAL::1.00000
+DEAL::1.00000
+DEAL::
+DEAL::open direction = 1
+DEAL::1.00000
+DEAL::2.00000
+DEAL::2.00000
+DEAL::
+DEAL::open direction = 2
+DEAL::1.00000
+DEAL::3.00000
+DEAL::3.00000
+DEAL::