From 6b97567de9d45b6df8032b36cac2500f97045f74 Mon Sep 17 00:00:00 2001 From: Jihuan Tian Date: Fri, 22 Apr 2016 07:05:08 +0800 Subject: [PATCH] 1. Bug fix for wrong sign of curl integrators in integrators/maxwell.h: nitsche_curl_matrix, curl_curl_matrix and curl_matrix. 2. Add test code 'tests/integrators/maxwell_curl.cc' and related output file for verifying curl integrators in integrators/maxwell.h. 3. Add basis function support checking on faces in the function nitsche_curl_matrix in integrators/maxwell.h. --- doc/news/changes.h | 10 ++- include/deal.II/integrators/maxwell.h | 46 ++++++----- tests/integrators/maxwell_curl.cc | 109 ++++++++++++++++++++++++++ tests/integrators/maxwell_curl.output | 76 ++++++++++++++++++ 4 files changed, 219 insertions(+), 22 deletions(-) create mode 100644 tests/integrators/maxwell_curl.cc create mode 100644 tests/integrators/maxwell_curl.output diff --git a/doc/news/changes.h b/doc/news/changes.h index c62ea5a10d..c3780ac9c3 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -241,6 +241,14 @@ inconvenience this causes.

Specific improvements

    +
  1. Fixed: Corrected the sign of curl calculated in the functions: + LocalIntegrators::curl_curl_matrix, LocalIntegrators::curl_matrix, + LocalIntegrators::nitsche_curl_matrix and LocalIntegrators::ip_curl_matrix in + integrators/maxwell.h. +
    + (Jihuan Tian, 2016/05/09) +
  2. +
  3. Fixed: Bug in the RelaxationBlock class function do_step. Before, the corrections were not added together, which leads to a wrong update whenever the Jacobi blocks are overlapping. For SOR, SSOR and non-overlapping Jacobi this was @@ -287,7 +295,7 @@ inconvenience this causes.
    (Martin Kronbichler, Daniel Jodlbauer, 2016/04/21)
  4. - +
  5. New: Added an optional string parameter to the ParameterHandler::read_input () and ParameterHandler::read_input_from_string() functions. When a line which equals this string is encountered, the parsing of parameters diff --git a/include/deal.II/integrators/maxwell.h b/include/deal.II/integrators/maxwell.h index 64432a9ee8..e7e08589c5 100644 --- a/include/deal.II/integrators/maxwell.h +++ b/include/deal.II/integrators/maxwell.h @@ -37,9 +37,9 @@ namespace LocalIntegrators * * @f[ * \nabla\times \mathbf u = \begin{pmatrix} - * \partial_3 u_2 - \partial_2 u_3 \\ - * \partial_1 u_3 - \partial_3 u_1 \\ - * \partial_2 u_1 - \partial_1 u_2 + * \partial_2 u_3 - \partial_3 u_2 \\ + * \partial_3 u_1 - \partial_1 u_3 \\ + * \partial_1 u_2 - \partial_2 u_1 * \end{pmatrix}. * @f] * @@ -49,7 +49,7 @@ namespace LocalIntegrators * function and the vector curl of a scalar function. The current * implementation exchanges the sign and we have: * @f[ - * \nabla \times \mathbf u = \partial_1 u_2 - \partial 2 u_1, + * \nabla \times \mathbf u = \partial_1 u_2 - \partial_2 u_1, * \qquad * \nabla \times p = \begin{pmatrix} * \partial_2 p \\ -\partial_1 p @@ -197,8 +197,8 @@ namespace LocalIntegrators const unsigned int d1 = (d+1)%dim; const unsigned int d2 = (d+2)%dim; - const double cv = fe.shape_grad_component(i,k,d1)[d2] - fe.shape_grad_component(i,k,d2)[d1]; - const double cu = fe.shape_grad_component(j,k,d1)[d2] - fe.shape_grad_component(j,k,d2)[d1]; + const double cv = fe.shape_grad_component(i,k,d2)[d1] - fe.shape_grad_component(i,k,d1)[d2]; + const double cu = fe.shape_grad_component(j,k,d2)[d1] - fe.shape_grad_component(j,k,d1)[d2]; M(i,j) += dx * cu * cv; } @@ -246,7 +246,7 @@ namespace LocalIntegrators const unsigned int d2 = (d+2)%dim; const double vv = fetest.shape_value_component(i,k,d); - const double cu = fe.shape_grad_component(j,k,d1)[d2] - fe.shape_grad_component(j,k,d2)[d1]; + const double cu = fe.shape_grad_component(j,k,d2)[d1] - fe.shape_grad_component(j,k,d1)[d2]; M(i,j) += dx * cu * vv; } } @@ -272,6 +272,7 @@ namespace LocalIntegrators void nitsche_curl_matrix ( FullMatrix &M, const FEValuesBase &fe, + const unsigned int face_no, double penalty, double factor = 1.) { @@ -299,17 +300,20 @@ namespace LocalIntegrators const Tensor<1,dim> n = fe.normal_vector(k); for (unsigned int i=0; i +// +// 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. +// +// --------------------------------------------------------------------- + + +// Test curl related functions in integrators/maxwell.h + +#include +#include +#include +#include + +#include +#include +#include +#include + +#include +#include +#include + +#include +#include +#include + +#include + +#include + +using namespace dealii; +using namespace LocalIntegrators::Maxwell; + +template +void make_grid(Triangulation &tr) +{ + GridGenerator::subdivided_hyper_cube(tr, 1, -1., 1.); +} + +template +void TestMaxwellCurl(Triangulation &tr) +{ + int element_order = 1; + int quadrature_order = 3; + + DoFHandler dof_handler(tr); + FESystem fe(FE_Q(element_order), dim); + FEValues fe_values(fe, QGauss(quadrature_order), update_values | update_JxW_values | update_gradients | update_quadrature_points); + FEFaceValues fe_face_values(fe, QGauss(quadrature_order), update_values | update_JxW_values | update_gradients | update_normal_vectors | update_quadrature_points); + + dof_handler.distribute_dofs(fe); + + const unsigned int dofs_per_cell = fe.dofs_per_cell; + + FullMatrix curl_curl_check(dofs_per_cell, dofs_per_cell); + FullMatrix curl_check(dofs_per_cell, dofs_per_cell); + FullMatrix nitsche_curl_check(dofs_per_cell, dofs_per_cell); + + curl_curl_check = 0; + curl_check = 0; + nitsche_curl_check = 0; + + typename::DoFHandler::cell_iterator cell = dof_handler.begin(0); + + fe_values.reinit(cell); + + curl_curl_matrix(curl_curl_check, fe_values, 1.); + deallog << "curl_curl_matrix" << std::endl; + curl_curl_check.print(deallog, 10); + + curl_matrix(curl_check, fe_values, fe_values, 1.); + deallog << "curl_matrix" << std::endl; + curl_check.print(deallog, 10); + + for (unsigned int face = 0; face < GeometryInfo::faces_per_cell; ++face) + { + fe_face_values.reinit(cell, face); + nitsche_curl_matrix(nitsche_curl_check, fe_face_values, face, 1., 1.); + } + + deallog << "nitsche_curl_matrix" << std::endl; + nitsche_curl_check.print(deallog, 10); + + dof_handler.clear(); +} + +int main() +{ + const std::string logname = "output"; + std::ofstream logfile(logname.c_str()); + deallog.attach(logfile); + deallog.threshold_double(1.e-10); + + Triangulation<3> tr; + make_grid(tr); + TestMaxwellCurl(tr); +} diff --git a/tests/integrators/maxwell_curl.output b/tests/integrators/maxwell_curl.output new file mode 100644 index 0000000000..30f6ffd82c --- /dev/null +++ b/tests/integrators/maxwell_curl.output @@ -0,0 +1,76 @@ + +DEAL::curl_curl_matrix +DEAL::0.44 -0.17 -0.17 0.22 0.17 0.17 -0.11 -0.17 -0.083 -0.056 0.17 0.083 -0.11 -0.083 -0.17 -0.056 0.083 0.17 -0.22 -0.083 -0.083 -0.11 0.083 0.083 +DEAL::-0.17 0.44 -0.17 -0.17 -0.11 -0.083 0.17 0.22 0.17 0.17 -0.056 0.083 -0.083 -0.11 -0.17 -0.083 -0.22 -0.083 0.083 -0.056 0.17 0.083 -0.11 0.083 +DEAL::-0.17 -0.17 0.44 -0.17 -0.083 -0.11 -0.083 -0.17 -0.11 -0.083 -0.083 -0.22 0.17 0.17 0.22 0.17 0.083 -0.056 0.083 0.17 -0.056 0.083 0.083 -0.11 +DEAL::0.22 -0.17 -0.17 0.44 0.17 0.17 -0.056 -0.17 -0.083 -0.11 0.17 0.083 -0.056 -0.083 -0.17 -0.11 0.083 0.17 -0.11 -0.083 -0.083 -0.22 0.083 0.083 +DEAL::0.17 -0.11 -0.083 0.17 0.44 -0.17 -0.17 -0.056 0.083 -0.17 0.22 0.17 0.083 -0.22 -0.083 0.083 -0.11 -0.17 -0.083 -0.11 0.083 -0.083 -0.056 0.17 +DEAL::0.17 -0.083 -0.11 0.17 -0.17 0.44 0.083 -0.083 -0.22 0.083 -0.17 -0.11 -0.17 0.083 -0.056 -0.17 0.17 0.22 -0.083 0.083 -0.11 -0.083 0.17 -0.056 +DEAL::-0.11 0.17 -0.083 -0.056 -0.17 0.083 0.44 0.17 -0.17 0.22 -0.17 0.17 -0.22 0.083 -0.083 -0.11 -0.083 0.083 -0.11 0.083 -0.17 -0.056 -0.083 0.17 +DEAL::-0.17 0.22 -0.17 -0.17 -0.056 -0.083 0.17 0.44 0.17 0.17 -0.11 0.083 -0.083 -0.056 -0.17 -0.083 -0.11 -0.083 0.083 -0.11 0.17 0.083 -0.22 0.083 +DEAL::-0.083 0.17 -0.11 -0.083 0.083 -0.22 -0.17 0.17 0.44 -0.17 0.083 -0.11 0.083 -0.17 -0.056 0.083 -0.083 -0.11 0.17 -0.17 0.22 0.17 -0.083 -0.056 +DEAL::-0.056 0.17 -0.083 -0.11 -0.17 0.083 0.22 0.17 -0.17 0.44 -0.17 0.17 -0.11 0.083 -0.083 -0.22 -0.083 0.083 -0.056 0.083 -0.17 -0.11 -0.083 0.17 +DEAL::0.17 -0.056 -0.083 0.17 0.22 -0.17 -0.17 -0.11 0.083 -0.17 0.44 0.17 0.083 -0.11 -0.083 0.083 -0.056 -0.17 -0.083 -0.22 0.083 -0.083 -0.11 0.17 +DEAL::0.083 0.083 -0.22 0.083 0.17 -0.11 0.17 0.083 -0.11 0.17 0.17 0.44 -0.083 -0.083 -0.11 -0.083 -0.17 -0.056 -0.17 -0.083 -0.056 -0.17 -0.17 0.22 +DEAL::-0.11 -0.083 0.17 -0.056 0.083 -0.17 -0.22 -0.083 0.083 -0.11 0.083 -0.083 0.44 -0.17 0.17 0.22 0.17 -0.17 -0.11 -0.17 0.083 -0.056 0.17 -0.083 +DEAL::-0.083 -0.11 0.17 -0.083 -0.22 0.083 0.083 -0.056 -0.17 0.083 -0.11 -0.083 -0.17 0.44 0.17 -0.17 -0.11 0.083 0.17 0.22 -0.17 0.17 -0.056 -0.083 +DEAL::-0.17 -0.17 0.22 -0.17 -0.083 -0.056 -0.083 -0.17 -0.056 -0.083 -0.083 -0.11 0.17 0.17 0.44 0.17 0.083 -0.11 0.083 0.17 -0.11 0.083 0.083 -0.22 +DEAL::-0.056 -0.083 0.17 -0.11 0.083 -0.17 -0.11 -0.083 0.083 -0.22 0.083 -0.083 0.22 -0.17 0.17 0.44 0.17 -0.17 -0.056 -0.17 0.083 -0.11 0.17 -0.083 +DEAL::0.083 -0.22 0.083 0.083 -0.11 0.17 -0.083 -0.11 -0.083 -0.083 -0.056 -0.17 0.17 -0.11 0.083 0.17 0.44 0.17 -0.17 -0.056 -0.083 -0.17 0.22 -0.17 +DEAL::0.17 -0.083 -0.056 0.17 -0.17 0.22 0.083 -0.083 -0.11 0.083 -0.17 -0.056 -0.17 0.083 -0.11 -0.17 0.17 0.44 -0.083 0.083 -0.22 -0.083 0.17 -0.11 +DEAL::-0.22 0.083 0.083 -0.11 -0.083 -0.083 -0.11 0.083 0.17 -0.056 -0.083 -0.17 -0.11 0.17 0.083 -0.056 -0.17 -0.083 0.44 0.17 0.17 0.22 -0.17 -0.17 +DEAL::-0.083 -0.056 0.17 -0.083 -0.11 0.083 0.083 -0.11 -0.17 0.083 -0.22 -0.083 -0.17 0.22 0.17 -0.17 -0.056 0.083 0.17 0.44 -0.17 0.17 -0.11 -0.083 +DEAL::-0.083 0.17 -0.056 -0.083 0.083 -0.11 -0.17 0.17 0.22 -0.17 0.083 -0.056 0.083 -0.17 -0.11 0.083 -0.083 -0.22 0.17 -0.17 0.44 0.17 -0.083 -0.11 +DEAL::-0.11 0.083 0.083 -0.22 -0.083 -0.083 -0.056 0.083 0.17 -0.11 -0.083 -0.17 -0.056 0.17 0.083 -0.11 -0.17 -0.083 0.22 0.17 0.17 0.44 -0.17 -0.17 +DEAL::0.083 -0.11 0.083 0.083 -0.056 0.17 -0.083 -0.22 -0.083 -0.083 -0.11 -0.17 0.17 -0.056 0.083 0.17 0.22 0.17 -0.17 -0.11 -0.083 -0.17 0.44 -0.17 +DEAL::0.083 0.083 -0.11 0.083 0.17 -0.056 0.17 0.083 -0.056 0.17 0.17 0.22 -0.083 -0.083 -0.22 -0.083 -0.17 -0.11 -0.17 -0.083 -0.11 -0.17 -0.17 0.44 +DEAL::curl_matrix +DEAL::0 0.22 -0.22 0 0.11 -0.11 0 0.11 0.22 0 0.056 0.11 0 -0.22 -0.11 0 -0.11 -0.056 0 -0.11 0.11 0 -0.056 0.056 +DEAL::-0.22 0 0.22 -0.11 0 -0.22 -0.11 0 0.11 -0.056 0 -0.11 0.22 0 0.11 0.11 0 -0.11 0.11 0 0.056 0.056 0 -0.056 +DEAL::0.22 -0.22 0 0.11 0.22 0 -0.22 -0.11 0 -0.11 0.11 0 0.11 -0.11 0 0.056 0.11 0 -0.11 -0.056 0 -0.056 0.056 0 +DEAL::0 0.11 -0.11 0 0.22 -0.22 0 0.056 0.11 0 0.11 0.22 0 -0.11 -0.056 0 -0.22 -0.11 0 -0.056 0.056 0 -0.11 0.11 +DEAL::-0.11 0 0.22 -0.22 0 -0.22 -0.056 0 0.11 -0.11 0 -0.11 0.11 0 0.11 0.22 0 -0.11 0.056 0 0.056 0.11 0 -0.056 +DEAL::0.11 -0.22 0 0.22 0.22 0 -0.11 -0.11 0 -0.22 0.11 0 0.056 -0.11 0 0.11 0.11 0 -0.056 -0.056 0 -0.11 0.056 0 +DEAL::0 0.11 -0.22 0 0.056 -0.11 0 0.22 0.22 0 0.11 0.11 0 -0.11 -0.11 0 -0.056 -0.056 0 -0.22 0.11 0 -0.11 0.056 +DEAL::-0.11 0 0.11 -0.056 0 -0.11 -0.22 0 0.22 -0.11 0 -0.22 0.11 0 0.056 0.056 0 -0.056 0.22 0 0.11 0.11 0 -0.11 +DEAL::0.22 -0.11 0 0.11 0.11 0 -0.22 -0.22 0 -0.11 0.22 0 0.11 -0.056 0 0.056 0.056 0 -0.11 -0.11 0 -0.056 0.11 0 +DEAL::0 0.056 -0.11 0 0.11 -0.22 0 0.11 0.11 0 0.22 0.22 0 -0.056 -0.056 0 -0.11 -0.11 0 -0.11 0.056 0 -0.22 0.11 +DEAL::-0.056 0 0.11 -0.11 0 -0.11 -0.11 0 0.22 -0.22 0 -0.22 0.056 0 0.056 0.11 0 -0.056 0.11 0 0.11 0.22 0 -0.11 +DEAL::0.11 -0.11 0 0.22 0.11 0 -0.11 -0.22 0 -0.22 0.22 0 0.056 -0.056 0 0.11 0.056 0 -0.056 -0.11 0 -0.11 0.11 0 +DEAL::0 0.22 -0.11 0 0.11 -0.056 0 0.11 0.11 0 0.056 0.056 0 -0.22 -0.22 0 -0.11 -0.11 0 -0.11 0.22 0 -0.056 0.11 +DEAL::-0.22 0 0.11 -0.11 0 -0.11 -0.11 0 0.056 -0.056 0 -0.056 0.22 0 0.22 0.11 0 -0.22 0.11 0 0.11 0.056 0 -0.11 +DEAL::0.11 -0.11 0 0.056 0.11 0 -0.11 -0.056 0 -0.056 0.056 0 0.22 -0.22 0 0.11 0.22 0 -0.22 -0.11 0 -0.11 0.11 0 +DEAL::0 0.11 -0.056 0 0.22 -0.11 0 0.056 0.056 0 0.11 0.11 0 -0.11 -0.11 0 -0.22 -0.22 0 -0.056 0.11 0 -0.11 0.22 +DEAL::-0.11 0 0.11 -0.22 0 -0.11 -0.056 0 0.056 -0.11 0 -0.056 0.11 0 0.22 0.22 0 -0.22 0.056 0 0.11 0.11 0 -0.11 +DEAL::0.056 -0.11 0 0.11 0.11 0 -0.056 -0.056 0 -0.11 0.056 0 0.11 -0.22 0 0.22 0.22 0 -0.11 -0.11 0 -0.22 0.11 0 +DEAL::0 0.11 -0.11 0 0.056 -0.056 0 0.22 0.11 0 0.11 0.056 0 -0.11 -0.22 0 -0.056 -0.11 0 -0.22 0.22 0 -0.11 0.11 +DEAL::-0.11 0 0.056 -0.056 0 -0.056 -0.22 0 0.11 -0.11 0 -0.11 0.11 0 0.11 0.056 0 -0.11 0.22 0 0.22 0.11 0 -0.22 +DEAL::0.11 -0.056 0 0.056 0.056 0 -0.11 -0.11 0 -0.056 0.11 0 0.22 -0.11 0 0.11 0.11 0 -0.22 -0.22 0 -0.11 0.22 0 +DEAL::0 0.056 -0.056 0 0.11 -0.11 0 0.11 0.056 0 0.22 0.11 0 -0.056 -0.11 0 -0.11 -0.22 0 -0.11 0.11 0 -0.22 0.22 +DEAL::-0.056 0 0.056 -0.11 0 -0.056 -0.11 0 0.11 -0.22 0 -0.11 0.056 0 0.11 0.11 0 -0.11 0.11 0 0.22 0.22 0 -0.22 +DEAL::0.056 -0.056 0 0.11 0.056 0 -0.056 -0.11 0 -0.11 0.11 0 0.11 -0.11 0 0.22 0.11 0 -0.11 -0.22 0 -0.22 0.22 0 +DEAL::nitsche_curl_matrix +DEAL::2.7 -0.67 -0.67 1.3 0.33 0.33 0.67 -0.33 -0.33 0.33 0 0.17 0.67 -0.33 -0.33 0.33 0.17 0 0 -0.17 -0.17 0 0 0 +DEAL::-0.67 2.7 -0.67 -0.33 0.67 -0.33 0.33 1.3 0.33 0 0.33 0.17 -0.33 0.67 -0.33 -0.17 0 -0.17 0.17 0.33 0 0 0 0 +DEAL::-0.67 -0.67 2.7 -0.33 -0.33 0.67 -0.33 -0.33 0.67 -0.17 -0.17 0 0.33 0.33 1.3 0 0.17 0.33 0.17 0 0.33 0 0 0 +DEAL::1.3 -0.33 -0.33 2.7 0.67 0.67 0.33 0 -0.17 0.67 0.33 0.33 0.33 -0.17 0 0.67 0.33 0.33 0 0 0 0 0.17 0.17 +DEAL::0.33 0.67 -0.33 0.67 2.7 -0.67 0 0.33 0.17 -0.33 1.3 0.33 0.17 0 -0.17 0.33 0.67 -0.33 0 0 0 -0.17 0.33 0 +DEAL::0.33 -0.33 0.67 0.67 -0.67 2.7 0.17 -0.17 0 0.33 -0.33 0.67 0 0.17 0.33 -0.33 0.33 1.3 0 0 0 -0.17 0 0.33 +DEAL::0.67 0.33 -0.33 0.33 0 0.17 2.7 0.67 -0.67 1.3 -0.33 0.33 0 0.17 -0.17 0 0 0 0.67 0.33 -0.33 0.33 -0.17 0 +DEAL::-0.33 1.3 -0.33 0 0.33 -0.17 0.67 2.7 0.67 0.33 0.67 0.33 -0.17 0.33 0 0 0 0 0.33 0.67 0.33 0.17 0 0.17 +DEAL::-0.33 0.33 0.67 -0.17 0.17 0 -0.67 0.67 2.7 -0.33 0.33 0.67 0.17 0 0.33 0 0 0 0.33 -0.33 1.3 0 -0.17 0.33 +DEAL::0.33 0 -0.17 0.67 -0.33 0.33 1.3 0.33 -0.33 2.7 -0.67 0.67 0 0 0 0 -0.17 0.17 0.33 0.17 0 0.67 -0.33 0.33 +DEAL::0 0.33 -0.17 0.33 1.3 -0.33 -0.33 0.67 0.33 -0.67 2.7 0.67 0 0 0 0.17 0.33 0 -0.17 0 0.17 -0.33 0.67 0.33 +DEAL::0.17 0.17 0 0.33 0.33 0.67 0.33 0.33 0.67 0.67 0.67 2.7 0 0 0 -0.17 0 0.33 0 -0.17 0.33 -0.33 -0.33 1.3 +DEAL::0.67 -0.33 0.33 0.33 0.17 0 0 -0.17 0.17 0 0 0 2.7 -0.67 0.67 1.3 0.33 -0.33 0.67 -0.33 0.33 0.33 0 -0.17 +DEAL::-0.33 0.67 0.33 -0.17 0 0.17 0.17 0.33 0 0 0 0 -0.67 2.7 0.67 -0.33 0.67 0.33 0.33 1.3 -0.33 0 0.33 -0.17 +DEAL::-0.33 -0.33 1.3 0 -0.17 0.33 -0.17 0 0.33 0 0 0 0.67 0.67 2.7 0.33 0.33 0.67 0.33 0.33 0.67 0.17 0.17 0 +DEAL::0.33 -0.17 0 0.67 0.33 -0.33 0 0 0 0 0.17 -0.17 1.3 -0.33 0.33 2.7 0.67 -0.67 0.33 0 0.17 0.67 0.33 -0.33 +DEAL::0.17 0 0.17 0.33 0.67 0.33 0 0 0 -0.17 0.33 0 0.33 0.67 0.33 0.67 2.7 0.67 0 0.33 -0.17 -0.33 1.3 -0.33 +DEAL::0 -0.17 0.33 0.33 -0.33 1.3 0 0 0 0.17 0 0.33 -0.33 0.33 0.67 -0.67 0.67 2.7 -0.17 0.17 0 -0.33 0.33 0.67 +DEAL::0 0.17 0.17 0 0 0 0.67 0.33 0.33 0.33 -0.17 0 0.67 0.33 0.33 0.33 0 -0.17 2.7 0.67 0.67 1.3 -0.33 -0.33 +DEAL::-0.17 0.33 0 0 0 0 0.33 0.67 -0.33 0.17 0 -0.17 -0.33 1.3 0.33 0 0.33 0.17 0.67 2.7 -0.67 0.33 0.67 -0.33 +DEAL::-0.17 0 0.33 0 0 0 -0.33 0.33 1.3 0 0.17 0.33 0.33 -0.33 0.67 0.17 -0.17 0 0.67 -0.67 2.7 0.33 -0.33 0.67 +DEAL::0 0 0 0 -0.17 -0.17 0.33 0.17 0 0.67 -0.33 -0.33 0.33 0 0.17 0.67 -0.33 -0.33 1.3 0.33 0.33 2.7 -0.67 -0.67 +DEAL::0 0 0 0.17 0.33 0 -0.17 0 -0.17 -0.33 0.67 -0.33 0 0.33 0.17 0.33 1.3 0.33 -0.33 0.67 -0.33 -0.67 2.7 -0.67 +DEAL::0 0 0 0.17 0 0.33 0 0.17 0.33 0.33 0.33 1.3 -0.17 -0.17 0 -0.33 -0.33 0.67 -0.33 -0.33 0.67 -0.67 -0.67 2.7 -- 2.39.5