From 8ad86a52b716835c06a2716df3299a421d8d3dc4 Mon Sep 17 00:00:00 2001 From: Matthias Maier Date: Tue, 25 Aug 2015 13:40:01 -0500 Subject: [PATCH] LinearOperator: Also test for complex number support --- tests/lac/linear_operator_01.cc | 16 +- tests/lac/linear_operator_01a.cc | 242 ++++++++++++++++++ .../linear_operator_01a.with_cxx11=on.output | 20 ++ tests/lac/linear_operator_02a.cc | 159 ++++++++++++ .../linear_operator_02a.with_cxx11=on.output | 27 ++ 5 files changed, 456 insertions(+), 8 deletions(-) create mode 100644 tests/lac/linear_operator_01a.cc create mode 100644 tests/lac/linear_operator_01a.with_cxx11=on.output create mode 100644 tests/lac/linear_operator_02a.cc create mode 100644 tests/lac/linear_operator_02a.with_cxx11=on.output diff --git a/tests/lac/linear_operator_01.cc b/tests/lac/linear_operator_01.cc index 0f47a0d139..0f4f99385c 100644 --- a/tests/lac/linear_operator_01.cc +++ b/tests/lac/linear_operator_01.cc @@ -29,19 +29,19 @@ using namespace dealii; struct LeftVector { typedef double value_type; - double value; + value_type value; - LeftVector & operator = (double new_value) + LeftVector & operator = (value_type new_value) { value = new_value; return *this; } - LeftVector & operator *= (double scale) + LeftVector & operator *= (value_type scale) { value *= scale; return *this; } - LeftVector & operator /= (double scale) + LeftVector & operator /= (value_type scale) { value /= scale; return *this; @@ -64,19 +64,19 @@ struct LeftVector struct RightVector { typedef double value_type; - double value; + value_type value; - RightVector & operator = (double new_value) + RightVector & operator = (value_type new_value) { value = new_value; return *this; } - RightVector & operator *= (double scale) + RightVector & operator *= (value_type scale) { value *= scale; return *this; } - RightVector & operator /= (double scale) + RightVector & operator /= (value_type scale) { value /= scale; return *this; diff --git a/tests/lac/linear_operator_01a.cc b/tests/lac/linear_operator_01a.cc new file mode 100644 index 0000000000..e1511c96e3 --- /dev/null +++ b/tests/lac/linear_operator_01a.cc @@ -0,0 +1,242 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +// Variant of linear_operator_01 that uses std::complex as +// value_type: Test the LinearOperator template on a trivial vector +// implementation :: RightVector -> LeftVector with complex numbers. + +#include "../tests.h" + +#include + +#include +#include + + +using namespace dealii; + +// Dummy vectors with different, non convertible types: + +struct LeftVector +{ + typedef std::complex value_type; + value_type value; + + LeftVector & operator = (value_type new_value) + { + value = new_value; + return *this; + } + LeftVector & operator *= (value_type scale) + { + value *= scale; + return *this; + } + LeftVector & operator /= (value_type scale) + { + value /= scale; + return *this; + } + LeftVector & operator += (const LeftVector &u) + { + value += u.value; + return *this; + } + int size() const + { + return 1; + } + std::size_t memory_consumption () const + { + return 1; + } +}; + +struct RightVector +{ + typedef std::complex value_type; + value_type value; + + RightVector & operator = (value_type new_value) + { + value = new_value; + return *this; + } + RightVector & operator *= (value_type scale) + { + value *= scale; + return *this; + } + RightVector & operator /= (value_type scale) + { + value /= scale; + return *this; + } + RightVector & operator += (const RightVector &u) + { + value += u.value; + return *this; + } + int size() const + { + return 1; + } + std::size_t memory_consumption () const + { + return 1; + } +}; + +int main() +{ + initlog(); + + // Create to distinct linear operators: + + LinearOperator multiply2; + multiply2.vmult = [](LeftVector &v, const RightVector &u) + { + v.value = 2. * u.value; + }; + multiply2.vmult_add = [](LeftVector &v, const RightVector &u) + { + v.value += 2. * u.value; + }; + multiply2.Tvmult = [](RightVector &v, const LeftVector &u) + { + v.value = 2. * u.value; + }; + multiply2.Tvmult_add = [](RightVector &v, const LeftVector &u) + { + v.value += 2. * u.value; + }; + multiply2.reinit_range_vector = [](LeftVector &, bool) + { + }; + multiply2.reinit_domain_vector = [](RightVector &, bool) + { + }; + + auto multiply4 = multiply2; + multiply4.vmult = [](LeftVector &v, const RightVector &u) + { + v.value = 4. * u.value; + }; + multiply4.vmult_add = [](LeftVector &v, const RightVector &u) + { + v.value += 4. * u.value; + }; + multiply4.Tvmult = [](RightVector &v, const LeftVector &u) + { + v.value = 4. * u.value; + }; + multiply4.Tvmult_add = [](RightVector &v, const LeftVector &u) + { + v.value += 4. * u.value; + }; + + + // Small unit tests for all functions: + + RightVector u = { 4. }; + LeftVector v = { 0. }; + + // vmult, vmult_add + + multiply2.vmult(v, u); + deallog << "2 * " << u.value << " = " << v.value << std::endl; + + multiply4.vmult_add(v, u); + deallog << "... + 4 * " << u.value << " = " << v.value << std::endl; + + multiply4.vmult(v, u); + deallog << "4 * " << u.value << " = " << v.value << std::endl; + + multiply2.vmult_add(v, u); + deallog << "... + 2 * " << u.value << " = " << v.value << std::endl; + + // Tvmult, Tvmult_add + + v.value = 4.; + + multiply2.Tvmult(u, v); + deallog << "2 * " << v.value << " = " << u.value << std::endl; + + multiply4.Tvmult_add(u, v); + deallog << "... + 4 * " << v.value << " = " << u.value << std::endl; + + multiply4.Tvmult(u, v); + deallog << "4 * " << v.value << " = " << u.value << std::endl; + + multiply2.Tvmult_add(u, v); + deallog << "... + 2 * " << v.value << " = " << u.value << std::endl; + + // operator+, operator-, operator+=, operator-= + + auto test = multiply2 + multiply4; + test.vmult(v, u); + deallog << "(2 + 4) * " << u.value << " = " << v.value << std::endl; + + test = multiply2 - multiply4; + test.vmult(v, u); + deallog << "(2 - 4) * " << u.value << " = " << v.value << std::endl; + + test += multiply2; + test.vmult(v, u); + deallog << "(2 - 4 + 2) * " << u.value << " = " << v.value << std::endl; + + test -= multiply4; + test.vmult(v, u); + deallog << "(2 - 4 + 2 - 4) * " << u.value << " = " << v.value << std::endl; + + // operator* with scalar + + test = 4. * multiply4; + test.vmult(v, u); + deallog << "(4 * 4) * " << u.value << " = " << v.value << std::endl; + + test = multiply4; + test *= std::complex(4.); + test.vmult(v, u); + deallog << "(4 * 4) * " << u.value << " = " << v.value << std::endl; + + test = multiply4 * 4.; + test.vmult(v, u); + deallog << "(4 * 4) * " << u.value << " = " << v.value << std::endl; + + // operator* and transpose + + auto test2 = transpose_operator(multiply2) * multiply4; + RightVector w = { 0. }; + test2.vmult(w, u); + deallog << "(2 * 4) * " << u.value << " = " << w.value << std::endl; + + test2 *= identity_operator(test2.reinit_range_vector); + test2.vmult(w, u); + deallog << "(2 * 4) * 1 * " << u.value << " = " << w.value << std::endl; + + // identity + + auto test3 = identity_operator(test2.reinit_range_vector) + test2; + test3.vmult(w, u); + deallog << "(1 + 2 * 4) * " << u.value << " = " << w.value << std::endl; + + // null operator + + auto test4 = null_operator(test2); + test4.vmult(w, u); + deallog << " 0 * " << u.value << " = " << w.value << std::endl; + +} diff --git a/tests/lac/linear_operator_01a.with_cxx11=on.output b/tests/lac/linear_operator_01a.with_cxx11=on.output new file mode 100644 index 0000000000..2f70fcac14 --- /dev/null +++ b/tests/lac/linear_operator_01a.with_cxx11=on.output @@ -0,0 +1,20 @@ + +DEAL::2 * (4.00000,0.00000) = (8.00000,0.00000) +DEAL::... + 4 * (4.00000,0.00000) = (24.0000,0.00000) +DEAL::4 * (4.00000,0.00000) = (16.0000,0.00000) +DEAL::... + 2 * (4.00000,0.00000) = (24.0000,0.00000) +DEAL::2 * (4.00000,0.00000) = (8.00000,0.00000) +DEAL::... + 4 * (4.00000,0.00000) = (24.0000,0.00000) +DEAL::4 * (4.00000,0.00000) = (16.0000,0.00000) +DEAL::... + 2 * (4.00000,0.00000) = (24.0000,0.00000) +DEAL::(2 + 4) * (24.0000,0.00000) = (144.000,0.00000) +DEAL::(2 - 4) * (24.0000,0.00000) = (-48.0000,0.00000) +DEAL::(2 - 4 + 2) * (24.0000,0.00000) = (0.00000,0.00000) +DEAL::(2 - 4 + 2 - 4) * (24.0000,0.00000) = (-96.0000,0.00000) +DEAL::(4 * 4) * (24.0000,0.00000) = (384.000,0.00000) +DEAL::(4 * 4) * (24.0000,0.00000) = (384.000,0.00000) +DEAL::(4 * 4) * (24.0000,0.00000) = (384.000,0.00000) +DEAL::(2 * 4) * (24.0000,0.00000) = (192.000,0.00000) +DEAL::(2 * 4) * 1 * (24.0000,0.00000) = (192.000,0.00000) +DEAL::(1 + 2 * 4) * (24.0000,0.00000) = (216.000,0.00000) +DEAL:: 0 * (24.0000,0.00000) = (0.00000,0.00000) diff --git a/tests/lac/linear_operator_02a.cc b/tests/lac/linear_operator_02a.cc new file mode 100644 index 0000000000..4cbfbc6a87 --- /dev/null +++ b/tests/lac/linear_operator_02a.cc @@ -0,0 +1,159 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +// Variant of linear_operator_02 that uses std::complex as +// value_type: Tests for the LinearOperator template with +// dealii::Vector +// dealii::SparseMatrix +// dealii::FullMatrix + +#include "../tests.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + +using namespace dealii; + +int main() +{ + initlog(); + deallog << std::setprecision(10); + + static const int dim = 2; + + Triangulation triangulation; + GridGenerator::hyper_cube (triangulation); + triangulation.refine_global(2); + + MappingQ1 mapping_q1; + FE_Q q1(1); + DoFHandler dof_handler(triangulation); + dof_handler.distribute_dofs(q1); + + DynamicSparsityPattern dsp(dof_handler.n_dofs(), dof_handler.n_dofs()); + DoFTools::make_sparsity_pattern(dof_handler, dsp); + dsp.compress(); + SparsityPattern sparsity_pattern; + sparsity_pattern.copy_from(dsp); + sparsity_pattern.compress(); + + // use a real valued matrix and rely on the template overloads of vmult, + // etc. to also work with std::complex + SparseMatrix a (sparsity_pattern); + SparseMatrix b (sparsity_pattern); + + QGauss quadrature(4); + MatrixCreator::create_laplace_matrix(mapping_q1, dof_handler, quadrature, a); + MatrixCreator::create_mass_matrix(mapping_q1, dof_handler, quadrature, b); + + + // Constructors and assignment: + + auto op_a = linear_operator>>(a); + auto op_b = linear_operator>>(b); + + { + LinearOperator>> op_x (a); + op_a = a; + op_b = b; + } + + // vmult: + + Vector> u; + op_a.reinit_domain_vector(u, true); + for (unsigned int i = 0; i < u.size(); ++i) { + u[i] = (std::complex)(i+1); + } + + deallog << "u: " << u << std::endl; + + Vector> v; + op_a.reinit_domain_vector(v, false); + Vector> w; + op_a.reinit_domain_vector(w, false); + Vector> x; + op_a.reinit_domain_vector(x, false); + + op_a.vmult(v, u); + deallog << "Au: " << v << std::endl; + + op_b.vmult(w, u); + deallog << "Bu: " << w << std::endl; + + // operator+, operator-, operator+=, operator-=: + + x = v; + x += w; + deallog << "Au+Bu: " << x << std::endl; + + (op_a + op_b).vmult(x, u); + deallog << "(A+B)u: " << x << std::endl; + + auto op_x = op_a; + op_x += op_b; + op_x.vmult(x, u); + deallog << "(A+=B)u: " << x << std::endl; + + x = v; + x -= w; + deallog << "Au-Bu: " << x << std::endl; + + (op_a - op_b).vmult(x, u); + deallog << "(A-B)u: " << x << std::endl; + + op_x = op_a; + op_x -= op_b; + op_x.vmult(x, u); + deallog << "(A-=B)u: " << x << std::endl; + + // operator*, operator*= + + op_b.vmult(v,u); + op_a.vmult(w,v); + deallog << "(A(Bu)): " << w << std::endl; + + (op_a * op_b).vmult(x, u); + deallog << "(A*B)u: " << x << std::endl; + + op_x = op_a; + op_x *= op_b; + op_x.vmult(x, u); + deallog << "(A*=B)u: " << x << std::endl; + + op_x *= std::complex(4.); + op_x.vmult(x, u); + deallog << "(A*=B*=4.)u: " << x << std::endl; + + // solver interface: + + // TODO: Well, update as soon as SolverCG can work with complex data + // types. +} diff --git a/tests/lac/linear_operator_02a.with_cxx11=on.output b/tests/lac/linear_operator_02a.with_cxx11=on.output new file mode 100644 index 0000000000..4b35304bf5 --- /dev/null +++ b/tests/lac/linear_operator_02a.with_cxx11=on.output @@ -0,0 +1,27 @@ + +DEAL::u: (1.000000000,0.000000000) (2.000000000,0.000000000) (3.000000000,0.000000000) (4.000000000,0.000000000) (5.000000000,0.000000000) (6.000000000,0.000000000) (7.000000000,0.000000000) (8.000000000,0.000000000) (9.000000000,0.000000000) (10.00000000,0.000000000) (11.00000000,0.000000000) (12.00000000,0.000000000) (13.00000000,0.000000000) (14.00000000,0.000000000) (15.00000000,0.000000000) (16.00000000,0.000000000) (17.00000000,0.000000000) (18.00000000,0.000000000) (19.00000000,0.000000000) (20.00000000,0.000000000) (21.00000000,0.000000000) (22.00000000,0.000000000) (23.00000000,0.000000000) (24.00000000,0.000000000) (25.00000000,0.000000000) +DEAL:: +DEAL::Au: (-1.500000000,0.000000000) (-2.666666667,0.000000000) (-2.000000000,0.000000000) (-3.000000000,0.000000000) (-2.333333333,0.000000000) (-5.000000000,0.000000000) (-3.500000000,0.000000000) (-5.333333333,0.000000000) (-9.333333333,0.000000000) (0.5000000000,0.000000000) (1.333333333,0.000000000) (0.5000000000,0.000000000) (1.166666667,0.000000000) (-1.666666667,0.000000000) (-1.666666667,0.000000000) (2.000000000,0.000000000) (6.000000000,0.000000000) (3.000000000,0.000000000) (1.000000000,0.000000000) (3.000000000,0.000000000) (1.666666667,0.000000000) (9.000000000,0.000000000) (4.000000000,0.000000000) (3.333333333,0.000000000) (1.500000000,0.000000000) +DEAL:: +DEAL::Bu: (0.03125000000,0.000000000) (0.09201388889,0.000000000) (0.1145833333,0.000000000) (0.2812500000,0.000000000) (0.1788194444,0.000000000) (0.4270833333,0.000000000) (0.2552083333,0.000000000) (0.5538194444,0.000000000) (0.6631944444,0.000000000) (0.3072916667,0.000000000) (0.6753472222,0.000000000) (0.1822916667,0.000000000) (0.3923611111,0.000000000) (0.8888888889,0.000000000) (0.4878472222,0.000000000) (0.4791666667,0.000000000) (1.000000000,0.000000000) (1.093750000,0.000000000) (0.2864583333,0.000000000) (0.5937500000,0.000000000) (0.6371527778,0.000000000) (1.281250000,0.000000000) (0.6770833333,0.000000000) (0.7170138889,0.000000000) (0.3750000000,0.000000000) +DEAL:: +DEAL::Au+Bu: (-1.468750000,0.000000000) (-2.574652778,0.000000000) (-1.885416667,0.000000000) (-2.718750000,0.000000000) (-2.154513889,0.000000000) (-4.572916667,0.000000000) (-3.244791667,0.000000000) (-4.779513889,0.000000000) (-8.670138889,0.000000000) (0.8072916667,0.000000000) (2.008680556,0.000000000) (0.6822916667,0.000000000) (1.559027778,0.000000000) (-0.7777777778,0.000000000) (-1.178819444,0.000000000) (2.479166667,0.000000000) (7.000000000,0.000000000) (4.093750000,0.000000000) (1.286458333,0.000000000) (3.593750000,0.000000000) (2.303819444,0.000000000) (10.28125000,0.000000000) (4.677083333,0.000000000) (4.050347222,0.000000000) (1.875000000,0.000000000) +DEAL:: +DEAL::(A+B)u: (-1.468750000,0.000000000) (-2.574652778,0.000000000) (-1.885416667,0.000000000) (-2.718750000,0.000000000) (-2.154513889,0.000000000) (-4.572916667,0.000000000) (-3.244791667,0.000000000) (-4.779513889,0.000000000) (-8.670138889,0.000000000) (0.8072916667,0.000000000) (2.008680556,0.000000000) (0.6822916667,0.000000000) (1.559027778,0.000000000) (-0.7777777778,0.000000000) (-1.178819444,0.000000000) (2.479166667,0.000000000) (7.000000000,0.000000000) (4.093750000,0.000000000) (1.286458333,0.000000000) (3.593750000,0.000000000) (2.303819444,0.000000000) (10.28125000,0.000000000) (4.677083333,0.000000000) (4.050347222,0.000000000) (1.875000000,0.000000000) +DEAL:: +DEAL::(A+=B)u: (-1.468750000,0.000000000) (-2.574652778,0.000000000) (-1.885416667,0.000000000) (-2.718750000,0.000000000) (-2.154513889,0.000000000) (-4.572916667,0.000000000) (-3.244791667,0.000000000) (-4.779513889,0.000000000) (-8.670138889,0.000000000) (0.8072916667,0.000000000) (2.008680556,0.000000000) (0.6822916667,0.000000000) (1.559027778,0.000000000) (-0.7777777778,0.000000000) (-1.178819444,0.000000000) (2.479166667,0.000000000) (7.000000000,0.000000000) (4.093750000,0.000000000) (1.286458333,0.000000000) (3.593750000,0.000000000) (2.303819444,0.000000000) (10.28125000,0.000000000) (4.677083333,0.000000000) (4.050347222,0.000000000) (1.875000000,0.000000000) +DEAL:: +DEAL::Au-Bu: (-1.531250000,0.000000000) (-2.758680556,0.000000000) (-2.114583333,0.000000000) (-3.281250000,0.000000000) (-2.512152778,0.000000000) (-5.427083333,0.000000000) (-3.755208333,0.000000000) (-5.887152778,0.000000000) (-9.996527778,0.000000000) (0.1927083333,0.000000000) (0.6579861111,0.000000000) (0.3177083333,0.000000000) (0.7743055556,0.000000000) (-2.555555556,0.000000000) (-2.154513889,0.000000000) (1.520833333,0.000000000) (5.000000000,0.000000000) (1.906250000,0.000000000) (0.7135416667,0.000000000) (2.406250000,0.000000000) (1.029513889,0.000000000) (7.718750000,0.000000000) (3.322916667,0.000000000) (2.616319444,0.000000000) (1.125000000,0.000000000) +DEAL:: +DEAL::(A-B)u: (-1.531250000,0.000000000) (-2.758680556,0.000000000) (-2.114583333,0.000000000) (-3.281250000,0.000000000) (-2.512152778,0.000000000) (-5.427083333,0.000000000) (-3.755208333,0.000000000) (-5.887152778,0.000000000) (-9.996527778,0.000000000) (0.1927083333,0.000000000) (0.6579861111,0.000000000) (0.3177083333,0.000000000) (0.7743055556,0.000000000) (-2.555555556,0.000000000) (-2.154513889,0.000000000) (1.520833333,0.000000000) (5.000000000,0.000000000) (1.906250000,0.000000000) (0.7135416667,0.000000000) (2.406250000,0.000000000) (1.029513889,0.000000000) (7.718750000,0.000000000) (3.322916667,0.000000000) (2.616319444,0.000000000) (1.125000000,0.000000000) +DEAL:: +DEAL::(A-=B)u: (-1.531250000,0.000000000) (-2.758680556,0.000000000) (-2.114583333,0.000000000) (-3.281250000,0.000000000) (-2.512152778,0.000000000) (-5.427083333,0.000000000) (-3.755208333,0.000000000) (-5.887152778,0.000000000) (-9.996527778,0.000000000) (0.1927083333,0.000000000) (0.6579861111,0.000000000) (0.3177083333,0.000000000) (0.7743055556,0.000000000) (-2.555555556,0.000000000) (-2.154513889,0.000000000) (1.520833333,0.000000000) (5.000000000,0.000000000) (1.906250000,0.000000000) (0.7135416667,0.000000000) (2.406250000,0.000000000) (1.029513889,0.000000000) (7.718750000,0.000000000) (3.322916667,0.000000000) (2.616319444,0.000000000) (1.125000000,0.000000000) +DEAL:: +DEAL::(A(Bu)): (-0.1073495370,0.000000000) (-0.1866319444,0.000000000) (-0.2039930556,0.000000000) (-0.02199074074,0.000000000) (-0.2893518519,0.000000000) (-0.07465277778,0.000000000) (-0.3703703704,0.000000000) (0.03877314815,0.000000000) (-0.2986111111,0.000000000) (-0.1487268519,0.000000000) (0.6250000000,0.000000000) (-0.2201967593,0.000000000) (-0.2123842593,0.000000000) (0.4710648148,0.000000000) (-0.4762731481,0.000000000) (-0.1672453704,0.000000000) (1.145833333,0.000000000) (0.8049768519,0.000000000) (-0.3211805556,0.000000000) (-0.2199074074,0.000000000) (-0.4939236111,0.000000000) (1.570023148,0.000000000) (-0.2034143519,0.000000000) (-0.2300347222,0.000000000) (-0.4094328704,0.000000000) +DEAL:: +DEAL::(A*B)u: (-0.1073495370,0.000000000) (-0.1866319444,0.000000000) (-0.2039930556,0.000000000) (-0.02199074074,0.000000000) (-0.2893518519,0.000000000) (-0.07465277778,0.000000000) (-0.3703703704,0.000000000) (0.03877314815,0.000000000) (-0.2986111111,0.000000000) (-0.1487268519,0.000000000) (0.6250000000,0.000000000) (-0.2201967593,0.000000000) (-0.2123842593,0.000000000) (0.4710648148,0.000000000) (-0.4762731481,0.000000000) (-0.1672453704,0.000000000) (1.145833333,0.000000000) (0.8049768519,0.000000000) (-0.3211805556,0.000000000) (-0.2199074074,0.000000000) (-0.4939236111,0.000000000) (1.570023148,0.000000000) (-0.2034143519,0.000000000) (-0.2300347222,0.000000000) (-0.4094328704,0.000000000) +DEAL:: +DEAL::(A*=B)u: (-0.1073495370,0.000000000) (-0.1866319444,0.000000000) (-0.2039930556,0.000000000) (-0.02199074074,0.000000000) (-0.2893518519,0.000000000) (-0.07465277778,0.000000000) (-0.3703703704,0.000000000) (0.03877314815,0.000000000) (-0.2986111111,0.000000000) (-0.1487268519,0.000000000) (0.6250000000,0.000000000) (-0.2201967593,0.000000000) (-0.2123842593,0.000000000) (0.4710648148,0.000000000) (-0.4762731481,0.000000000) (-0.1672453704,0.000000000) (1.145833333,0.000000000) (0.8049768519,0.000000000) (-0.3211805556,0.000000000) (-0.2199074074,0.000000000) (-0.4939236111,0.000000000) (1.570023148,0.000000000) (-0.2034143519,0.000000000) (-0.2300347222,0.000000000) (-0.4094328704,0.000000000) +DEAL:: +DEAL::(A*=B*=4.)u: (-0.4293981481,0.000000000) (-0.7465277778,0.000000000) (-0.8159722222,0.000000000) (-0.08796296296,0.000000000) (-1.157407407,0.000000000) (-0.2986111111,0.000000000) (-1.481481481,0.000000000) (0.1550925926,0.000000000) (-1.194444444,0.000000000) (-0.5949074074,0.000000000) (2.500000000,0.000000000) (-0.8807870370,0.000000000) (-0.8495370370,0.000000000) (1.884259259,0.000000000) (-1.905092593,0.000000000) (-0.6689814815,0.000000000) (4.583333333,0.000000000) (3.219907407,0.000000000) (-1.284722222,0.000000000) (-0.8796296296,0.000000000) (-1.975694444,0.000000000) (6.280092593,0.000000000) (-0.8136574074,0.000000000) (-0.9201388889,0.000000000) (-1.637731481,0.000000000) +DEAL:: -- 2.39.5