<h3>Specific improvements</h3>
<ol>
+<li> New: There is now a second DoFTools::map_dofs_to_support_points
+function that also works for parallel::distributed::Triangulation
+triangulations.
+<br>
+(Wolfgang Bangerth, 2012/04/26)
+
<li> New: There is now a second DoFTools::extract_boundary_dofs
function that also works for parallel::distributed::Triangulation
triangulations.
* 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 <int dim, int spacedim>
void
std::vector<Point<spacedim> > &support_points);
/**
- * Same as above for the hp case.
+ * Same as the previous function but for the hp case.
*/
-
template <int dim, int spacedim>
void
map_dofs_to_support_points (const dealii::hp::MappingCollection<dim,spacedim> &mapping,
const hp::DoFHandler<dim,spacedim> &dof_handler,
std::vector<Point<spacedim> > &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 <int dim, int spacedim>
+ void
+ map_dofs_to_support_points (const Mapping<dim,spacedim> &mapping,
+ const DoFHandler<dim,spacedim> &dof_handler,
+ std::map<unsigned int, Point<spacedim> > &support_points);
+
+ /**
+ * Same as the previous function but for the hp case.
+ */
+ template <int dim, int spacedim>
+ void
+ map_dofs_to_support_points (const dealii::hp::MappingCollection<dim,spacedim> &mapping,
+ const hp::DoFHandler<dim,spacedim> &dof_handler,
+ std::map<unsigned int, Point<spacedim> > &support_points);
+
/**
* This is the opposite function
* to the one above. It generates
{
namespace
{
- template<class DH>
- void map_dofs_to_support_points_internal(
- const hp::MappingCollection<DH::dimension, DH::space_dimension> & mapping,
- const DH &dof_handler,
- std::vector<Point<DH::space_dimension> > &support_points)
+ template <class DH>
+ void
+ map_dofs_to_support_points(const hp::MappingCollection<DH::dimension, DH::space_dimension> & mapping,
+ const DH &dof_handler,
+ std::map<unsigned int,Point<DH::space_dimension> > &support_points)
{
- AssertDimension(support_points.size(), dof_handler.n_dofs());
-
const unsigned int dim = DH::dimension;
const unsigned int spacedim = DH::space_dimension;
// check whether every fe in the collection
// has support points
Assert(fe_collection[fe_index].has_support_points(),
- typename FiniteElement<dim>::ExcFEHasNoSupportPoints());
+ typename FiniteElement<dim>::ExcFEHasNoSupportPoints());
q_coll_dummy.push_back(
Quadrature<dim> (
fe_collection[fe_index].get_unit_support_points()));
// rule have been set to invalid values
// by the used constructor.
hp::FEValues<dim, spacedim> 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<unsigned int> local_dof_indices;
for (; cell != endc; ++cell)
- {
- hp_fe_values.reinit(cell);
- const FEValues<dim, spacedim> &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<dim, spacedim> &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<Point<spacedim> > & 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<Point<spacedim> > & 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 <class DH>
+ void
+ map_dofs_to_support_points(const hp::MappingCollection<DH::dimension, DH::space_dimension> & mapping,
+ const DH &dof_handler,
+ std::vector<Point<DH::space_dimension> > &support_points)
+ {
+ // get the data in the form of the map as above
+ std::map<unsigned int,Point<DH::space_dimension> > 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<dof_handler.n_dofs(); ++i)
+ {
+ Assert (x_support_points.find(i) != x_support_points.end(),
+ ExcInternalError());
+ support_points[i] = x_support_points[i];
+ }
+ }
}
}
std::vector<Point<spacedim> > &support_points)
{
AssertDimension(support_points.size(), dof_handler.n_dofs());
+ Assert ((dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim>*>
+ (&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<dim, spacedim> 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<int dim, int spacedim>
void
- map_dofs_to_support_points(
- const hp::MappingCollection<dim, spacedim> &mapping,
- const hp::DoFHandler<dim, spacedim> &dof_handler,
- std::vector<Point<spacedim> > &support_points)
+ map_dofs_to_support_points(const hp::MappingCollection<dim, spacedim> &mapping,
+ const hp::DoFHandler<dim, spacedim> &dof_handler,
+ std::vector<Point<spacedim> > &support_points)
{
AssertDimension(support_points.size(), dof_handler.n_dofs());
+ Assert ((dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim>*>
+ (&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 <int dim, int spacedim>
+ void
+ map_dofs_to_support_points (const Mapping<dim,spacedim> &mapping,
+ const DoFHandler<dim,spacedim> &dof_handler,
+ std::map<unsigned int, Point<spacedim> > &support_points)
+ {
+ support_points.clear();
+
+ //Let the internal function do all the work,
+ //just make sure that it gets a MappingCollection
+ const hp::MappingCollection<dim, spacedim> mapping_collection(mapping);
+
+ internal::map_dofs_to_support_points (mapping_collection,
+ dof_handler,
+ support_points);
+ }
+
+
+ template<int dim, int spacedim>
+ void
+ map_dofs_to_support_points(const hp::MappingCollection<dim, spacedim> &mapping,
+ const hp::DoFHandler<dim, spacedim> &dof_handler,
+ std::map<unsigned int, Point<spacedim> > &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);
}
const hp::DoFHandler<deal_II_dimension>&,
std::vector<Point<deal_II_dimension> >&);
+
+template
+void
+DoFTools::map_dofs_to_support_points<deal_II_dimension>
+(const Mapping<deal_II_dimension,deal_II_dimension>&,
+ const DoFHandler<deal_II_dimension>&,
+ std::map<unsigned int, Point<deal_II_dimension> >&);
+
+
+template
+void
+DoFTools::map_dofs_to_support_points<deal_II_dimension>
+(const hp::MappingCollection<deal_II_dimension,deal_II_dimension>&,
+ const hp::DoFHandler<deal_II_dimension>&,
+ std::map<unsigned int, Point<deal_II_dimension> >&);
+
#if deal_II_dimension < 3
template
const hp::DoFHandler<deal_II_dimension, deal_II_dimension+1>&,
std::vector<Point<deal_II_dimension+1> >&);
+template
+void
+DoFTools::map_dofs_to_support_points<deal_II_dimension,deal_II_dimension+1>
+(const Mapping<deal_II_dimension,deal_II_dimension+1>&,
+ const DoFHandler<deal_II_dimension, deal_II_dimension+1>&,
+ std::map<unsigned int, Point<deal_II_dimension+1> >&);
+
+template
+void
+DoFTools::map_dofs_to_support_points<deal_II_dimension,deal_II_dimension+1>
+(const hp::MappingCollection<deal_II_dimension,deal_II_dimension+1>&,
+ const hp::DoFHandler<deal_II_dimension, deal_II_dimension+1>&,
+ std::map<unsigned int, Point<deal_II_dimension+1> >&);
+
+
template
void
DoFTools::count_dofs_per_block<DoFHandler<deal_II_dimension,deal_II_dimension+1> > (
--- /dev/null
+//---------------------------------------------------------------------------
+// $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 <deal.II/base/logstream.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/fe/fe_q.h>
+
+#include <fstream>
+
+
+template<int dim>
+void test()
+{
+ parallel::distributed::Triangulation<dim> tr(MPI_COMM_WORLD);
+
+ GridGenerator::hyper_ball(tr);
+ tr.refine_global (1);
+
+ const FE_Q<dim> fe(1);
+ DoFHandler<dim> dofh(tr);
+ dofh.distribute_dofs (fe);
+
+ std::map<unsigned int, Point<dim> > points;
+ DoFTools::map_dofs_to_support_points (MappingQ1<dim>(),
+ dofh,
+ points);
+ if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0)
+ {
+ for (typename std::map<unsigned int, Point<dim> >::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<dofh.n_dofs(); ++i)
+ if (relevant_set.is_element(i))
+ Assert (points.find(i) != points.end(),
+ ExcInternalError())
+ else
+ Assert (points.find(i) == points.end(),
+ ExcInternalError());
+}
+
+
+int main(int argc, char *argv[])
+{
+#ifdef DEAL_II_COMPILER_SUPPORTS_MPI
+ MPI_Init (&argc,&argv);
+#else
+ (void)argc;
+ (void)argv;
+#endif
+
+ unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+
+
+ deallog.push(Utilities::int_to_string(myid));
+
+ if (myid == 0)
+ {
+ std::ofstream logfile(output_file_for_mpi("map_dofs_to_support_points").c_str());
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+ deallog.threshold_double(1.e-10);
+
+ deallog.push("2d");
+ test<2>();
+ 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
+}
--- /dev/null
+
+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
--- /dev/null
+
+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