From e6beda99df04a380b0c60be306763c7d8d2bdb65 Mon Sep 17 00:00:00 2001 From: bangerth Date: Thu, 26 Apr 2012 16:11:48 +0000 Subject: [PATCH] Provide a second function that computes the support point locations even for distributed triangulations. git-svn-id: https://svn.dealii.org/trunk@25455 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 6 + deal.II/include/deal.II/dofs/dof_tools.h | 56 +++++++- deal.II/source/dofs/dof_tools.cc | 126 +++++++++++++----- deal.II/source/dofs/dof_tools.inst.in | 31 +++++ tests/mpi/map_dofs_to_support_points.cc | 115 ++++++++++++++++ .../ncpu_10/cmp/generic | 69 ++++++++++ .../ncpu_4/cmp/generic | 87 ++++++++++++ 7 files changed, 457 insertions(+), 33 deletions(-) create mode 100644 tests/mpi/map_dofs_to_support_points.cc create mode 100644 tests/mpi/map_dofs_to_support_points/ncpu_10/cmp/generic create mode 100644 tests/mpi/map_dofs_to_support_points/ncpu_4/cmp/generic diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index fbbccfd87e..eb9b381013 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -175,6 +175,12 @@ enabled due to a missing include file in file

Specific improvements

    +
  1. New: There is now a second DoFTools::map_dofs_to_support_points +function that also works for parallel::distributed::Triangulation +triangulations. +
    +(Wolfgang Bangerth, 2012/04/26) +
  2. New: There is now a second DoFTools::extract_boundary_dofs function that also works for parallel::distributed::Triangulation triangulations. diff --git a/deal.II/include/deal.II/dofs/dof_tools.h b/deal.II/include/deal.II/dofs/dof_tools.h index 27014eeb90..3db11e7135 100644 --- a/deal.II/include/deal.II/dofs/dof_tools.h +++ b/deal.II/include/deal.II/dofs/dof_tools.h @@ -2251,9 +2251,19 @@ namespace DoFTools * or the like. Otherwise, an * exception is thrown. * - * The given array must have a + * @pre The given array must have a * length of as many elements as * there are degrees of freedom. + * + * @note The precondition to this function + * that the output argument needs to have + * size equal to the total number of degrees + * of freedom makes this function + * unsuitable for the case that the given + * DoFHandler object derives from a + * parallel::distributed::Triangulation object. + * Consequently, this function will produce an + * error if called with such a DoFHandler. */ template void @@ -2262,15 +2272,55 @@ namespace DoFTools std::vector > &support_points); /** - * Same as above for the hp case. + * Same as the previous function but for the hp case. */ - template void map_dofs_to_support_points (const dealii::hp::MappingCollection &mapping, const hp::DoFHandler &dof_handler, std::vector > &support_points); + /** + * This function is a version of the above map_dofs_to_support_points + * function that doesn't simply return a vector of support_points with one + * entry for each global degree of freedom, but instead a map that + * maps from the DoFs index to its location. The point of this + * function is that it is also usable in cases where the DoFHandler + * is based on a parallel::distributed::Triangulation object. In such cases, + * each processor will not be able to determine the support point location + * of all DoFs, and worse no processor may be able to hold a vector that + * would contain the locations of all DoFs even if they were known. As + * a consequence, this function constructs a map from those DoFs for which + * we can know the locations (namely, those DoFs that are + * locally relevant (see @ref GlossLocallyRelevantDof "locally relevant DoFs") + * to their locations. + * + * For non-distributed triangulations, the map returned as @p support_points + * is of course dense, i.e., every DoF is to be found in it. + * + * @param mapping The mapping from the reference cell to the real cell on + * which DoFs are defined. + * @param dof_handler The object that describes which DoF indices live on + * which cell of the triangulation. + * @param support_points A map that for every locally relevant DoF index + * contains the corresponding location in real space coordinates. + * Previous content of this object is deleted in this function. + */ + template + void + map_dofs_to_support_points (const Mapping &mapping, + const DoFHandler &dof_handler, + std::map > &support_points); + + /** + * Same as the previous function but for the hp case. + */ + template + void + map_dofs_to_support_points (const dealii::hp::MappingCollection &mapping, + const hp::DoFHandler &dof_handler, + std::map > &support_points); + /** * This is the opposite function * to the one above. It generates diff --git a/deal.II/source/dofs/dof_tools.cc b/deal.II/source/dofs/dof_tools.cc index b5eb53fe81..b6efd20146 100644 --- a/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/source/dofs/dof_tools.cc @@ -5743,14 +5743,12 @@ namespace DoFTools { namespace { - template - void map_dofs_to_support_points_internal( - const hp::MappingCollection & mapping, - const DH &dof_handler, - std::vector > &support_points) + template + void + map_dofs_to_support_points(const hp::MappingCollection & mapping, + const DH &dof_handler, + std::map > &support_points) { - AssertDimension(support_points.size(), dof_handler.n_dofs()); - const unsigned int dim = DH::dimension; const unsigned int spacedim = DH::space_dimension; @@ -5762,7 +5760,7 @@ namespace DoFTools // check whether every fe in the collection // has support points Assert(fe_collection[fe_index].has_support_points(), - typename FiniteElement::ExcFEHasNoSupportPoints()); + typename FiniteElement::ExcFEHasNoSupportPoints()); q_coll_dummy.push_back( Quadrature ( fe_collection[fe_index].get_unit_support_points())); @@ -5781,25 +5779,49 @@ namespace DoFTools // rule have been set to invalid values // by the used constructor. hp::FEValues hp_fe_values(mapping, fe_collection, - q_coll_dummy, update_quadrature_points); + q_coll_dummy, update_quadrature_points); typename DH::active_cell_iterator cell = dof_handler.begin_active(), endc = dof_handler.end(); std::vector local_dof_indices; for (; cell != endc; ++cell) - { - hp_fe_values.reinit(cell); - const FEValues &fe_values = hp_fe_values.get_present_fe_values(); + // only work on locally relevant cells + if (cell->is_artificial() == false) + { + hp_fe_values.reinit(cell); + const FEValues &fe_values = hp_fe_values.get_present_fe_values(); - local_dof_indices.resize(cell->get_fe().dofs_per_cell); - cell->get_dof_indices(local_dof_indices); + local_dof_indices.resize(cell->get_fe().dofs_per_cell); + cell->get_dof_indices(local_dof_indices); - const std::vector > & points = - fe_values.get_quadrature_points(); - for (unsigned int i = 0; i < cell->get_fe().dofs_per_cell; ++i) - support_points[local_dof_indices[i]] = points[i]; - }; + const std::vector > & points = + fe_values.get_quadrature_points(); + for (unsigned int i = 0; i < cell->get_fe().dofs_per_cell; ++i) + // insert the values into the map + support_points[local_dof_indices[i]] = points[i]; + } } + + + template + void + map_dofs_to_support_points(const hp::MappingCollection & mapping, + const DH &dof_handler, + std::vector > &support_points) + { + // get the data in the form of the map as above + std::map > x_support_points; + map_dofs_to_support_points(mapping, dof_handler, x_support_points); + + // now convert from the map to the linear vector. make sure every + // entry really appeared in the map + for (unsigned int i=0; i > &support_points) { AssertDimension(support_points.size(), dof_handler.n_dofs()); + Assert ((dynamic_cast*> + (&dof_handler.get_tria()) + == + 0), + ExcMessage ("This function can not be used with distributed triangulations." + "See the documentation for more information.")); //Let the internal function do all the work, //just make sure that it gets a MappingCollection const hp::MappingCollection mapping_collection(mapping); - internal::map_dofs_to_support_points_internal( - mapping_collection, - dof_handler, support_points); + internal::map_dofs_to_support_points (mapping_collection, + dof_handler, + support_points); } template void - map_dofs_to_support_points( - const hp::MappingCollection &mapping, - const hp::DoFHandler &dof_handler, - std::vector > &support_points) + map_dofs_to_support_points(const hp::MappingCollection &mapping, + const hp::DoFHandler &dof_handler, + std::vector > &support_points) { AssertDimension(support_points.size(), dof_handler.n_dofs()); + Assert ((dynamic_cast*> + (&dof_handler.get_tria()) + == + 0), + ExcMessage ("This function can not be used with distributed triangulations." + "See the documentation for more information.")); + + //Let the internal function do all the work, + //just make sure that it gets a MappingCollection + internal::map_dofs_to_support_points (mapping, + dof_handler, + support_points); + } + + + template + void + map_dofs_to_support_points (const Mapping &mapping, + const DoFHandler &dof_handler, + std::map > &support_points) + { + support_points.clear(); + + //Let the internal function do all the work, + //just make sure that it gets a MappingCollection + const hp::MappingCollection mapping_collection(mapping); + + internal::map_dofs_to_support_points (mapping_collection, + dof_handler, + support_points); + } + + + template + void + map_dofs_to_support_points(const hp::MappingCollection &mapping, + const hp::DoFHandler &dof_handler, + std::map > &support_points) + { + support_points.clear(); //Let the internal function do all the work, //just make sure that it gets a MappingCollection - internal::map_dofs_to_support_points_internal( - mapping, - dof_handler, - support_points); + internal::map_dofs_to_support_points (mapping, + dof_handler, + support_points); } diff --git a/deal.II/source/dofs/dof_tools.inst.in b/deal.II/source/dofs/dof_tools.inst.in index fe1259fe03..b4c14494df 100644 --- a/deal.II/source/dofs/dof_tools.inst.in +++ b/deal.II/source/dofs/dof_tools.inst.in @@ -812,6 +812,22 @@ DoFTools::map_dofs_to_support_points const hp::DoFHandler&, std::vector >&); + +template +void +DoFTools::map_dofs_to_support_points +(const Mapping&, + const DoFHandler&, + std::map >&); + + +template +void +DoFTools::map_dofs_to_support_points +(const hp::MappingCollection&, + const hp::DoFHandler&, + std::map >&); + #if deal_II_dimension < 3 template @@ -828,6 +844,21 @@ DoFTools::map_dofs_to_support_points const hp::DoFHandler&, std::vector >&); +template +void +DoFTools::map_dofs_to_support_points +(const Mapping&, + const DoFHandler&, + std::map >&); + +template +void +DoFTools::map_dofs_to_support_points +(const hp::MappingCollection&, + const hp::DoFHandler&, + std::map >&); + + template void DoFTools::count_dofs_per_block > ( diff --git a/tests/mpi/map_dofs_to_support_points.cc b/tests/mpi/map_dofs_to_support_points.cc new file mode 100644 index 0000000000..ea51ad6295 --- /dev/null +++ b/tests/mpi/map_dofs_to_support_points.cc @@ -0,0 +1,115 @@ +//--------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2009, 2010, 2011, 2012 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 DoFTools::map_dofs_to_support_points for parallel DoFHandlers + +#include "../tests.h" +#include "coarse_grid_common.h" +#include +#include +#include +#include +#include +#include + +#include + + +template +void test() +{ + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + GridGenerator::hyper_ball(tr); + tr.refine_global (1); + + const FE_Q fe(1); + DoFHandler dofh(tr); + dofh.distribute_dofs (fe); + + std::map > points; + DoFTools::map_dofs_to_support_points (MappingQ1(), + dofh, + points); + if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0) + { + for (typename std::map >::const_iterator + p = points.begin(); + p != points.end(); + ++p) + deallog << p->first << " -> " << p->second + << std::endl; + } + + // the result of the call above is + // supposed to be a map that + // contains exactly the locally + // relevant dofs, so test this + IndexSet relevant_set; + DoFTools::extract_locally_relevant_dofs (dofh, relevant_set); + + for (unsigned int i=0; i(); + deallog.pop(); + + deallog.push("3d"); + test<3>(); + deallog.pop(); + } + else + { + deallog.push("2d"); + test<2>(); + deallog.pop(); + + deallog.push("3d"); + test<3>(); + deallog.pop(); + } + +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + MPI_Finalize(); +#endif +} diff --git a/tests/mpi/map_dofs_to_support_points/ncpu_10/cmp/generic b/tests/mpi/map_dofs_to_support_points/ncpu_10/cmp/generic new file mode 100644 index 0000000000..072424fe61 --- /dev/null +++ b/tests/mpi/map_dofs_to_support_points/ncpu_10/cmp/generic @@ -0,0 +1,69 @@ + +DEAL:0:2d::0 -> -0.707107 -0.707107 +DEAL:0:2d::1 -> 0.00000 -0.707107 +DEAL:0:2d::2 -> -0.500000 -0.500000 +DEAL:0:2d::3 -> 0.00000 -0.500000 +DEAL:0:2d::4 -> 0.707107 -0.707107 +DEAL:0:2d::5 -> 0.500000 -0.500000 +DEAL:0:2d::6 -> -0.292893 -0.292893 +DEAL:0:2d::7 -> 0.00000 -0.292893 +DEAL:0:2d::8 -> 0.292893 -0.292893 +DEAL:0:2d::9 -> -0.707107 0.00000 +DEAL:0:2d::10 -> -0.500000 0.00000 +DEAL:0:2d::11 -> -0.292893 0.00000 +DEAL:0:2d::15 -> 0.707107 0.00000 +DEAL:0:2d::16 -> 0.500000 0.00000 +DEAL:0:2d::19 -> 0.292893 0.00000 +DEAL:0:2d::21 -> 0.00000 0.00000 +DEAL:0:3d::0 -> -0.577350 -0.577350 -0.577350 +DEAL:0:3d::1 -> 0.00000 -0.577350 -0.577350 +DEAL:0:3d::2 -> -0.577350 0.00000 -0.577350 +DEAL:0:3d::3 -> 0.00000 0.00000 -0.577350 +DEAL:0:3d::4 -> -0.394338 -0.394338 -0.394338 +DEAL:0:3d::5 -> 1.73472e-18 -0.394338 -0.394338 +DEAL:0:3d::6 -> -0.394338 0.00000 -0.394338 +DEAL:0:3d::7 -> 2.89121e-19 0.00000 -0.394338 +DEAL:0:3d::8 -> 0.577350 -0.577350 -0.577350 +DEAL:0:3d::9 -> 0.577350 0.00000 -0.577350 +DEAL:0:3d::10 -> 0.394338 -0.394338 -0.394338 +DEAL:0:3d::11 -> 0.394338 0.00000 -0.394338 +DEAL:0:3d::12 -> -0.577350 0.577350 -0.577350 +DEAL:0:3d::13 -> 0.00000 0.577350 -0.577350 +DEAL:0:3d::14 -> -0.394338 0.394338 -0.394338 +DEAL:0:3d::15 -> 1.73472e-18 0.394338 -0.394338 +DEAL:0:3d::16 -> 0.577350 0.577350 -0.577350 +DEAL:0:3d::17 -> 0.394338 0.394338 -0.394338 +DEAL:0:3d::18 -> -0.211325 -0.211325 -0.211325 +DEAL:0:3d::19 -> 0.00000 -0.211325 -0.211325 +DEAL:0:3d::20 -> -0.211325 0.00000 -0.211325 +DEAL:0:3d::21 -> 0.00000 0.00000 -0.211325 +DEAL:0:3d::22 -> 0.211325 -0.211325 -0.211325 +DEAL:0:3d::23 -> 0.211325 0.00000 -0.211325 +DEAL:0:3d::24 -> -0.211325 0.211325 -0.211325 +DEAL:0:3d::25 -> 0.00000 0.211325 -0.211325 +DEAL:0:3d::26 -> 0.211325 0.211325 -0.211325 +DEAL:0:3d::27 -> 0.577350 -0.577350 0.00000 +DEAL:0:3d::28 -> 0.577350 0.00000 0.00000 +DEAL:0:3d::29 -> 0.394338 -0.394338 1.73472e-18 +DEAL:0:3d::30 -> 0.394338 0.00000 6.93889e-18 +DEAL:0:3d::31 -> 0.577350 0.577350 0.00000 +DEAL:0:3d::32 -> 0.394338 0.394338 1.73472e-18 +DEAL:0:3d::33 -> 0.211325 -0.211325 0.00000 +DEAL:0:3d::34 -> 0.211325 0.00000 0.00000 +DEAL:0:3d::35 -> 0.211325 0.211325 0.00000 +DEAL:0:3d::45 -> -0.577350 -0.577350 0.00000 +DEAL:0:3d::46 -> -0.394338 -0.394338 0.00000 +DEAL:0:3d::47 -> -0.577350 0.00000 0.00000 +DEAL:0:3d::48 -> -0.394338 1.44560e-19 -6.93889e-18 +DEAL:0:3d::49 -> -0.211325 -0.211325 0.00000 +DEAL:0:3d::50 -> -0.211325 0.00000 0.00000 +DEAL:0:3d::51 -> -0.577350 0.577350 0.00000 +DEAL:0:3d::52 -> -0.394338 0.394338 0.00000 +DEAL:0:3d::53 -> -0.211325 0.211325 0.00000 +DEAL:0:3d::63 -> 0.00000 -0.577350 0.00000 +DEAL:0:3d::64 -> 1.44560e-19 -0.394338 6.93889e-18 +DEAL:0:3d::65 -> 0.00000 -0.211325 0.00000 +DEAL:0:3d::69 -> 0.00000 0.577350 0.00000 +DEAL:0:3d::70 -> 1.44560e-19 0.394338 0.00000 +DEAL:0:3d::71 -> 0.00000 0.211325 0.00000 +DEAL:0:3d::75 -> 0.00000 0.00000 0.00000 diff --git a/tests/mpi/map_dofs_to_support_points/ncpu_4/cmp/generic b/tests/mpi/map_dofs_to_support_points/ncpu_4/cmp/generic new file mode 100644 index 0000000000..bde99dec8b --- /dev/null +++ b/tests/mpi/map_dofs_to_support_points/ncpu_4/cmp/generic @@ -0,0 +1,87 @@ + +DEAL:0:2d::0 -> -0.707107 -0.707107 +DEAL:0:2d::1 -> 0.00000 -0.707107 +DEAL:0:2d::2 -> -0.500000 -0.500000 +DEAL:0:2d::3 -> 0.00000 -0.500000 +DEAL:0:2d::4 -> 0.707107 -0.707107 +DEAL:0:2d::5 -> 0.500000 -0.500000 +DEAL:0:2d::6 -> -0.292893 -0.292893 +DEAL:0:2d::7 -> 0.00000 -0.292893 +DEAL:0:2d::8 -> 0.292893 -0.292893 +DEAL:0:2d::9 -> -0.707107 0.00000 +DEAL:0:2d::10 -> -0.500000 0.00000 +DEAL:0:2d::11 -> -0.292893 0.00000 +DEAL:0:2d::15 -> 0.707107 0.00000 +DEAL:0:2d::16 -> 0.500000 0.00000 +DEAL:0:2d::19 -> 0.292893 0.00000 +DEAL:0:2d::21 -> 0.00000 0.00000 +DEAL:0:3d::0 -> -0.577350 -0.577350 -0.577350 +DEAL:0:3d::1 -> 0.00000 -0.577350 -0.577350 +DEAL:0:3d::2 -> -0.577350 0.00000 -0.577350 +DEAL:0:3d::3 -> 0.00000 0.00000 -0.577350 +DEAL:0:3d::4 -> -0.394338 -0.394338 -0.394338 +DEAL:0:3d::5 -> 1.73472e-18 -0.394338 -0.394338 +DEAL:0:3d::6 -> -0.394338 0.00000 -0.394338 +DEAL:0:3d::7 -> 2.89121e-19 0.00000 -0.394338 +DEAL:0:3d::8 -> 0.577350 -0.577350 -0.577350 +DEAL:0:3d::9 -> 0.577350 0.00000 -0.577350 +DEAL:0:3d::10 -> 0.394338 -0.394338 -0.394338 +DEAL:0:3d::11 -> 0.394338 0.00000 -0.394338 +DEAL:0:3d::12 -> -0.577350 0.577350 -0.577350 +DEAL:0:3d::13 -> 0.00000 0.577350 -0.577350 +DEAL:0:3d::14 -> -0.394338 0.394338 -0.394338 +DEAL:0:3d::15 -> 1.73472e-18 0.394338 -0.394338 +DEAL:0:3d::16 -> 0.577350 0.577350 -0.577350 +DEAL:0:3d::17 -> 0.394338 0.394338 -0.394338 +DEAL:0:3d::18 -> -0.211325 -0.211325 -0.211325 +DEAL:0:3d::19 -> 0.00000 -0.211325 -0.211325 +DEAL:0:3d::20 -> -0.211325 0.00000 -0.211325 +DEAL:0:3d::21 -> 0.00000 0.00000 -0.211325 +DEAL:0:3d::22 -> 0.211325 -0.211325 -0.211325 +DEAL:0:3d::23 -> 0.211325 0.00000 -0.211325 +DEAL:0:3d::24 -> -0.211325 0.211325 -0.211325 +DEAL:0:3d::25 -> 0.00000 0.211325 -0.211325 +DEAL:0:3d::26 -> 0.211325 0.211325 -0.211325 +DEAL:0:3d::27 -> 0.577350 -0.577350 0.00000 +DEAL:0:3d::28 -> 0.577350 0.00000 0.00000 +DEAL:0:3d::29 -> 0.394338 -0.394338 1.73472e-18 +DEAL:0:3d::30 -> 0.394338 0.00000 6.93889e-18 +DEAL:0:3d::31 -> 0.577350 0.577350 0.00000 +DEAL:0:3d::32 -> 0.394338 0.394338 1.73472e-18 +DEAL:0:3d::33 -> 0.211325 -0.211325 0.00000 +DEAL:0:3d::34 -> 0.211325 0.00000 0.00000 +DEAL:0:3d::35 -> 0.211325 0.211325 0.00000 +DEAL:0:3d::36 -> 0.577350 -0.577350 0.577350 +DEAL:0:3d::37 -> 0.577350 0.00000 0.577350 +DEAL:0:3d::38 -> 0.394338 -0.394338 0.394338 +DEAL:0:3d::39 -> 0.394338 0.00000 0.394338 +DEAL:0:3d::40 -> 0.577350 0.577350 0.577350 +DEAL:0:3d::41 -> 0.394338 0.394338 0.394338 +DEAL:0:3d::42 -> 0.211325 -0.211325 0.211325 +DEAL:0:3d::43 -> 0.211325 0.00000 0.211325 +DEAL:0:3d::44 -> 0.211325 0.211325 0.211325 +DEAL:0:3d::45 -> -0.577350 -0.577350 0.00000 +DEAL:0:3d::46 -> -0.394338 -0.394338 0.00000 +DEAL:0:3d::47 -> -0.577350 0.00000 0.00000 +DEAL:0:3d::48 -> -0.394338 1.44560e-19 -6.93889e-18 +DEAL:0:3d::49 -> -0.211325 -0.211325 0.00000 +DEAL:0:3d::50 -> -0.211325 0.00000 0.00000 +DEAL:0:3d::51 -> -0.577350 0.577350 0.00000 +DEAL:0:3d::52 -> -0.394338 0.394338 0.00000 +DEAL:0:3d::53 -> -0.211325 0.211325 0.00000 +DEAL:0:3d::63 -> 0.00000 -0.577350 0.00000 +DEAL:0:3d::64 -> 1.44560e-19 -0.394338 6.93889e-18 +DEAL:0:3d::65 -> 0.00000 -0.211325 0.00000 +DEAL:0:3d::66 -> 0.00000 -0.577350 0.577350 +DEAL:0:3d::67 -> 0.00000 -0.394338 0.394338 +DEAL:0:3d::68 -> 0.00000 -0.211325 0.211325 +DEAL:0:3d::69 -> 0.00000 0.577350 0.00000 +DEAL:0:3d::70 -> 1.44560e-19 0.394338 0.00000 +DEAL:0:3d::71 -> 0.00000 0.211325 0.00000 +DEAL:0:3d::72 -> 0.00000 0.577350 0.577350 +DEAL:0:3d::73 -> 0.00000 0.394338 0.394338 +DEAL:0:3d::74 -> 0.00000 0.211325 0.211325 +DEAL:0:3d::75 -> 0.00000 0.00000 0.00000 +DEAL:0:3d::76 -> 0.00000 0.00000 0.211325 +DEAL:0:3d::77 -> 0.00000 0.00000 0.577350 +DEAL:0:3d::78 -> 0.00000 6.93889e-18 0.394338 -- 2.39.5