Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F79254497
mesh_io_msh.cc
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Sun, Aug 25, 12:32
Size
32 KB
Mime Type
text/x-c
Expires
Tue, Aug 27, 12:32 (2 d)
Engine
blob
Format
Raw Data
Handle
20061848
Attached To
rAKA akantu
mesh_io_msh.cc
View Options
/**
* @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: Mon Dec 07 2015
*
* @brief Read/Write for MSH files generated by gmsh
*
* @section LICENSE
*
* Copyright (©) 2010-2012, 2014, 2015 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>
/* -------------------------------------------------------------------------- */
__BEGIN_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
[
_kirchhoff_shell
]
=
_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
()
{}
/* -------------------------------------------------------------------------- */
/* 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.getSize());
// 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_DEBUG_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().getSize());
// 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
)
{
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
;
if
(
!
infile
.
good
())
{
AKANTU_DEBUG_ERROR
(
"Cannot open file "
<<
filename
);
}
while
(
infile
.
good
())
{
my_getline
(
infile
,
line
);
current_line
++
;
/// read the header
if
(
line
==
"$MeshFormat"
)
{
my_getline
(
infile
,
line
);
/// the format line
std
::
stringstream
sstr
(
line
);
std
::
string
version
;
sstr
>>
version
;
Int
format
;
sstr
>>
format
;
if
(
format
!=
0
)
AKANTU_DEBUG_ERROR
(
"This reader can only read ASCII files."
);
my_getline
(
infile
,
line
);
/// the end of block line
current_line
+=
2
;
file_format
=
2
;
}
/// read the physical names
if
(
line
==
"$PhysicalNames"
)
{
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
;
}
}
/// read all nodes
if
(
line
==
"$Nodes"
||
line
==
"$NOD"
)
{
UInt
nb_nodes
;
my_getline
(
infile
,
line
);
std
::
stringstream
sstr
(
line
);
sstr
>>
nb_nodes
;
current_line
++
;
Array
<
Real
>
&
nodes
=
const_cast
<
Array
<
Real
>
&>
(
mesh
.
getNodes
());
nodes
.
resize
(
nb_nodes
);
mesh
.
nb_global_nodes
=
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
}
/// read all elements
if
(
line
==
"$Elements"
||
line
==
"$ELM"
)
{
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
=
NULL
;
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
.
getConnectivityPointer
(
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
}
UInt
*
local_connect
=
new
UInt
[
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
);
delete
[]
local_connect
;
}
my_getline
(
infile
,
line
);
/// the end of block line
}
if
((
line
[
0
]
==
'$'
)
&&
(
line
.
find
(
"End"
)
==
std
::
string
::
npos
))
{
AKANTU_DEBUG_WARNING
(
"Unsuported block_kind "
<<
line
<<
" at line "
<<
current_line
);
}
}
mesh
.
updateTypesOffsets
(
_not_ghost
);
infile
.
close
();
this
->
constructPhysicalNames
(
"tag_0"
,
mesh
);
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"
<<
std
::
endl
;
outfile
<<
"2.1 0 8"
<<
std
::
endl
;
outfile
<<
"$EndMeshFormat"
<<
std
::
endl
;
outfile
<<
std
::
setprecision
(
std
::
numeric_limits
<
Real
>::
digits10
);
outfile
<<
"$Nodes"
<<
std
::
endl
;
outfile
<<
nodes
.
getSize
()
<<
std
::
endl
;
outfile
<<
std
::
uppercase
;
for
(
UInt
i
=
0
;
i
<
nodes
.
getSize
();
++
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
<<
std
::
endl
;
;
}
outfile
<<
std
::
nouppercase
;
outfile
<<
"$EndNodes"
<<
std
::
endl
;
;
outfile
<<
"$Elements"
<<
std
::
endl
;
;
Mesh
::
type_iterator
it
=
mesh
.
firstType
(
_all_dimensions
,
_not_ghost
,
_ek_not_defined
);
Mesh
::
type_iterator
end
=
mesh
.
lastType
(
_all_dimensions
,
_not_ghost
,
_ek_not_defined
);
Int
nb_elements
=
0
;
for
(;
it
!=
end
;
++
it
)
{
const
Array
<
UInt
>
&
connectivity
=
mesh
.
getConnectivity
(
*
it
,
_not_ghost
);
nb_elements
+=
connectivity
.
getSize
();
}
outfile
<<
nb_elements
<<
std
::
endl
;
UInt
element_idx
=
1
;
for
(
it
=
mesh
.
firstType
(
_all_dimensions
,
_not_ghost
,
_ek_not_defined
);
it
!=
end
;
++
it
)
{
ElementType
type
=
*
it
;
const
Array
<
UInt
>
&
connectivity
=
mesh
.
getConnectivity
(
type
,
_not_ghost
);
UInt
*
tag
[
2
]
=
{
NULL
,
NULL
};
try
{
const
Array
<
UInt
>
&
data_tag_0
=
mesh
.
getData
<
UInt
>
(
"tag_0"
,
type
,
_not_ghost
);
tag
[
0
]
=
data_tag_0
.
storage
();
}
catch
(...)
{
tag
[
0
]
=
NULL
;
}
try
{
const
Array
<
UInt
>
&
data_tag_1
=
mesh
.
getData
<
UInt
>
(
"tag_1"
,
type
,
_not_ghost
);
tag
[
1
]
=
data_tag_1
.
storage
();
}
catch
(...)
{
tag
[
1
]
=
NULL
;
}
for
(
UInt
i
=
0
;
i
<
connectivity
.
getSize
();
++
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
<<
std
::
endl
;
element_idx
++
;
}
}
outfile
<<
"$EndElements"
<<
std
::
endl
;
;
outfile
.
close
();
}
/* -------------------------------------------------------------------------- */
__END_AKANTU__
Event Timeline
Log In to Comment