diff --git a/src/io/mesh_io/mesh_io_msh.cc b/src/io/mesh_io/mesh_io_msh.cc
index 80941e3ea..40de37cc6 100644
--- a/src/io/mesh_io/mesh_io_msh.cc
+++ b/src/io/mesh_io/mesh_io_msh.cc
@@ -1,1120 +1,731 @@
/**
* @file mesh_io_msh.cc
*
* @author Dana Christen <dana.christen@gmail.com>
* @author Mauro Corrado <mauro.corrado@epfl.ch>
* @author David Simon Kammer <david.kammer@epfl.ch>
* @author Nicolas Richart <nicolas.richart@epfl.ch>
*
* @date creation: Fri Jun 18 2010
* @date last modification: Tue Feb 20 2018
*
* @brief Read/Write for MSH files generated by gmsh
*
* @section LICENSE
*
* Copyright (©) 2010-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* Akantu is free software: you can 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 3 of the License, or (at your option) any
* later version.
*
* Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
* A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Akantu. If not, see <http://www.gnu.org/licenses/>.
*
*/
/* -----------------------------------------------------------------------------
Version (Legacy) 1.0
$NOD
number-of-nodes
node-number x-coord y-coord z-coord
...
$ENDNOD
$ELM
number-of-elements
elm-number elm-type reg-phys reg-elem number-of-nodes node-number-list
...
$ENDELM
-----------------------------------------------------------------------------
Version 2.1
$MeshFormat
version-number file-type data-size
$EndMeshFormat
$Nodes
number-of-nodes
node-number x-coord y-coord z-coord
...
$EndNodes
$Elements
number-of-elements
elm-number elm-type number-of-tags < tag > ... node-number-list
...
$EndElements
$PhysicalNames
number-of-names
physical-dimension physical-number "physical-name"
...
$EndPhysicalNames
$NodeData
number-of-string-tags
< "string-tag" >
...
number-of-real-tags
< real-tag >
...
number-of-integer-tags
< integer-tag >
...
node-number value ...
...
$EndNodeData
$ElementData
number-of-string-tags
< "string-tag" >
...
number-of-real-tags
< real-tag >
...
number-of-integer-tags
< integer-tag >
...
elm-number value ...
...
$EndElementData
$ElementNodeData
number-of-string-tags
< "string-tag" >
...
number-of-real-tags
< real-tag >
...
number-of-integer-tags
< integer-tag >
...
elm-number number-of-nodes-per-element value ...
...
$ElementEndNodeData
-----------------------------------------------------------------------------
elem-type
1: 2-node line.
2: 3-node triangle.
3: 4-node quadrangle.
4: 4-node tetrahedron.
5: 8-node hexahedron.
6: 6-node prism.
7: 5-node pyramid.
8: 3-node second order line
9: 6-node second order triangle
10: 9-node second order quadrangle
11: 10-node second order tetrahedron
12: 27-node second order hexahedron
13: 18-node second order prism
14: 14-node second order pyramid
15: 1-node point.
16: 8-node second order quadrangle
17: 20-node second order hexahedron
18: 15-node second order prism
19: 13-node second order pyramid
20: 9-node third order incomplete triangle
21: 10-node third order triangle
22: 12-node fourth order incomplete triangle
23: 15-node fourth order triangle
24: 15-node fifth order incomplete triangle
25: 21-node fifth order complete triangle
26: 4-node third order edge
27: 5-node fourth order edge
28: 6-node fifth order edge
29: 20-node third order tetrahedron
30: 35-node fourth order tetrahedron
31: 56-node fifth order tetrahedron
-------------------------------------------------------------------------- */
/* -------------------------------------------------------------------------- */
#include <fstream>
/* -------------------------------------------------------------------------- */
#include "mesh_io.hh"
#include "mesh_utils.hh"
/* -------------------------------------------------------------------------- */
-/* -------------------------------------------------------------------------- */
-// The boost spirit is a work on the way it does not compile so I kept the
-// current code. The current code does not handle files generated on Windows
-// <CRLF>
-
-// #include <boost/config/warning_disable.hpp>
-// #include <boost/spirit/include/qi.hpp>
-// #include <boost/spirit/include/phoenix_core.hpp>
-// #include <boost/spirit/include/phoenix_fusion.hpp>
-// #include <boost/spirit/include/phoenix_object.hpp>
-// #include <boost/spirit/include/phoenix_container.hpp>
-// #include <boost/spirit/include/phoenix_operator.hpp>
-// #include <boost/spirit/include/phoenix_bind.hpp>
-// #include <boost/spirit/include/phoenix_stl.hpp>
-/* -------------------------------------------------------------------------- */
-
namespace akantu {
/* -------------------------------------------------------------------------- */
/* Methods Implentations */
/* -------------------------------------------------------------------------- */
MeshIOMSH::MeshIOMSH() {
canReadSurface = true;
canReadExtendedData = true;
_msh_nodes_per_elem[_msh_not_defined] = 0;
_msh_nodes_per_elem[_msh_segment_2] = 2;
_msh_nodes_per_elem[_msh_triangle_3] = 3;
_msh_nodes_per_elem[_msh_quadrangle_4] = 4;
_msh_nodes_per_elem[_msh_tetrahedron_4] = 4;
_msh_nodes_per_elem[_msh_hexahedron_8] = 8;
_msh_nodes_per_elem[_msh_prism_1] = 6;
_msh_nodes_per_elem[_msh_pyramid_1] = 1;
_msh_nodes_per_elem[_msh_segment_3] = 3;
_msh_nodes_per_elem[_msh_triangle_6] = 6;
_msh_nodes_per_elem[_msh_quadrangle_9] = 9;
_msh_nodes_per_elem[_msh_tetrahedron_10] = 10;
_msh_nodes_per_elem[_msh_hexahedron_27] = 27;
_msh_nodes_per_elem[_msh_hexahedron_20] = 20;
_msh_nodes_per_elem[_msh_prism_18] = 18;
_msh_nodes_per_elem[_msh_prism_15] = 15;
_msh_nodes_per_elem[_msh_pyramid_14] = 14;
_msh_nodes_per_elem[_msh_point] = 1;
_msh_nodes_per_elem[_msh_quadrangle_8] = 8;
_msh_to_akantu_element_types[_msh_not_defined] = _not_defined;
_msh_to_akantu_element_types[_msh_segment_2] = _segment_2;
_msh_to_akantu_element_types[_msh_triangle_3] = _triangle_3;
_msh_to_akantu_element_types[_msh_quadrangle_4] = _quadrangle_4;
_msh_to_akantu_element_types[_msh_tetrahedron_4] = _tetrahedron_4;
_msh_to_akantu_element_types[_msh_hexahedron_8] = _hexahedron_8;
_msh_to_akantu_element_types[_msh_prism_1] = _pentahedron_6;
_msh_to_akantu_element_types[_msh_pyramid_1] = _not_defined;
_msh_to_akantu_element_types[_msh_segment_3] = _segment_3;
_msh_to_akantu_element_types[_msh_triangle_6] = _triangle_6;
_msh_to_akantu_element_types[_msh_quadrangle_9] = _not_defined;
_msh_to_akantu_element_types[_msh_tetrahedron_10] = _tetrahedron_10;
_msh_to_akantu_element_types[_msh_hexahedron_27] = _not_defined;
_msh_to_akantu_element_types[_msh_hexahedron_20] = _hexahedron_20;
_msh_to_akantu_element_types[_msh_prism_18] = _not_defined;
_msh_to_akantu_element_types[_msh_prism_15] = _pentahedron_15;
_msh_to_akantu_element_types[_msh_pyramid_14] = _not_defined;
_msh_to_akantu_element_types[_msh_point] = _point_1;
_msh_to_akantu_element_types[_msh_quadrangle_8] = _quadrangle_8;
_akantu_to_msh_element_types[_not_defined] = _msh_not_defined;
_akantu_to_msh_element_types[_segment_2] = _msh_segment_2;
_akantu_to_msh_element_types[_segment_3] = _msh_segment_3;
_akantu_to_msh_element_types[_triangle_3] = _msh_triangle_3;
_akantu_to_msh_element_types[_triangle_6] = _msh_triangle_6;
_akantu_to_msh_element_types[_tetrahedron_4] = _msh_tetrahedron_4;
_akantu_to_msh_element_types[_tetrahedron_10] = _msh_tetrahedron_10;
_akantu_to_msh_element_types[_quadrangle_4] = _msh_quadrangle_4;
_akantu_to_msh_element_types[_quadrangle_8] = _msh_quadrangle_8;
_akantu_to_msh_element_types[_hexahedron_8] = _msh_hexahedron_8;
_akantu_to_msh_element_types[_hexahedron_20] = _msh_hexahedron_20;
_akantu_to_msh_element_types[_pentahedron_6] = _msh_prism_1;
_akantu_to_msh_element_types[_pentahedron_15] = _msh_prism_15;
_akantu_to_msh_element_types[_point_1] = _msh_point;
#if defined(AKANTU_STRUCTURAL_MECHANICS)
_akantu_to_msh_element_types[_bernoulli_beam_2] = _msh_segment_2;
_akantu_to_msh_element_types[_bernoulli_beam_3] = _msh_segment_2;
_akantu_to_msh_element_types[_discrete_kirchhoff_triangle_18] =
_msh_triangle_3;
#endif
std::map<ElementType, MSHElementType>::iterator it;
for (it = _akantu_to_msh_element_types.begin();
it != _akantu_to_msh_element_types.end(); ++it) {
UInt nb_nodes = _msh_nodes_per_elem[it->second];
std::vector<UInt> tmp(nb_nodes);
for (UInt i = 0; i < nb_nodes; ++i) {
tmp[i] = i;
}
switch (it->first) {
case _tetrahedron_10:
tmp[8] = 9;
tmp[9] = 8;
break;
case _pentahedron_6:
tmp[0] = 2;
tmp[1] = 0;
tmp[2] = 1;
tmp[3] = 5;
tmp[4] = 3;
tmp[5] = 4;
break;
case _pentahedron_15:
tmp[0] = 2;
tmp[1] = 0;
tmp[2] = 1;
tmp[3] = 5;
tmp[4] = 3;
tmp[5] = 4;
tmp[6] = 8;
tmp[8] = 11;
tmp[9] = 6;
tmp[10] = 9;
tmp[11] = 10;
tmp[12] = 14;
tmp[14] = 12;
break;
case _hexahedron_20:
tmp[9] = 11;
tmp[10] = 12;
tmp[11] = 9;
tmp[12] = 13;
tmp[13] = 10;
tmp[17] = 19;
tmp[18] = 17;
tmp[19] = 18;
break;
default:
// nothing to change
break;
}
_read_order[it->first] = tmp;
}
}
/* -------------------------------------------------------------------------- */
MeshIOMSH::~MeshIOMSH() = default;
/* -------------------------------------------------------------------------- */
-/* Spirit stuff */
-/* -------------------------------------------------------------------------- */
-// namespace _parser {
-// namespace spirit = ::boost::spirit;
-// namespace qi = ::boost::spirit::qi;
-// namespace ascii = ::boost::spirit::ascii;
-// namespace lbs = ::boost::spirit::qi::labels;
-// namespace phx = ::boost::phoenix;
-
-// /* ------------------------------------------------------------------------
-// */
-// /* Lazy functors */
-// /* ------------------------------------------------------------------------
-// */
-// struct _Element {
-// int index;
-// std::vector<int> tags;
-// std::vector<int> connectivity;
-// ElementType type;
-// };
-
-// /* ------------------------------------------------------------------------
-// */
-// struct lazy_get_nb_nodes_ {
-// template <class elem_type> struct result { typedef int type; };
-// template <class elem_type> bool operator()(elem_type et) {
-// return MeshIOMSH::_msh_nodes_per_elem[et];
-// }
-// };
-
-// /* ------------------------------------------------------------------------
-// */
-// struct lazy_element_ {
-// template <class id_t, class tags_t, class elem_type, class conn_t>
-// struct result {
-// typedef _Element type;
-// };
-// template <class id_t, class tags_t, class elem_type, class conn_t>
-// _Element operator()(id_t id, const elem_type & et, const tags_t & t,
-// const conn_t & c) {
-// _Element tmp_el;
-// tmp_el.index = id;
-// tmp_el.tags = t;
-// tmp_el.connectivity = c;
-// tmp_el.type = et;
-// return tmp_el;
-// }
-// };
-
-// /* ------------------------------------------------------------------------
-// */
-// struct lazy_check_value_ {
-// template <class T> struct result { typedef void type; };
-// template <class T> void operator()(T v1, T v2) {
-// if (v1 != v2) {
-// AKANTU_EXCEPTION("The msh parser expected a "
-// << v2 << " in the header bug got a " << v1);
-// }
-// }
-// };
-
-// /* ------------------------------------------------------------------------
-// */
-// struct lazy_node_read_ {
-// template <class Mesh, class ID, class V, class size, class Map>
-// struct result {
-// typedef bool type;
-// };
-// template <class Mesh, class ID, class V, class size, class Map>
-// bool operator()(Mesh & mesh, const ID & id, const V & pos, size max,
-// Map & nodes_mapping) const {
-// Vector<Real> tmp_pos(mesh.getSpatialDimension());
-// UInt i = 0;
-// for (typename V::const_iterator it = pos.begin();
-// it != pos.end() || i < mesh.getSpatialDimension(); ++it)
-// tmp_pos[i++] = *it;
-
-// nodes_mapping[id] = mesh.getNbNodes();
-// mesh.getNodes().push_back(tmp_pos);
-// return (mesh.getNbNodes() < UInt(max));
-// }
-// };
-
-// /* ------------------------------------------------------------------------
-// */
-// struct lazy_element_read_ {
-// template <class Mesh, class EL, class size, class NodeMap, class ElemMap>
-// struct result {
-// typedef bool type;
-// };
-
-// template <class Mesh, class EL, class size, class NodeMap, class ElemMap>
-// bool operator()(Mesh & mesh, const EL & element, size max,
-// const NodeMap & nodes_mapping,
-// ElemMap & elements_mapping) const {
-// Vector<UInt> tmp_conn(Mesh::getNbNodesPerElement(element.type));
-
-// AKANTU_DEBUG_ASSERT(element.connectivity.size() == tmp_conn.size(),
-// "The element "
-// << element.index
-// << "in the MSH file has too many nodes.");
-
-// mesh.addConnectivityType(element.type);
-// Array<UInt> & connectivity = mesh.getConnectivity(element.type);
-
-// UInt i = 0;
-// for (std::vector<int>::const_iterator it =
-// element.connectivity.begin();
-// it != element.connectivity.end(); ++it) {
-// typename NodeMap::const_iterator nit = nodes_mapping.find(*it);
-// AKANTU_DEBUG_ASSERT(nit != nodes_mapping.end(),
-// "There is an unknown node in the connectivity.");
-// tmp_conn[i++] = nit->second;
-// }
-
-// ::akantu::Element el(element.type, connectivity.size());
-// elements_mapping[element.index] = el;
-
-// connectivity.push_back(tmp_conn);
-
-// for (UInt i = 0; i < element.tags.size(); ++i) {
-// std::stringstream tag_sstr;
-// tag_sstr << "tag_" << i;
-// Array<UInt> * data =
-// mesh.template getDataPointer<UInt>(tag_sstr.str(), element.type,
-// _not_ghost);
-// data->push_back(element.tags[i]);
-// }
-
-// return (mesh.getNbElement() < UInt(max));
-// }
-// };
-
-// /* ------------------------------------------------------------------------
-// */
-// template <class Iterator, typename Skipper = ascii::space_type>
-// struct MshMeshGrammar : qi::grammar<Iterator, void(), Skipper> {
-// public:
-// MshMeshGrammar(Mesh & mesh)
-// : MshMeshGrammar::base_type(start, "msh_mesh_reader"), mesh(mesh) {
-// phx::function<lazy_element_> lazy_element;
-// phx::function<lazy_get_nb_nodes_> lazy_get_nb_nodes;
-// phx::function<lazy_check_value_> lazy_check_value;
-// phx::function<lazy_node_read_> lazy_node_read;
-// phx::function<lazy_element_read_> lazy_element_read;
-
-// clang-format off
-
-// start
-// = *( known_section | unknown_section
-// )
-// ;
-
-// known_section
-// = qi::char_("$")
-// >> sections [ lbs::_a = lbs::_1 ]
-// >> qi::lazy(*lbs::_a)
-// >> qi::lit("$End")
-// //>> qi::lit(phx::ref(lbs::_a))
-// ;
-
-// unknown_section
-// = qi::char_("$")
-// >> qi::char_("") [ lbs::_a = lbs::_1 ]
-// >> ignore_section
-// >> qi::lit("$End")
-// >> qi::lit(phx::val(lbs::_a))
-// ;
-
-// mesh_format // version followed by 0 or 1 for ascii or binary
-// = version >> (
-// ( qi::char_("0")
-// >> qi::int_ [ lazy_check_value(lbs::_1, 8) ]
-// )
-// | ( qi::char_("1")
-// >> qi::int_ [ lazy_check_value(lbs::_1, 8) ]
-// >> qi::dword [ lazy_check_value(lbs::_1, 1) ]
-// )
-// )
-// ;
-
-// nodes
-// = nodes_
-// ;
-
-// nodes_
-// = qi::int_ [ lbs::_a = lbs::_1 ]
-// > *(
-// ( qi::int_ >> position ) [ lazy_node_read(phx::ref(mesh),
-// lbs::_1,
-// phx::cref(lbs::_2),
-// lbs::_a,
-// phx::ref(this->msh_nodes_to_akantu)) ]
-// )
-// ;
-
-// element
-// = elements_
-// ;
-
-// elements_
-// = qi::int_ [ lbs::_a = lbs::_1 ]
-// > qi::repeat(phx::ref(lbs::_a))[ element [ lazy_element_read(phx::ref(mesh),
-// lbs::_1,
-// phx::cref(lbs::_a),
-// phx::cref(this->msh_nodes_to_akantu),
-// phx::ref(this->msh_elemt_to_akantu)) ]]
-// ;
-
-// ignore_section
-// = *(qi::char_ - qi::char_('$'))
-// ;
-
-// interpolation_scheme = ignore_section;
-// element_data = ignore_section;
-// node_data = ignore_section;
-
-// version
-// = qi::int_ [ phx::push_back(lbs::_val, lbs::_1) ]
-// >> *( qi::char_(".") >> qi::int_ [ phx::push_back(lbs::_val, lbs::_1) ] )
-// ;
-
-// position
-// = real [ phx::push_back(lbs::_val, lbs::_1) ]
-// > real [ phx::push_back(lbs::_val, lbs::_1) ]
-// > real [ phx::push_back(lbs::_val, lbs::_1) ]
-// ;
-
-// tags
-// = qi::int_ [ lbs::_a = lbs::_1 ]
-// > qi::repeat(phx::val(lbs::_a))[ qi::int_ [ phx::push_back(lbs::_val,
-// lbs::_1) ] ]
-// ;
-
-// element
-// = ( qi::int_ [ lbs::_a = lbs::_1 ]
-// > msh_element_type [ lbs::_b = lazy_get_nb_nodes(lbs::_1) ]
-// > tags [ lbs::_c = lbs::_1 ]
-// > connectivity(lbs::_a) [ lbs::_d = lbs::_1 ]
-// ) [ lbs::_val = lazy_element(lbs::_a,
-// phx::cref(lbs::_b),
-// phx::cref(lbs::_c),
-// phx::cref(lbs::_d)) ]
-// ;
-// connectivity
-// = qi::repeat(lbs::_r1)[ qi::int_ [ phx::push_back(lbs::_val,
-// lbs::_1) ] ]
-// ;
-
-// sections.add
-// ("MeshFormat", &mesh_format)
-// ("Nodes", &nodes)
-// ("Elements", &elements)
-// ("PhysicalNames", &physical_names)
-// ("InterpolationScheme", &interpolation_scheme)
-// ("ElementData", &element_data)
-// ("NodeData", &node_data);
-
-// msh_element_type.add
-// ("0" , _not_defined )
-// ("1" , _segment_2 )
-// ("2" , _triangle_3 )
-// ("3" , _quadrangle_4 )
-// ("4" , _tetrahedron_4 )
-// ("5" , _hexahedron_8 )
-// ("6" , _pentahedron_6 )
-// ("7" , _not_defined )
-// ("8" , _segment_3 )
-// ("9" , _triangle_6 )
-// ("10", _not_defined )
-// ("11", _tetrahedron_10)
-// ("12", _not_defined )
-// ("13", _not_defined )
-// ("14", _hexahedron_20 )
-// ("15", _pentahedron_15)
-// ("16", _not_defined )
-// ("17", _point_1 )
-// ("18", _quadrangle_8 );
-
-// mesh_format .name("MeshFormat" );
-// nodes .name("Nodes" );
-// elements .name("Elements" );
-// physical_names .name("PhysicalNames" );
-// interpolation_scheme.name("InterpolationScheme");
-// element_data .name("ElementData" );
-// node_data .name("NodeData" );
-
-// clang-format on
-// }
-
-// /* ----------------------------------------------------------------------
-// */
-// /* Rules */
-// /* ----------------------------------------------------------------------
-// */
-// private:
-// qi::symbols<char, ElementType> msh_element_type;
-// qi::symbols<char, qi::rule<Iterator, void(), Skipper> *> sections;
-// qi::rule<Iterator, void(), Skipper> start;
-// qi::rule<Iterator, void(), Skipper, qi::locals<std::string> >
-// unknown_section;
-// qi::rule<Iterator, void(), Skipper, qi::locals<qi::rule<Iterator,
-// Skipper> *> > known_section;
-// qi::rule<Iterator, void(), Skipper> mesh_format, nodes, elements,
-// physical_names, ignore_section,
-// interpolation_scheme, element_data, node_data, any_line;
-// qi::rule<Iterator, void(), Skipper, qi::locals<int> > nodes_;
-// qi::rule<Iterator, void(), Skipper, qi::locals< int, int, vector<int>,
-// vector<int> > > elements_;
-
-// qi::rule<Iterator, std::vector<int>(), Skipper> version;
-// qi::rule<Iterator, _Element(), Skipper, qi::locals<ElementType> >
-// element;
-// qi::rule<Iterator, std::vector<int>(), Skipper, qi::locals<int> > tags;
-// qi::rule<Iterator, std::vector<int>(int), Skipper> connectivity;
-// qi::rule<Iterator, std::vector<Real>(), Skipper> position;
-
-// qi::real_parser<Real, qi::real_policies<Real> > real;
-
-// /* ----------------------------------------------------------------------
-// */
-// /* Members */
-// /* ----------------------------------------------------------------------
-// */
-// private:
-// /// reference to the mesh to read
-// Mesh & mesh;
-
-// /// correspondance between the numbering of nodes in the abaqus file and
-// in
-// /// the akantu mesh
-// std::map<UInt, UInt> msh_nodes_to_akantu;
-
-// /// correspondance between the element number in the abaqus file and the
-// /// Element in the akantu mesh
-// std::map<UInt, Element> msh_elemt_to_akantu;
-// };
-// }
-
-// /* --------------------------------------------------------------------------
-// */
-// void MeshIOAbaqus::read(const std::string& filename, Mesh& mesh) {
-// namespace qi = boost::spirit::qi;
-// namespace ascii = boost::spirit::ascii;
-
-// std::ifstream infile;
-// infile.open(filename.c_str());
-
-// if(!infile.good()) {
-// AKANTU_ERROR("Cannot open file " << filename);
-// }
-
-// std::string storage; // We will read the contents here.
-// infile.unsetf(std::ios::skipws); // No white space skipping!
-// std::copy(std::istream_iterator<char>(infile),
-// std::istream_iterator<char>(),
-// std::back_inserter(storage));
-
-// typedef std::string::const_iterator iterator_t;
-// typedef ascii::space_type skipper;
-// typedef _parser::MshMeshGrammar<iterator_t, skipper> grammar;
-
-// grammar g(mesh);
-// skipper ws;
-
-// iterator_t iter = storage.begin();
-// iterator_t end = storage.end();
-
-// qi::phrase_parse(iter, end, g, ws);
-
-// this->setNbGlobalNodes(mesh, mesh.getNodes().size());
-// MeshUtils::fillElementToSubElementsData(mesh);
-// }
-
static void my_getline(std::ifstream & infile, std::string & str) {
std::string tmp_str;
std::getline(infile, tmp_str);
str = trim(tmp_str);
}
/* -------------------------------------------------------------------------- */
void MeshIOMSH::read(const std::string & filename, Mesh & mesh) {
MeshAccessor mesh_accessor(mesh);
std::ifstream infile;
infile.open(filename.c_str());
std::string line;
UInt first_node_number = std::numeric_limits<UInt>::max(),
last_node_number = 0, file_format = 1, current_line = 0;
bool has_physical_names = false;
if (!infile.good()) {
AKANTU_EXCEPTION("Cannot open file " << filename);
}
std::map<std::string, std::function<void(const std::string &)>> readers;
readers["$MeshFormat"] = [&](const std::string &) {
my_getline(infile, line); /// the format line
std::stringstream sstr(line);
int version;
sstr >> version;
if (version > 2)
AKANTU_ERROR("This reader can not read version "
<< version << " of the MSH file format");
Int format;
sstr >> format;
if (format != 0)
AKANTU_ERROR("This reader can only read ASCII files.");
my_getline(infile, line); /// the end of block line
current_line += 2;
file_format = 2;
};
readers["$NOD"] = readers["$Nodes"] = [&](const std::string &) {
UInt nb_nodes;
my_getline(infile, line);
std::stringstream sstr(line);
sstr >> nb_nodes;
current_line++;
Array<Real> & nodes = mesh_accessor.getNodes();
nodes.resize(nb_nodes);
mesh_accessor.setNbGlobalNodes(nb_nodes);
UInt index;
Real coord[3];
UInt spatial_dimension = nodes.getNbComponent();
/// for each node, read the coordinates
for (UInt i = 0; i < nb_nodes; ++i) {
UInt offset = i * spatial_dimension;
my_getline(infile, line);
std::stringstream sstr_node(line);
sstr_node >> index >> coord[0] >> coord[1] >> coord[2];
current_line++;
first_node_number = std::min(first_node_number, index);
last_node_number = std::max(last_node_number, index);
/// read the coordinates
for (UInt j = 0; j < spatial_dimension; ++j)
nodes.storage()[offset + j] = coord[j];
}
my_getline(infile, line); /// the end of block line
};
readers["$ELM"] = readers["$Elements"] = [&](const std::string &) {
UInt nb_elements;
std::vector<UInt> read_order;
my_getline(infile, line);
std::stringstream sstr(line);
sstr >> nb_elements;
current_line++;
Int index;
UInt msh_type;
ElementType akantu_type, akantu_type_old = _not_defined;
Array<UInt> * connectivity = nullptr;
UInt node_per_element = 0;
for (UInt i = 0; i < nb_elements; ++i) {
my_getline(infile, line);
std::stringstream sstr_elem(line);
current_line++;
sstr_elem >> index;
sstr_elem >> msh_type;
/// get the connectivity vector depending on the element type
akantu_type =
this->_msh_to_akantu_element_types[(MSHElementType)msh_type];
if (akantu_type == _not_defined) {
AKANTU_DEBUG_WARNING("Unsuported element kind "
<< msh_type << " at line " << current_line);
continue;
}
if (akantu_type != akantu_type_old) {
connectivity = &mesh_accessor.getConnectivity(akantu_type);
// connectivity->resize(0);
node_per_element = connectivity->getNbComponent();
akantu_type_old = akantu_type;
read_order = this->_read_order[akantu_type];
}
/// read tags informations
if (file_format == 2) {
UInt nb_tags;
sstr_elem >> nb_tags;
for (UInt j = 0; j < nb_tags; ++j) {
Int tag;
sstr_elem >> tag;
std::stringstream sstr_tag_name;
sstr_tag_name << "tag_" << j;
Array<UInt> & data = mesh.getDataPointer<UInt>(
sstr_tag_name.str(), akantu_type, _not_ghost);
data.push_back(tag);
}
} else if (file_format == 1) {
Int tag;
sstr_elem >> tag; // reg-phys
std::string tag_name = "tag_0";
Array<UInt> * data =
&mesh.getDataPointer<UInt>(tag_name, akantu_type, _not_ghost);
data->push_back(tag);
sstr_elem >> tag; // reg-elem
tag_name = "tag_1";
data = &mesh.getDataPointer<UInt>(tag_name, akantu_type, _not_ghost);
data->push_back(tag);
sstr_elem >> tag; // number-of-nodes
}
Vector<UInt> local_connect(node_per_element);
for (UInt j = 0; j < node_per_element; ++j) {
UInt node_index;
sstr_elem >> node_index;
AKANTU_DEBUG_ASSERT(node_index <= last_node_number,
"Node number not in range : line " << current_line);
node_index -= first_node_number;
local_connect(read_order[j]) = node_index;
}
connectivity->push_back(local_connect);
}
my_getline(infile, line); /// the end of block line
};
readers["$PhysicalNames"] = [&](const std::string &) {
has_physical_names = true;
my_getline(infile, line); /// the format line
std::stringstream sstr(line);
UInt num_of_phys_names;
sstr >> num_of_phys_names;
for (UInt k(0); k < num_of_phys_names; k++) {
my_getline(infile, line);
std::stringstream sstr_phys_name(line);
UInt phys_name_id;
UInt phys_dim;
sstr_phys_name >> phys_dim >> phys_name_id;
std::size_t b = line.find('\"');
std::size_t e = line.rfind('\"');
std::string phys_name = line.substr(b + 1, e - b - 1);
phys_name_map[phys_name_id] = phys_name;
}
my_getline(infile, line); /// the end of block line
};
readers["$Periodic"] = [&](const std::string &) {
UInt nb_periodic_entities;
my_getline(infile, line);
std::stringstream sstr(line);
sstr >> nb_periodic_entities;
mesh_accessor.getNodesFlags().resize(mesh.getNbNodes(),
NodeFlag::_normal);
for (UInt p = 0; p < nb_periodic_entities; ++p) {
// dimension slave-tag master-tag
my_getline(infile, line);
UInt dimension;
{
std::stringstream sstr(line);
sstr >> dimension;
}
// transformation
my_getline(infile, line);
// nb nodes
my_getline(infile, line);
UInt nb_nodes;
{
std::stringstream sstr(line);
sstr >> nb_nodes;
}
for (UInt n = 0; n < nb_nodes; ++n) {
// slave master
my_getline(infile, line);
// The info in the mesh seem inconsistent so they are ignored for know.
continue;
if (dimension == mesh.getSpatialDimension() - 1) {
UInt slave, master;
std::stringstream sstr(line);
sstr >> slave;
sstr >> master;
mesh_accessor.addPeriodicSlave(slave, master);
}
}
}
// mesh_accessor.markMeshPeriodic();
my_getline(infile, line);
};
readers["$NodeData"] = [&](const std::string &) {
/* $NodeData
number-of-string-tags
< "string-tag" >
…
number-of-real-tags
< real-tag >
…
number-of-integer-tags
< integer-tag >
…
node-number value …
…
$EndNodeData */
auto read_data_tags = [&](auto && tags) {
my_getline(infile, line); /// number of tags
UInt number_of_tags{0};
std::stringstream sstr(line);
sstr >> number_of_tags;
tags.resize(number_of_tags);
for (auto && tag : tags) {
my_getline(infile, line);
std::stringstream sstr(line);
sstr >> tag;
}
return tags;
};
auto && string_tags = read_data_tags(std::vector<std::string>());
auto && real_tags[[gnu::unused]] = read_data_tags(std::vector<double>());
auto && int_tags = read_data_tags(std::vector<int>());
auto && data = mesh.registerNodalData<double>(trim(string_tags[0], '"'), int_tags[1]);
data.resize(mesh.getNbNodes(), 0.);
- for (auto _[[gnu::unused]] : arange(int_tags[2])) {
+ for (auto n [[gnu::unused]] : arange(int_tags[2])) {
my_getline(infile, line);
std::stringstream sstr(line);
int node;
double value;
sstr >> node;
- sstr >> value;
-
- data[node - first_node_number] = value;
+ for (auto c : arange(int_tags[1])) {
+ sstr >> value;
+ data(node - first_node_number, c) = value;
+ }
}
my_getline(infile, line);
};
readers["Unsupported"] = [&](const std::string & block) {
AKANTU_DEBUG_WARNING("Unsuported block_kind " << line << " at line "
<< current_line);
auto && end_block = "$End" + block;
while (line != end_block) {
my_getline(infile, line);
current_line++;
}
};
while (infile.good()) {
my_getline(infile, line);
current_line++;
auto && it = readers.find(line);
if (it != readers.end()) {
it->second(line);
} else if(line.size() != 0) {
readers["Unsupported"](line);
}
}
// mesh.updateTypesOffsets(_not_ghost);
infile.close();
this->constructPhysicalNames("tag_0", mesh);
if (has_physical_names)
mesh.createGroupsFromMeshData<std::string>("physical_names");
MeshUtils::fillElementToSubElementsData(mesh);
}
/* -------------------------------------------------------------------------- */
void MeshIOMSH::write(const std::string & filename, const Mesh & mesh) {
std::ofstream outfile;
const Array<Real> & nodes = mesh.getNodes();
outfile.open(filename.c_str());
outfile << "$MeshFormat"
<< "\n";
outfile << "2.1 0 8"
<< "\n";
outfile << "$EndMeshFormat"
<< "\n";
outfile << std::setprecision(std::numeric_limits<Real>::digits10);
outfile << "$Nodes"
<< "\n";
outfile << nodes.size() << "\n";
outfile << std::uppercase;
for (UInt i = 0; i < nodes.size(); ++i) {
Int offset = i * nodes.getNbComponent();
outfile << i + 1;
for (UInt j = 0; j < nodes.getNbComponent(); ++j) {
outfile << " " << nodes.storage()[offset + j];
}
for (UInt p = nodes.getNbComponent(); p < 3; ++p)
outfile << " " << 0.;
outfile << "\n";
;
}
outfile << std::nouppercase;
outfile << "$EndNodes"
<< "\n";
outfile << "$Elements"
<< "\n";
Int nb_elements = 0;
for (auto && type :
mesh.elementTypes(_all_dimensions, _not_ghost, _ek_not_defined)) {
const Array<UInt> & connectivity = mesh.getConnectivity(type, _not_ghost);
nb_elements += connectivity.size();
}
outfile << nb_elements << "\n";
UInt element_idx = 1;
for (auto && type :
mesh.elementTypes(_all_dimensions, _not_ghost, _ek_not_defined)) {
const auto & connectivity = mesh.getConnectivity(type, _not_ghost);
UInt * tag[2] = {nullptr, nullptr};
if (mesh.hasData<UInt>("tag_0", type, _not_ghost)) {
const auto & data_tag_0 = mesh.getData<UInt>("tag_0", type, _not_ghost);
tag[0] = data_tag_0.storage();
}
if (mesh.hasData<UInt>("tag_1", type, _not_ghost)) {
const auto & data_tag_1 = mesh.getData<UInt>("tag_1", type, _not_ghost);
tag[1] = data_tag_1.storage();
}
for (UInt i = 0; i < connectivity.size(); ++i) {
UInt offset = i * connectivity.getNbComponent();
outfile << element_idx << " " << _akantu_to_msh_element_types[type]
<< " 2";
/// \todo write the real data in the file
for (UInt t = 0; t < 2; ++t)
if (tag[t])
outfile << " " << tag[t][i];
else
outfile << " 0";
for (UInt j = 0; j < connectivity.getNbComponent(); ++j) {
outfile << " " << connectivity.storage()[offset + j] + 1;
}
outfile << "\n";
element_idx++;
}
}
outfile << "$EndElements"
<< "\n";
if (mesh.hasData(MeshDataType::_nodal)) {
auto && tags = mesh.getTagNames();
for (auto && tag : tags) {
if (mesh.getTypeCode(tag, MeshDataType::_nodal) != _tc_real)
continue;
auto && data = mesh.getNodalData<double>(tag);
outfile << "$NodeData"
<< "\n";
outfile << "1"
<< "\n";
outfile << "\"" << tag << "\"\n";
outfile << "1\n0.0"
<< "\n";
outfile << "3\n0"
<< "\n";
outfile << data.getNbComponent() << "\n";
outfile << data.size() << "\n";
for (auto && d : enumerate(make_view(data, data.getNbComponent()))) {
outfile << std::get<0>(d) + 1;
for (auto && v : std::get<1>(d)) {
outfile << " " << v;
}
outfile << "\n";
}
outfile << "$EndNodeData"
<< "\n";
}
}
outfile.close();
}
/* -------------------------------------------------------------------------- */
} // namespace akantu
diff --git a/src/synchronizer/node_info_per_processor.cc b/src/synchronizer/node_info_per_processor.cc
index c3bae37a2..71fbdbf55 100644
--- a/src/synchronizer/node_info_per_processor.cc
+++ b/src/synchronizer/node_info_per_processor.cc
@@ -1,868 +1,866 @@
/**
* @file node_info_per_processor.cc
*
* @author Nicolas Richart <nicolas.richart@epfl.ch>
*
* @date creation: Wed Mar 16 2016
* @date last modification: Wed Nov 08 2017
*
* @brief Please type the brief for file: Helper classes to create the
* distributed synchronizer and distribute a mesh
*
* @section LICENSE
*
* Copyright (©) 2016-2018 EPFL (Ecole Polytechnique Fédérale de Lausanne)
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* Akantu is free software: you can 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 3 of the License, or (at your option) any
* later version.
*
* Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
* A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Akantu. If not, see <http://www.gnu.org/licenses/>.
*
*/
/* -------------------------------------------------------------------------- */
#include "node_info_per_processor.hh"
#include "communicator.hh"
#include "node_group.hh"
#include "node_synchronizer.hh"
/* -------------------------------------------------------------------------- */
#include <algorithm>
/* -------------------------------------------------------------------------- */
namespace akantu {
/* -------------------------------------------------------------------------- */
NodeInfoPerProc::NodeInfoPerProc(NodeSynchronizer & synchronizer,
UInt message_cnt, UInt root)
: MeshAccessor(synchronizer.getMesh()), synchronizer(synchronizer),
comm(synchronizer.getCommunicator()), rank(comm.whoAmI()),
nb_proc(comm.getNbProc()), root(root), mesh(synchronizer.getMesh()),
spatial_dimension(synchronizer.getMesh().getSpatialDimension()),
message_count(message_cnt) {}
/* -------------------------------------------------------------------------- */
void NodeInfoPerProc::synchronize() {
synchronizeNodes();
synchronizeTypes();
synchronizeGroups();
synchronizePeriodicity();
synchronizeTags();
}
/* -------------------------------------------------------------------------- */
template <class CommunicationBuffer>
void NodeInfoPerProc::fillNodeGroupsFromBuffer(CommunicationBuffer & buffer) {
AKANTU_DEBUG_IN();
std::vector<std::vector<std::string>> node_to_group;
buffer >> node_to_group;
AKANTU_DEBUG_ASSERT(node_to_group.size() == mesh.getNbGlobalNodes(),
"Not the good amount of nodes where transmitted");
const auto & global_nodes = mesh.getGlobalNodesIds();
for (auto && data : enumerate(global_nodes)) {
for (const auto & node : node_to_group[std::get<1>(data)]) {
mesh.getNodeGroup(node).add(std::get<0>(data), false);
}
}
for (auto && ng_data : mesh.iterateNodeGroups()) {
- ng_data.optimize();
+ ng_data.optimize();
}
AKANTU_DEBUG_OUT();
}
/* -------------------------------------------------------------------------- */
void NodeInfoPerProc::fillNodesType() {
AKANTU_DEBUG_IN();
UInt nb_nodes = mesh.getNbNodes();
auto & nodes_flags = this->getNodesFlags();
Array<UInt> nodes_set(nb_nodes);
nodes_set.set(0);
enum NodeSet {
NORMAL_SET = 1,
GHOST_SET = 2,
};
Array<bool> already_seen(nb_nodes, 1, false);
for (auto gt : ghost_types) {
UInt set = NORMAL_SET;
if (gt == _ghost)
set = GHOST_SET;
already_seen.set(false);
for (auto && type :
mesh.elementTypes(_all_dimensions, gt, _ek_not_defined)) {
const auto & connectivity = mesh.getConnectivity(type, gt);
for (auto & conn :
make_view(connectivity, connectivity.getNbComponent())) {
for (UInt n = 0; n < conn.size(); ++n) {
AKANTU_DEBUG_ASSERT(conn(n) < nb_nodes,
"Node " << conn(n)
<< " bigger than number of nodes "
<< nb_nodes);
if (!already_seen(conn(n))) {
nodes_set(conn(n)) += set;
already_seen(conn(n)) = true;
}
}
}
}
}
nodes_flags.resize(nb_nodes);
for (UInt i = 0; i < nb_nodes; ++i) {
if (nodes_set(i) == NORMAL_SET)
nodes_flags(i) = NodeFlag::_normal;
else if (nodes_set(i) == GHOST_SET)
nodes_flags(i) = NodeFlag::_pure_ghost;
else if (nodes_set(i) == (GHOST_SET + NORMAL_SET))
nodes_flags(i) = NodeFlag::_master;
else {
AKANTU_EXCEPTION("Gni ?");
}
}
AKANTU_DEBUG_OUT();
}
/* -------------------------------------------------------------------------- */
void NodeInfoPerProc::fillCommunicationScheme(const Array<UInt> & master_info) {
AKANTU_DEBUG_IN();
Communications<UInt> & communications =
this->synchronizer.getCommunications();
{ // send schemes
auto it = master_info.begin_reinterpret(2, master_info.size() / 2);
auto end = master_info.end_reinterpret(2, master_info.size() / 2);
std::map<UInt, Array<UInt>> send_array_per_proc;
for (; it != end; ++it) {
const Vector<UInt> & send_info = *it;
send_array_per_proc[send_info(0)].push_back(send_info(1));
}
for (auto & send_schemes : send_array_per_proc) {
auto & scheme = communications.createSendScheme(send_schemes.first);
auto & sends = send_schemes.second;
std::sort(sends.begin(), sends.end());
std::transform(sends.begin(), sends.end(), sends.begin(),
[this](UInt g) -> UInt { return mesh.getNodeLocalId(g); });
scheme.copy(sends);
}
}
{ // receive schemes
std::map<UInt, Array<UInt>> recv_array_per_proc;
for (auto node : arange(mesh.getNbNodes())) {
if (mesh.isSlaveNode(node)) {
recv_array_per_proc[mesh.getNodePrank(node)].push_back(
mesh.getNodeGlobalId(node));
}
}
for (auto & recv_schemes : recv_array_per_proc) {
auto & scheme = communications.createRecvScheme(recv_schemes.first);
auto & recvs = recv_schemes.second;
std::sort(recvs.begin(), recvs.end());
std::transform(recvs.begin(), recvs.end(), recvs.begin(),
[this](UInt g) -> UInt { return mesh.getNodeLocalId(g); });
scheme.copy(recvs);
}
}
AKANTU_DEBUG_OUT();
}
/* -------------------------------------------------------------------------- */
void NodeInfoPerProc::fillPeriodicPairs(const Array<UInt> & global_pairs,
std::vector<UInt> & missing_nodes) {
this->wipePeriodicInfo();
auto & nodes_flags = this->getNodesFlags();
auto checkIsLocal = [&](auto && global_node) {
auto && node = mesh.getNodeLocalId(global_node);
if (node == UInt(-1)) {
auto & global_nodes = this->getNodesGlobalIds();
node = global_nodes.size();
global_nodes.push_back(global_node);
nodes_flags.push_back(NodeFlag::_pure_ghost);
missing_nodes.push_back(global_node);
std::cout << "Missing node " << node << std::endl;
}
return node;
};
for (auto && pairs : make_view(global_pairs, 2)) {
UInt slave = checkIsLocal(pairs(0));
UInt master = checkIsLocal(pairs(1));
this->addPeriodicSlave(slave, master);
}
this->markMeshPeriodic();
}
/* -------------------------------------------------------------------------- */
void NodeInfoPerProc::receiveMissingPeriodic(
DynamicCommunicationBuffer & buffer) {
auto & nodes = this->getNodes();
Communications<UInt> & communications =
this->synchronizer.getCommunications();
std::size_t nb_nodes;
buffer >> nb_nodes;
for (auto _[[gnu::unused]] : arange(nb_nodes)) {
Vector<Real> pos(spatial_dimension);
Int prank;
buffer >> pos;
buffer >> prank;
UInt node = nodes.size();
this->setNodePrank(node, prank);
nodes.push_back(pos);
auto & scheme = communications.createRecvScheme(prank);
scheme.push_back(node);
}
while (buffer.getLeftToUnpack() != 0) {
Int prank;
UInt gnode;
buffer >> gnode;
buffer >> prank;
auto node = mesh.getNodeLocalId(gnode);
AKANTU_DEBUG_ASSERT(node != UInt(-1),
"I cannot send the node "
<< gnode << " to proc " << prank
<< " because it is note a local node");
auto & scheme = communications.createSendScheme(prank);
scheme.push_back(node);
}
}
/* -------------------------------------------------------------------------- */
void NodeInfoPerProc::fillNodalData(DynamicCommunicationBuffer & buffer,
std::string tag_name) {
#define AKANTU_DISTRIBUTED_SYNHRONIZER_TAG_DATA(r, _, elem) \
case BOOST_PP_TUPLE_ELEM(2, 0, elem): { \
auto & nodal_data = \
- mesh.getNodalData<BOOST_PP_TUPLE_ELEM(2, 1, elem)>(tag_name); \
+ mesh.getNodalData<BOOST_PP_TUPLE_ELEM(2, 1, elem)>(tag_name); \
nodal_data.resize(mesh.getNbNodes()); \
for (auto && data : make_view(nodal_data)) { \
- for (auto i[[gnu::unused]] : arange(nodal_data.getNbComponent())) { \
- buffer >> data; \
- } \
+ buffer >> data; \
} \
break; \
}
MeshDataTypeCode data_type_code =
mesh.getTypeCode(tag_name, MeshDataType::_nodal);
switch (data_type_code) {
BOOST_PP_SEQ_FOR_EACH(AKANTU_DISTRIBUTED_SYNHRONIZER_TAG_DATA, ,
AKANTU_MESH_DATA_TYPES)
default:
AKANTU_ERROR("Could not obtain the type of tag" << tag_name << "!");
break;
}
#undef AKANTU_DISTRIBUTED_SYNHRONIZER_TAG_DATA
}
/* -------------------------------------------------------------------------- */
/* -------------------------------------------------------------------------- */
/* -------------------------------------------------------------------------- */
MasterNodeInfoPerProc::MasterNodeInfoPerProc(NodeSynchronizer & synchronizer,
UInt message_cnt, UInt root)
: NodeInfoPerProc(synchronizer, message_cnt, root),
all_nodes(0, synchronizer.getMesh().getSpatialDimension()) {
UInt nb_global_nodes = this->mesh.getNbGlobalNodes();
this->comm.broadcast(nb_global_nodes, this->root);
}
/* -------------------------------------------------------------------------- */
void MasterNodeInfoPerProc::synchronizeNodes() {
this->nodes_per_proc.resize(nb_proc);
this->nb_nodes_per_proc.resize(nb_proc);
Array<Real> local_nodes(0, spatial_dimension);
Array<Real> & nodes = this->getNodes();
all_nodes.copy(nodes);
nodes_pranks.resize(nodes.size(), UInt(-1));
for (UInt p = 0; p < nb_proc; ++p) {
UInt nb_nodes = 0;
// UInt * buffer;
Array<Real> * nodes_to_send;
Array<UInt> & nodespp = nodes_per_proc[p];
if (p != root) {
nodes_to_send = new Array<Real>(0, spatial_dimension);
AKANTU_DEBUG_INFO("Receiving number of nodes from proc "
<< p << " " << Tag::genTag(p, 0, Tag::_NB_NODES));
comm.receive(nb_nodes, p, Tag::genTag(p, 0, Tag::_NB_NODES));
nodespp.resize(nb_nodes);
this->nb_nodes_per_proc(p) = nb_nodes;
AKANTU_DEBUG_INFO("Receiving list of nodes from proc "
<< p << " " << Tag::genTag(p, 0, Tag::_NODES));
comm.receive(nodespp, p, Tag::genTag(p, 0, Tag::_NODES));
} else {
Array<UInt> & local_ids = this->getNodesGlobalIds();
this->nb_nodes_per_proc(p) = local_ids.size();
nodespp.copy(local_ids);
nodes_to_send = &local_nodes;
}
Array<UInt>::const_scalar_iterator it = nodespp.begin();
Array<UInt>::const_scalar_iterator end = nodespp.end();
/// get the coordinates for the selected nodes
for (; it != end; ++it) {
Vector<Real> coord(nodes.storage() + spatial_dimension * *it,
spatial_dimension);
nodes_to_send->push_back(coord);
}
if (p != root) { /// send them for distant processors
AKANTU_DEBUG_INFO("Sending coordinates to proc "
<< p << " "
<< Tag::genTag(this->rank, 0, Tag::_COORDINATES));
comm.send(*nodes_to_send, p,
Tag::genTag(this->rank, 0, Tag::_COORDINATES));
delete nodes_to_send;
}
}
/// construct the local nodes coordinates
nodes.copy(local_nodes);
}
/* -------------------------------------------------------------------------- */
void MasterNodeInfoPerProc::synchronizeTypes() {
// <global_id, <proc, local_id> >
std::multimap<UInt, std::pair<UInt, UInt>> nodes_to_proc;
std::vector<Array<NodeFlag>> nodes_flags_per_proc(nb_proc);
std::vector<Array<Int>> nodes_prank_per_proc(nb_proc);
if (mesh.isPeriodic())
all_periodic_flags.copy(this->getNodesFlags());
// arrays containing pairs of (proc, node)
std::vector<Array<UInt>> nodes_to_send_per_proc(nb_proc);
for (UInt p = 0; p < nb_proc; ++p) {
nodes_flags_per_proc[p].resize(nb_nodes_per_proc(p), NodeFlag(0xFF));
nodes_prank_per_proc[p].resize(nb_nodes_per_proc(p), -1);
}
this->fillNodesType();
auto is_master = [](auto && flag) {
return (flag & NodeFlag::_shared_mask) == NodeFlag::_master;
};
auto is_local = [](auto && flag) {
return (flag & NodeFlag::_shared_mask) == NodeFlag::_normal;
};
for (auto p : arange(nb_proc)) {
auto & nodes_flags = nodes_flags_per_proc[p];
if (p != root) {
AKANTU_DEBUG_INFO(
"Receiving first nodes types from proc "
<< p << " "
<< Tag::genTag(this->rank, this->message_count, Tag::_NODES_TYPE));
comm.receive(nodes_flags, p, Tag::genTag(p, 0, Tag::_NODES_TYPE));
} else {
nodes_flags.copy(this->getNodesFlags());
}
// stack all processors claiming to be master for a node
for (auto local_node : arange(nb_nodes_per_proc(p))) {
auto global_node = nodes_per_proc[p](local_node);
if (is_master(nodes_flags(local_node))) {
nodes_to_proc.insert(
std::make_pair(global_node, std::make_pair(p, local_node)));
} else if (is_local(nodes_flags(local_node))) {
nodes_pranks[global_node] = p;
}
}
}
for (auto i : arange(mesh.getNbGlobalNodes())) {
auto it_range = nodes_to_proc.equal_range(i);
if (it_range.first == nodes_to_proc.end() || it_range.first->first != i)
continue;
// pick the first processor out of the multi-map as the actual master
UInt master_proc = (it_range.first)->second.first;
nodes_pranks[i] = master_proc;
for (auto && data : range(it_range.first, it_range.second)) {
auto proc = data.second.first;
auto node = data.second.second;
if (proc != master_proc) {
// store the info on all the slaves for a given master
nodes_flags_per_proc[proc](node) = NodeFlag::_slave;
nodes_to_send_per_proc[master_proc].push_back(proc);
nodes_to_send_per_proc[master_proc].push_back(i);
}
}
}
/// Fills the nodes prank per proc
for (auto && data : zip(arange(nb_proc), nodes_per_proc, nodes_prank_per_proc,
nodes_flags_per_proc)) {
for (auto && node_data :
zip(std::get<1>(data), std::get<2>(data), std::get<3>(data))) {
if (std::get<2>(node_data) == NodeFlag::_normal) {
std::get<1>(node_data) = std::get<0>(data);
} else {
std::get<1>(node_data) = nodes_pranks(std::get<0>(node_data));
}
}
}
std::vector<CommunicationRequest> requests_send_type;
std::vector<CommunicationRequest> requests_send_master_info;
for (UInt p = 0; p < nb_proc; ++p) {
if (p != root) {
AKANTU_DEBUG_INFO("Sending nodes types to proc "
<< p << " "
<< Tag::genTag(this->rank, 0, Tag::_NODES_TYPE));
requests_send_type.push_back(
comm.asyncSend(nodes_flags_per_proc[p], p,
Tag::genTag(this->rank, 0, Tag::_NODES_TYPE)));
requests_send_type.push_back(
comm.asyncSend(nodes_prank_per_proc[p], p,
Tag::genTag(this->rank, 2, Tag::_NODES_TYPE)));
auto & nodes_to_send = nodes_to_send_per_proc[p];
AKANTU_DEBUG_INFO("Sending nodes master info to proc "
<< p << " "
<< Tag::genTag(this->rank, 1, Tag::_NODES_TYPE));
requests_send_master_info.push_back(comm.asyncSend(
nodes_to_send, p, Tag::genTag(this->rank, 1, Tag::_NODES_TYPE)));
} else {
this->getNodesFlags().copy(nodes_flags_per_proc[p]);
for (auto && data : enumerate(nodes_prank_per_proc[p])) {
auto node = std::get<0>(data);
if (not(mesh.isMasterNode(node) or mesh.isLocalNode(node))) {
this->setNodePrank(node, std::get<1>(data));
}
}
this->fillCommunicationScheme(nodes_to_send_per_proc[root]);
}
}
comm.waitAll(requests_send_type);
comm.freeCommunicationRequest(requests_send_type);
requests_send_type.clear();
comm.waitAll(requests_send_master_info);
comm.freeCommunicationRequest(requests_send_master_info);
}
/* -------------------------------------------------------------------------- */
void MasterNodeInfoPerProc::synchronizeGroups() {
AKANTU_DEBUG_IN();
UInt nb_total_nodes = mesh.getNbGlobalNodes();
DynamicCommunicationBuffer buffer;
using NodeToGroup = std::vector<std::vector<std::string>>;
NodeToGroup node_to_group;
node_to_group.resize(nb_total_nodes);
GroupManager::const_node_group_iterator ngi = mesh.node_group_begin();
GroupManager::const_node_group_iterator nge = mesh.node_group_end();
for (; ngi != nge; ++ngi) {
NodeGroup & ng = *(ngi->second);
std::string name = ngi->first;
NodeGroup::const_node_iterator nit = ng.begin();
NodeGroup::const_node_iterator nend = ng.end();
for (; nit != nend; ++nit) {
node_to_group[*nit].push_back(name);
}
nit = ng.begin();
if (nit != nend)
ng.empty();
}
buffer << node_to_group;
std::vector<CommunicationRequest> requests;
for (UInt p = 0; p < nb_proc; ++p) {
if (p == this->rank)
continue;
AKANTU_DEBUG_INFO("Sending node groups to proc "
<< p << " "
<< Tag::genTag(this->rank, p, Tag::_NODE_GROUP));
requests.push_back(comm.asyncSend(
buffer, p, Tag::genTag(this->rank, p, Tag::_NODE_GROUP)));
}
this->fillNodeGroupsFromBuffer(buffer);
comm.waitAll(requests);
comm.freeCommunicationRequest(requests);
requests.clear();
AKANTU_DEBUG_OUT();
}
/* -------------------------------------------------------------------------- */
void MasterNodeInfoPerProc::synchronizePeriodicity() {
bool is_periodic = mesh.isPeriodic();
comm.broadcast(is_periodic, root);
if (not is_periodic)
return;
std::vector<CommunicationRequest> requests;
std::vector<Array<UInt>> periodic_info_to_send_per_proc;
for (auto p : arange(nb_proc)) {
periodic_info_to_send_per_proc.emplace_back(0, 2);
auto && periodic_info = periodic_info_to_send_per_proc.back();
for (UInt proc_local_node : arange(nb_nodes_per_proc(p))) {
UInt global_node = nodes_per_proc[p](proc_local_node);
if ((all_periodic_flags[global_node] & NodeFlag::_periodic_mask) ==
NodeFlag::_periodic_slave) {
periodic_info.push_back(
Vector<UInt>{global_node, mesh.getPeriodicMaster(global_node)});
}
}
if (p == root)
continue;
auto && tag = Tag::genTag(this->rank, p, Tag::_PERIODIC_SLAVES);
AKANTU_DEBUG_INFO("Sending periodic info to proc " << p << " " << tag);
requests.push_back(comm.asyncSend(periodic_info, p, tag));
}
CommunicationStatus status;
std::vector<DynamicCommunicationBuffer> buffers(nb_proc);
std::vector<std::vector<UInt>> proc_missings(nb_proc);
auto nodes_it = all_nodes.begin(spatial_dimension);
for (UInt p = 0; p < nb_proc; ++p) {
auto & proc_missing = proc_missings[p];
if (p != root) {
auto && tag = Tag::genTag(p, 0, Tag::_PERIODIC_NODES);
comm.probe<UInt>(p, tag, status);
proc_missing.resize(status.size());
comm.receive(proc_missing, p, tag);
} else {
fillPeriodicPairs(periodic_info_to_send_per_proc[root], proc_missing);
}
auto & buffer = buffers[p];
buffer.reserve((spatial_dimension * sizeof(Real) + sizeof(Int)) *
proc_missing.size());
buffer << proc_missing.size();
for (auto && node : proc_missing) {
buffer << *(nodes_it + node);
buffer << nodes_pranks(node);
}
}
for (UInt p = 0; p < nb_proc; ++p) {
for (auto && node : proc_missings[p]) {
auto & buffer = buffers[nodes_pranks(node)];
buffer << node;
buffer << p;
}
}
for (UInt p = 0; p < nb_proc; ++p) {
if (p != root) {
auto && tag_send = Tag::genTag(p, 1, Tag::_PERIODIC_NODES);
requests.push_back(comm.asyncSend(buffers[p], p, tag_send));
} else {
receiveMissingPeriodic(buffers[p]);
}
}
comm.waitAll(requests);
comm.freeCommunicationRequest(requests);
requests.clear();
}
/* -------------------------------------------------------------------------- */
void MasterNodeInfoPerProc::fillTagBuffers(
std::vector<DynamicCommunicationBuffer> & buffers,
const std::string & tag_name) {
#define AKANTU_DISTRIBUTED_SYNHRONIZER_TAG_DATA(r, _, elem) \
case BOOST_PP_TUPLE_ELEM(2, 0, elem): { \
auto & nodal_data = \
- mesh.getNodalData<BOOST_PP_TUPLE_ELEM(2, 1, elem)>(tag_name); \
+ mesh.getNodalData<BOOST_PP_TUPLE_ELEM(2, 1, elem)>(tag_name); \
for (auto && data : enumerate(nodes_per_proc)) { \
auto proc = std::get<0>(data); \
auto & nodes = std::get<1>(data); \
auto & buffer = buffers[proc]; \
for (auto & node : nodes) { \
for (auto i : arange(nodal_data.getNbComponent())) { \
buffer << nodal_data(node, i); \
} \
} \
} \
break; \
}
MeshDataTypeCode data_type_code =
mesh.getTypeCode(tag_name, MeshDataType::_nodal);
switch (data_type_code) {
BOOST_PP_SEQ_FOR_EACH(AKANTU_DISTRIBUTED_SYNHRONIZER_TAG_DATA, ,
AKANTU_MESH_DATA_TYPES)
default:
AKANTU_ERROR("Could not obtain the type of tag" << tag_name << "!");
break;
}
#undef AKANTU_DISTRIBUTED_SYNHRONIZER_TAG_DATA
} // namespace akantu
/* -------------------------------------------------------------------------- */
void MasterNodeInfoPerProc::synchronizeTags() {
/// tag info
auto tag_names = mesh.getTagNames();
DynamicCommunicationBuffer tags_buffer;
for (auto && tag_name : tag_names) {
tags_buffer << tag_name;
tags_buffer << mesh.getTypeCode(tag_name, MeshDataType::_nodal);
tags_buffer << mesh.getNbComponent(tag_name);
}
UInt tags_buffer_length = tags_buffer.size();
AKANTU_DEBUG_INFO(
"Broadcasting the size of the information about the mesh data tags: ("
<< tags_buffer_length << ").");
comm.broadcast(tags_buffer_length, root);
if (tags_buffer_length != 0)
comm.broadcast(tags_buffer, root);
for (auto && tag_data : enumerate(tag_names)) {
auto tag_count = std::get<0>(tag_data);
auto & tag_name = std::get<1>(tag_data);
std::vector<DynamicCommunicationBuffer> buffers;
std::vector<CommunicationRequest> requests;
buffers.resize(nb_proc);
fillTagBuffers(buffers, tag_name);
for (auto && data : enumerate(buffers)) {
auto && proc = std::get<0>(data);
auto & buffer = std::get<1>(data);
if (proc == root) {
fillNodalData(buffer, tag_name);
} else {
auto && tag = Tag::genTag(this->rank, tag_count, Tag::_MESH_DATA);
requests.push_back(comm.asyncSend(buffer, proc, tag));
}
}
comm.waitAll(requests);
comm.freeCommunicationRequest(requests);
}
}
/* -------------------------------------------------------------------------- */
/* -------------------------------------------------------------------------- */
/* -------------------------------------------------------------------------- */
SlaveNodeInfoPerProc::SlaveNodeInfoPerProc(NodeSynchronizer & synchronizer,
UInt message_cnt, UInt root)
: NodeInfoPerProc(synchronizer, message_cnt, root) {
UInt nb_global_nodes = 0;
comm.broadcast(nb_global_nodes, root);
this->setNbGlobalNodes(nb_global_nodes);
}
/* -------------------------------------------------------------------------- */
void SlaveNodeInfoPerProc::synchronizeNodes() {
AKANTU_DEBUG_INFO("Sending list of nodes to proc "
<< root << " " << Tag::genTag(this->rank, 0, Tag::_NB_NODES)
<< " " << Tag::genTag(this->rank, 0, Tag::_NODES));
Array<UInt> & local_ids = this->getNodesGlobalIds();
Array<Real> & nodes = this->getNodes();
UInt nb_nodes = local_ids.size();
comm.send(nb_nodes, root, Tag::genTag(this->rank, 0, Tag::_NB_NODES));
comm.send(local_ids, root, Tag::genTag(this->rank, 0, Tag::_NODES));
/* --------<<<<-COORDINATES---------------------------------------------- */
nodes.resize(nb_nodes);
AKANTU_DEBUG_INFO("Receiving coordinates from proc "
<< root << " " << Tag::genTag(root, 0, Tag::_COORDINATES));
comm.receive(nodes, root, Tag::genTag(root, 0, Tag::_COORDINATES));
}
/* -------------------------------------------------------------------------- */
void SlaveNodeInfoPerProc::synchronizeTypes() {
this->fillNodesType();
auto & nodes_types = this->getNodesFlags();
AKANTU_DEBUG_INFO("Sending first nodes types to proc "
<< root << ""
<< Tag::genTag(this->rank, 0, Tag::_NODES_TYPE));
comm.send(nodes_types, root, Tag::genTag(this->rank, 0, Tag::_NODES_TYPE));
AKANTU_DEBUG_INFO("Receiving nodes types from proc "
<< root << " " << Tag::genTag(root, 0, Tag::_NODES_TYPE));
comm.receive(nodes_types, root, Tag::genTag(root, 0, Tag::_NODES_TYPE));
Array<Int> nodes_prank(nodes_types.size());
comm.receive(nodes_prank, root, Tag::genTag(root, 2, Tag::_NODES_TYPE));
for (auto && data : enumerate(nodes_prank)) {
auto node = std::get<0>(data);
if (not(mesh.isMasterNode(node) or mesh.isLocalNode(node))) {
this->setNodePrank(node, std::get<1>(data));
}
}
AKANTU_DEBUG_INFO("Receiving nodes master info from proc "
<< root << " " << Tag::genTag(root, 1, Tag::_NODES_TYPE));
CommunicationStatus status;
comm.probe<UInt>(root, Tag::genTag(root, 1, Tag::_NODES_TYPE), status);
Array<UInt> nodes_master_info(status.size());
if (nodes_master_info.size() > 0)
comm.receive(nodes_master_info, root,
Tag::genTag(root, 1, Tag::_NODES_TYPE));
this->fillCommunicationScheme(nodes_master_info);
}
/* -------------------------------------------------------------------------- */
void SlaveNodeInfoPerProc::synchronizeGroups() {
AKANTU_DEBUG_IN();
AKANTU_DEBUG_INFO("Receiving node groups from proc "
<< root << " "
<< Tag::genTag(root, this->rank, Tag::_NODE_GROUP));
CommunicationStatus status;
comm.probe<char>(root, Tag::genTag(root, this->rank, Tag::_NODE_GROUP),
status);
CommunicationBuffer buffer(status.size());
comm.receive(buffer, root, Tag::genTag(root, this->rank, Tag::_NODE_GROUP));
this->fillNodeGroupsFromBuffer(buffer);
AKANTU_DEBUG_OUT();
}
/* -------------------------------------------------------------------------- */
void SlaveNodeInfoPerProc::synchronizePeriodicity() {
bool is_periodic;
comm.broadcast(is_periodic, root);
if (not is_periodic)
return;
CommunicationStatus status;
auto && tag = Tag::genTag(root, this->rank, Tag::_PERIODIC_SLAVES);
comm.probe<UInt>(root, tag, status);
Array<UInt> periodic_info(status.size() / 2, 2);
comm.receive(periodic_info, root, tag);
std::vector<UInt> proc_missing;
fillPeriodicPairs(periodic_info, proc_missing);
auto && tag_missing_request =
Tag::genTag(this->rank, 0, Tag::_PERIODIC_NODES);
comm.send(proc_missing, root, tag_missing_request);
DynamicCommunicationBuffer buffer;
auto && tag_missing = Tag::genTag(this->rank, 1, Tag::_PERIODIC_NODES);
comm.receive(buffer, root, tag_missing);
receiveMissingPeriodic(buffer);
}
/* -------------------------------------------------------------------------- */
void SlaveNodeInfoPerProc::synchronizeTags() {
UInt tags_buffer_length{0};
comm.broadcast(tags_buffer_length, root);
if (tags_buffer_length == 0) {
return;
}
CommunicationBuffer tags_buffer;
tags_buffer.resize(tags_buffer_length);
comm.broadcast(tags_buffer, root);
std::vector<std::string> tag_names;
while (tags_buffer.getLeftToUnpack() > 0) {
std::string name;
MeshDataTypeCode code;
UInt nb_components;
tags_buffer >> name;
tags_buffer >> code;
tags_buffer >> nb_components;
mesh.registerNodalData(name, nb_components, code);
tag_names.push_back(name);
}
for (auto && tag_data : enumerate(tag_names)) {
auto tag_count = std::get<0>(tag_data);
auto & tag_name = std::get<1>(tag_data);
DynamicCommunicationBuffer buffer;
auto && tag = Tag::genTag(this->root, tag_count, Tag::_MESH_DATA);
comm.receive(buffer, this->root, tag);
fillNodalData(buffer, tag_name);
}
}
-} // akantu
+} // namespace akantu