diff --git a/src/common/aka_common.hh b/src/common/aka_common.hh index 7a0bdf1e0..012b30183 100644 --- a/src/common/aka_common.hh +++ b/src/common/aka_common.hh @@ -1,519 +1,520 @@ /** * @file aka_common.hh * * @author Guillaume Anciaux * @author Nicolas Richart * * @date creation: Mon Jun 14 2010 * @date last modification: Mon Feb 12 2018 * * @brief common type descriptions for akantu * * @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 . * * @section DESCRIPTION * * All common things to be included in the projects files * */ /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_COMMON_HH__ #define __AKANTU_COMMON_HH__ #include "aka_compatibilty_with_cpp_standard.hh" /* -------------------------------------------------------------------------- */ #define __BEGIN_AKANTU_DUMPER__ namespace dumper { #define __END_AKANTU_DUMPER__ } /* -------------------------------------------------------------------------- */ #if defined(WIN32) #define __attribute__(x) #endif /* -------------------------------------------------------------------------- */ #include "aka_config.hh" #include "aka_error.hh" #include "aka_safe_enum.hh" /* -------------------------------------------------------------------------- */ #include #include #include #include #include #include /* -------------------------------------------------------------------------- */ namespace akantu { /* -------------------------------------------------------------------------- */ /* Common types */ /* -------------------------------------------------------------------------- */ using ID = std::string; #ifdef AKANTU_NDEBUG static const Real REAL_INIT_VALUE = Real(0.); #else static const Real REAL_INIT_VALUE = std::numeric_limits::quiet_NaN(); #endif /* -------------------------------------------------------------------------- */ /* Memory types */ /* -------------------------------------------------------------------------- */ using MemoryID = UInt; // using Surface = std::string; // using SurfacePair= std::pair; // using SurfacePairList = std::list; /* -------------------------------------------------------------------------- */ extern const UInt _all_dimensions; #define AKANTU_PP_ENUM(s, data, i, elem) \ BOOST_PP_TUPLE_REM() \ elem BOOST_PP_COMMA_IF(BOOST_PP_NOT_EQUAL(i, BOOST_PP_DEC(data))) } // namespace akantu #if (defined(__GNUC__) || defined(__GNUG__)) #define AKA_GCC_VERSION \ (__GNUC__ * 10000 + __GNUC_MINOR__ * 100 + __GNUC_PATCHLEVEL__) #if AKA_GCC_VERSION < 60000 #define AKANTU_ENUM_HASH(type_name) \ namespace std { \ template <> struct hash<::akantu::type_name> { \ using argument_type = ::akantu::type_name; \ size_t operator()(const argument_type & e) const noexcept { \ auto ue = underlying_type_t(e); \ return uh(ue); \ } \ \ private: \ const hash> uh{}; \ }; \ } #else #define AKANTU_ENUM_HASH(type_name) #endif // AKA_GCC_VERSION #endif // GNU #include "aka_element_classes_info.hh" namespace akantu { #define AKANTU_PP_CAT(s, data, elem) BOOST_PP_CAT(data, elem) #define AKANTU_PP_TYPE_TO_STR(s, data, elem) \ ({BOOST_PP_CAT(data::_, elem), BOOST_PP_STRINGIZE(elem)}) #define AKANTU_PP_STR_TO_TYPE(s, data, elem) \ ({BOOST_PP_STRINGIZE(elem), BOOST_PP_CAT(data::_, elem)}) #define AKANTU_ENUM_DECLARE(type_name, list) \ enum class type_name { \ BOOST_PP_SEQ_ENUM(BOOST_PP_SEQ_TRANSFORM(AKANTU_PP_CAT, _, list)) \ }; #define AKANTU_ENUM_OUTPUT_STREAM(type_name, list) \ } \ AKANTU_ENUM_HASH(type_name) \ namespace aka { \ inline std::string to_string(const ::akantu::type_name & type) { \ static std::unordered_map<::akantu::type_name, std::string> convert{ \ BOOST_PP_SEQ_FOR_EACH_I( \ AKANTU_PP_ENUM, BOOST_PP_SEQ_SIZE(list), \ BOOST_PP_SEQ_TRANSFORM(AKANTU_PP_TYPE_TO_STR, \ ::akantu::type_name, list))}; \ return convert.at(type); \ } \ } \ namespace akantu { \ inline std::ostream & operator<<(std::ostream & stream, \ const type_name & type) { \ stream << aka::to_string(type); \ return stream; \ } #define AKANTU_ENUM_INPUT_STREAM(type_name, list) \ inline std::istream & operator>>(std::istream & stream, type_name & type) { \ std::string str; \ stream >> str; \ static std::unordered_map convert{ \ BOOST_PP_SEQ_FOR_EACH_I( \ AKANTU_PP_ENUM, BOOST_PP_SEQ_SIZE(list), \ BOOST_PP_SEQ_TRANSFORM(AKANTU_PP_STR_TO_TYPE, type_name, list))}; \ type = convert.at(str); \ return stream; \ } /* -------------------------------------------------------------------------- */ /* Mesh/FEM/Model types */ /* -------------------------------------------------------------------------- */ /// small help to use names for directions enum SpacialDirection { _x = 0, _y = 1, _z = 2 }; /// enum MeshIOType type of mesh reader/writer enum MeshIOType { _miot_auto, ///< Auto guess of the reader to use based on the extension _miot_gmsh, ///< Gmsh files _miot_gmsh_struct, ///< Gsmh reader with reintpretation of elements has /// structures elements _miot_diana, ///< TNO Diana mesh format _miot_abaqus ///< Abaqus mesh format }; /// enum MeshEventHandlerPriority defines relative order of execution of events enum EventHandlerPriority { _ehp_highest = 0, _ehp_mesh = 5, _ehp_fe_engine = 9, _ehp_synchronizer = 10, _ehp_dof_manager = 20, _ehp_model = 94, _ehp_non_local_manager = 100, _ehp_lowest = 100 }; #ifndef SWIG // clang-format off #define AKANTU_MODEL_TYPES \ (model) \ (solid_mechanics_model) \ (solid_mechanics_model_cohesive) \ (heat_transfer_model) \ (structural_mechanics_model) \ (embedded_model) // clang-format on /// enum ModelType defines which type of physics is solved AKANTU_ENUM_DECLARE(ModelType, AKANTU_MODEL_TYPES) AKANTU_ENUM_OUTPUT_STREAM(ModelType, AKANTU_MODEL_TYPES) AKANTU_ENUM_INPUT_STREAM(ModelType, AKANTU_MODEL_TYPES) #else enum class ModelType { _model, _solid_mechanics_model, _solid_mechanics_model_cohesive, _heat_transfer_model, _structural_mechanics_model, _embedded_model }; #endif /// enum AnalysisMethod type of solving method used to solve the equation of /// motion enum AnalysisMethod { _static = 0, _implicit_dynamic = 1, _explicit_lumped_mass = 2, _explicit_lumped_capacity = 2, _explicit_consistent_mass = 3 }; /// enum DOFSupportType defines which kind of dof that can exists enum DOFSupportType { _dst_nodal, _dst_generic }; /// Type of non linear resolution available in akantu enum NonLinearSolverType { _nls_linear, ///< No non linear convergence loop _nls_newton_raphson, ///< Regular Newton-Raphson _nls_newton_raphson_modified, ///< Newton-Raphson with initial tangent _nls_lumped, ///< Case of lumped mass or equivalent matrix _nls_auto ///< This will take a default value that make sense in case of /// model::getNewSolver }; /// Define the node/dof type enum NodeType : Int { _nt_pure_ghost = -3, _nt_master = -2, _nt_normal = -1 }; /// Type of time stepping solver enum TimeStepSolverType { _tsst_static, ///< Static solution _tsst_dynamic, ///< Dynamic solver _tsst_dynamic_lumped, ///< Dynamic solver with lumped mass _tsst_not_defined, ///< For not defined cases }; /// Type of integration scheme enum IntegrationSchemeType { _ist_pseudo_time, ///< Pseudo Time _ist_forward_euler, ///< GeneralizedTrapezoidal(0) _ist_trapezoidal_rule_1, ///< GeneralizedTrapezoidal(1/2) _ist_backward_euler, ///< GeneralizedTrapezoidal(1) _ist_central_difference, ///< NewmarkBeta(0, 1/2) _ist_fox_goodwin, ///< NewmarkBeta(1/6, 1/2) _ist_trapezoidal_rule_2, ///< NewmarkBeta(1/2, 1/2) _ist_linear_acceleration, ///< NewmarkBeta(1/3, 1/2) _ist_newmark_beta, ///< generic NewmarkBeta with user defined /// alpha and beta _ist_generalized_trapezoidal ///< generic GeneralizedTrapezoidal with user /// defined alpha }; /// enum SolveConvergenceCriteria different convergence criteria enum SolveConvergenceCriteria { _scc_residual, ///< Use residual to test the convergence _scc_solution, ///< Use solution to test the convergence _scc_residual_mass_wgh ///< Use residual weighted by inv. nodal mass to testb }; /// enum CohesiveMethod type of insertion of cohesive elements enum CohesiveMethod { _intrinsic, _extrinsic }; /// @enum SparseMatrixType type of sparse matrix used enum MatrixType { _unsymmetric, _symmetric, _mt_not_defined }; /* -------------------------------------------------------------------------- */ /* Ghosts handling */ /* -------------------------------------------------------------------------- */ /// @enum CommunicatorType type of communication method to use enum CommunicatorType { _communicator_mpi, _communicator_dummy }; /// @enum SynchronizationTag type of synchronizations enum SynchronizationTag { //--- Generic tags --- _gst_whatever, _gst_update, _gst_ask_nodes, _gst_size, //--- SolidMechanicsModel tags --- _gst_smm_mass, ///< synchronization of the SolidMechanicsModel.mass _gst_smm_for_gradu, ///< synchronization of the /// SolidMechanicsModel.displacement _gst_smm_boundary, ///< synchronization of the boundary, forces, velocities /// and displacement _gst_smm_uv, ///< synchronization of the nodal velocities and displacement _gst_smm_res, ///< synchronization of the nodal residual _gst_smm_init_mat, ///< synchronization of the data to initialize materials _gst_smm_stress, ///< synchronization of the stresses to compute the internal /// forces _gst_smmc_facets, ///< synchronization of facet data to setup facet synch _gst_smmc_facets_conn, ///< synchronization of facet global connectivity _gst_smmc_facets_stress, ///< synchronization of facets' stress to setup facet /// synch _gst_smmc_damage, ///< synchronization of damage // --- GlobalIdsUpdater tags --- _gst_giu_global_conn, ///< synchronization of global connectivities // --- CohesiveElementInserter tags --- _gst_ce_groups, ///< synchronization of cohesive element insertion depending /// on facet groups // --- GroupManager tags --- _gst_gm_clusters, ///< synchronization of clusters // --- HeatTransfer tags --- _gst_htm_temperature, ///< synchronization of the nodal temperature _gst_htm_gradient_temperature, ///< synchronization of the element gradient /// temperature // --- LevelSet tags --- _gst_htm_phi, ///< synchronization of the nodal level set value phi _gst_htm_gradient_phi, ///< synchronization of the element gradient phi //--- Material non local --- _gst_mnl_for_average, ///< synchronization of data to average in non local /// material _gst_mnl_weight, ///< synchronization of data for the weight computations // --- NeighborhoodSynchronization tags --- _gst_nh_criterion, // --- General tags --- _gst_test, ///< Test tag _gst_user_1, ///< tag for user simulations _gst_user_2, ///< tag for user simulations _gst_material_id, ///< synchronization of the material ids _gst_for_dump, ///< everything that needs to be synch before dump // --- Contact & Friction --- _gst_cf_nodal, ///< synchronization of disp, velo, and current position _gst_cf_incr, ///< synchronization of increment // --- Solver tags --- _gst_solver_solution ///< synchronization of the solution obained with the /// PETSc solver }; /// standard output stream operator for SynchronizationTag inline std::ostream & operator<<(std::ostream & stream, SynchronizationTag type); /// @enum GhostType type of ghost enum GhostType { _not_ghost = 0, _ghost = 1, _casper // not used but a real cute ghost }; } // namespace akantu #ifndef SWIG AKANTU_ENUM_HASH(GhostType) #endif namespace akantu { /* -------------------------------------------------------------------------- */ struct GhostType_def { using type = GhostType; static const type _begin_ = _not_ghost; static const type _end_ = _casper; }; using ghost_type_t = safe_enum; extern ghost_type_t ghost_types; /// standard output stream operator for GhostType inline std::ostream & operator<<(std::ostream & stream, GhostType type); /* -------------------------------------------------------------------------- */ /* Global defines */ /* -------------------------------------------------------------------------- */ #define AKANTU_MIN_ALLOCATION 2000 #define AKANTU_INDENT " " #define AKANTU_INCLUDE_INLINE_IMPL /* -------------------------------------------------------------------------- */ /* Type traits */ /* -------------------------------------------------------------------------- */ struct TensorTrait {}; /* -------------------------------------------------------------------------- */ template using is_tensor = std::is_base_of; /* -------------------------------------------------------------------------- */ template using is_scalar = std::is_arithmetic; /* -------------------------------------------------------------------------- */ /* -------------------------------------------------------------------------- */ #define AKANTU_SET_MACRO(name, variable, type) \ inline void set##name(type variable) { this->variable = variable; } #define AKANTU_GET_MACRO(name, variable, type) \ inline type get##name() const { return variable; } #define AKANTU_GET_MACRO_NOT_CONST(name, variable, type) \ inline type get##name() { return variable; } #define AKANTU_GET_MACRO_BY_SUPPORT_TYPE(name, variable, type, support, con) \ inline con Array & get##name( \ const support & el_type, const GhostType & ghost_type = _not_ghost) \ con { \ return variable(el_type, ghost_type); \ } #define AKANTU_GET_MACRO_BY_ELEMENT_TYPE(name, variable, type) \ AKANTU_GET_MACRO_BY_SUPPORT_TYPE(name, variable, type, ElementType, ) #define AKANTU_GET_MACRO_BY_ELEMENT_TYPE_CONST(name, variable, type) \ AKANTU_GET_MACRO_BY_SUPPORT_TYPE(name, variable, type, ElementType, const) #define AKANTU_GET_MACRO_BY_GEOMETRIE_TYPE(name, variable, type) \ AKANTU_GET_MACRO_BY_SUPPORT_TYPE(name, variable, type, GeometricalType, ) #define AKANTU_GET_MACRO_BY_GEOMETRIE_TYPE_CONST(name, variable, type) \ AKANTU_GET_MACRO_BY_SUPPORT_TYPE(name, variable, type, GeometricalType, const) /* -------------------------------------------------------------------------- */ /// initialize the static part of akantu void initialize(int & argc, char **& argv); /// initialize the static part of akantu and read the global input_file void initialize(const std::string & input_file, int & argc, char **& argv); /* -------------------------------------------------------------------------- */ /// finilize correctly akantu and clean the memory void finalize(); /* -------------------------------------------------------------------------- */ /// Read an new input file void readInputFile(const std::string & input_file); /* -------------------------------------------------------------------------- */ /* * For intel compiler annoying remark */ // #if defined(__INTEL_COMPILER) // /// remark #981: operands are evaluated in unspecified order // #pragma warning(disable : 981) // /// remark #383: value copied to temporary, reference to temporary used // #pragma warning(disable : 383) // #endif // defined(__INTEL_COMPILER) /* -------------------------------------------------------------------------- */ /* string manipulation */ /* -------------------------------------------------------------------------- */ inline std::string to_lower(const std::string & str); /* -------------------------------------------------------------------------- */ inline std::string trim(const std::string & to_trim); +inline std::string trim(const std::string & to_trim, char c); /* -------------------------------------------------------------------------- */ /* -------------------------------------------------------------------------- */ /// give a string representation of the a human readable size in bit template std::string printMemorySize(UInt size); /* -------------------------------------------------------------------------- */ } // namespace akantu #include "aka_fwd.hh" namespace akantu { /// get access to the internal argument parser cppargparse::ArgumentParser & getStaticArgumentParser(); /// get access to the internal input file parser Parser & getStaticParser(); /// get access to the user part of the internal input file parser const ParserSection & getUserParser(); } // namespace akantu #include "aka_common_inline_impl.cc" /* -------------------------------------------------------------------------- */ #if AKANTU_INTEGER_SIZE == 4 #define AKANTU_HASH_COMBINE_MAGIC_NUMBER 0x9e3779b9 #elif AKANTU_INTEGER_SIZE == 8 #define AKANTU_HASH_COMBINE_MAGIC_NUMBER 0x9e3779b97f4a7c13LL #endif namespace std { /** * Hashing function for pairs based on hash_combine from boost The magic number * is coming from the golden number @f[\phi = \frac{1 + \sqrt5}{2}@f] * @f[\frac{2^32}{\phi} = 0x9e3779b9@f] * http://stackoverflow.com/questions/4948780/magic-number-in-boosthash-combine * http://burtleburtle.net/bob/hash/doobs.html */ template struct hash> { hash() = default; size_t operator()(const std::pair & p) const { size_t seed = ah(p.first); return bh(p.second) + AKANTU_HASH_COMBINE_MAGIC_NUMBER + (seed << 6) + (seed >> 2); } private: const hash ah{}; const hash bh{}; }; } // namespace std #endif /* __AKANTU_COMMON_HH__ */ diff --git a/src/common/aka_common_inline_impl.cc b/src/common/aka_common_inline_impl.cc index 3a0023066..9b4f1edee 100644 --- a/src/common/aka_common_inline_impl.cc +++ b/src/common/aka_common_inline_impl.cc @@ -1,458 +1,466 @@ /** * @file aka_common_inline_impl.cc * * @author Guillaume Anciaux * @author Nicolas Richart * * @date creation: Fri Jun 18 2010 * @date last modification: Tue Feb 20 2018 * * @brief inline implementations of common akantu type descriptions * * @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 . * * @section DESCRIPTION * * All common things to be included in the projects files * */ /* -------------------------------------------------------------------------- */ #include "aka_common.hh" #include #include #include #include /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_AKA_COMMON_INLINE_IMPL_CC__ #define __AKANTU_AKA_COMMON_INLINE_IMPL_CC__ namespace akantu { /* -------------------------------------------------------------------------- */ /// standard output stream operator for GhostType inline std::ostream & operator<<(std::ostream & stream, GhostType type) { switch (type) { case _not_ghost: stream << "not_ghost"; break; case _ghost: stream << "ghost"; break; case _casper: stream << "Casper the friendly ghost"; break; } return stream; } /* -------------------------------------------------------------------------- */ /// standard output stream operator for TimeStepSolverType inline std::ostream & operator<<(std::ostream & stream, const TimeStepSolverType & type) { switch (type) { case _tsst_static: stream << "static"; break; case _tsst_dynamic: stream << "dynamic"; break; case _tsst_dynamic_lumped: stream << "dynamic_lumped"; break; case _tsst_not_defined: stream << "not defined time step solver"; break; } return stream; } /* -------------------------------------------------------------------------- */ /// standard input stream operator for TimeStepSolverType inline std::istream & operator>>(std::istream & stream, TimeStepSolverType & type) { std::string str; stream >> str; if (str == "static") type = _tsst_static; else if (str == "dynamic") type = _tsst_dynamic; else if (str == "dynamic_lumped") type = _tsst_dynamic_lumped; else { AKANTU_ERROR("The type " << str << " is not a recognized TimeStepSolverType"); stream.setstate(std::ios::failbit); } return stream; } /* -------------------------------------------------------------------------- */ /// standard output stream operator for NonLinearSolverType inline std::ostream & operator<<(std::ostream & stream, const NonLinearSolverType & type) { switch (type) { case _nls_linear: stream << "linear"; break; case _nls_newton_raphson: stream << "newton_raphson"; break; case _nls_newton_raphson_modified: stream << "newton_raphson_modified"; break; case _nls_lumped: stream << "lumped"; break; case _nls_auto: stream << "auto"; break; } return stream; } /* -------------------------------------------------------------------------- */ /// standard input stream operator for NonLinearSolverType inline std::istream & operator>>(std::istream & stream, NonLinearSolverType & type) { std::string str; stream >> str; if (str == "linear") type = _nls_linear; else if (str == "newton_raphson") type = _nls_newton_raphson; else if (str == "newton_raphson_modified") type = _nls_newton_raphson_modified; else if (str == "lumped") type = _nls_lumped; else if (str == "auto") type = _nls_auto; else type = _nls_auto; return stream; } /* -------------------------------------------------------------------------- */ /// standard output stream operator for IntegrationSchemeType inline std::ostream & operator<<(std::ostream & stream, const IntegrationSchemeType & type) { switch (type) { case _ist_pseudo_time: stream << "pseudo_time"; break; case _ist_forward_euler: stream << "forward_euler"; break; case _ist_trapezoidal_rule_1: stream << "trapezoidal_rule_1"; break; case _ist_backward_euler: stream << "backward_euler"; break; case _ist_central_difference: stream << "central_difference"; break; case _ist_fox_goodwin: stream << "fox_goodwin"; break; case _ist_trapezoidal_rule_2: stream << "trapezoidal_rule_2"; break; case _ist_linear_acceleration: stream << "linear_acceleration"; break; case _ist_newmark_beta: stream << "newmark_beta"; break; case _ist_generalized_trapezoidal: stream << "generalized_trapezoidal"; break; } return stream; } /* -------------------------------------------------------------------------- */ /// standard input stream operator for IntegrationSchemeType inline std::istream & operator>>(std::istream & stream, IntegrationSchemeType & type) { std::string str; stream >> str; if (str == "pseudo_time") type = _ist_pseudo_time; else if (str == "forward_euler") type = _ist_forward_euler; else if (str == "trapezoidal_rule_1") type = _ist_trapezoidal_rule_1; else if (str == "backward_euler") type = _ist_backward_euler; else if (str == "central_difference") type = _ist_central_difference; else if (str == "fox_goodwin") type = _ist_fox_goodwin; else if (str == "trapezoidal_rule_2") type = _ist_trapezoidal_rule_2; else if (str == "linear_acceleration") type = _ist_linear_acceleration; else if (str == "newmark_beta") type = _ist_newmark_beta; else if (str == "generalized_trapezoidal") type = _ist_generalized_trapezoidal; else { AKANTU_ERROR("The type " << str << " is not a recognized IntegrationSchemeType"); stream.setstate(std::ios::failbit); } return stream; } /* -------------------------------------------------------------------------- */ /// standard output stream operator for SynchronizationTag inline std::ostream & operator<<(std::ostream & stream, SynchronizationTag type) { switch (type) { case _gst_whatever: stream << "_gst_whatever"; break; case _gst_ask_nodes: stream << "_gst_ask_nodes"; break; case _gst_update: stream << "_gst_update"; break; case _gst_size: stream << "_gst_size"; break; case _gst_smm_mass: stream << "_gst_smm_mass"; break; case _gst_smm_for_gradu: stream << "_gst_smm_for_gradu"; break; case _gst_smm_boundary: stream << "_gst_smm_boundary"; break; case _gst_smm_uv: stream << "_gst_smm_uv"; break; case _gst_smm_res: stream << "_gst_smm_res"; break; case _gst_smm_init_mat: stream << "_gst_smm_init_mat"; break; case _gst_smm_stress: stream << "_gst_smm_stress"; break; case _gst_smmc_facets: stream << "_gst_smmc_facets"; break; case _gst_smmc_facets_conn: stream << "_gst_smmc_facets_conn"; break; case _gst_smmc_facets_stress: stream << "_gst_smmc_facets_stress"; break; case _gst_smmc_damage: stream << "_gst_smmc_damage"; break; case _gst_giu_global_conn: stream << "_gst_giu_global_conn"; break; case _gst_ce_groups: stream << "_gst_ce_groups"; break; case _gst_gm_clusters: stream << "_gst_gm_clusters"; break; case _gst_htm_temperature: stream << "_gst_htm_temperature"; break; case _gst_htm_gradient_temperature: stream << "_gst_htm_gradient_temperature"; break; case _gst_htm_phi: stream << "_gst_htm_phi"; break; case _gst_htm_gradient_phi: stream << "_gst_htm_gradient_phi"; break; case _gst_mnl_for_average: stream << "_gst_mnl_for_average"; break; case _gst_mnl_weight: stream << "_gst_mnl_weight"; break; case _gst_nh_criterion: stream << "_gst_nh_criterion"; break; case _gst_test: stream << "_gst_test"; break; case _gst_user_1: stream << "_gst_user_1"; break; case _gst_user_2: stream << "_gst_user_2"; break; case _gst_material_id: stream << "_gst_material_id"; break; case _gst_for_dump: stream << "_gst_for_dump"; break; case _gst_cf_nodal: stream << "_gst_cf_nodal"; break; case _gst_cf_incr: stream << "_gst_cf_incr"; break; case _gst_solver_solution: stream << "_gst_solver_solution"; break; } return stream; } /* -------------------------------------------------------------------------- */ /// standard output stream operator for SolveConvergenceCriteria inline std::ostream & operator<<(std::ostream & stream, const SolveConvergenceCriteria & criteria) { switch (criteria) { case _scc_residual: stream << "_scc_residual"; break; case _scc_solution: stream << "_scc_solution"; break; case _scc_residual_mass_wgh: stream << "_scc_residual_mass_wgh"; break; } return stream; } inline std::istream & operator>>(std::istream & stream, SolveConvergenceCriteria & criteria) { std::string str; stream >> str; if (str == "residual") criteria = _scc_residual; else if (str == "solution") criteria = _scc_solution; else if (str == "residual_mass_wgh") criteria = _scc_residual_mass_wgh; else { stream.setstate(std::ios::failbit); } return stream; } /* -------------------------------------------------------------------------- */ inline std::string to_lower(const std::string & str) { std::string lstr = str; std::transform(lstr.begin(), lstr.end(), lstr.begin(), (int (*)(int))tolower); return lstr; } +namespace { + template + inline std::string trim_p(const std::string & to_trim, pred && p) { + std::string trimed = to_trim; + // left trim + trimed.erase(trimed.begin(), std::find_if(trimed.begin(), trimed.end(), std::not1(p))); + // right trim + trimed.erase(std::find_if(trimed.rbegin(), trimed.rend(), std::not1(p)).base(), + trimed.end()); + return trimed; + } + +} // namespace + /* -------------------------------------------------------------------------- */ inline std::string trim(const std::string & to_trim) { - std::string trimed = to_trim; - // left trim - trimed.erase(trimed.begin(), - std::find_if(trimed.begin(), trimed.end(), - std::not1(std::ptr_fun(isspace)))); - // right trim - trimed.erase(std::find_if(trimed.rbegin(), trimed.rend(), - std::not1(std::ptr_fun(isspace))) - .base(), - trimed.end()); - return trimed; + return trim_p(to_trim, [&](auto && a) { return std::isspace(a); }); +} + +inline std::string trim(const std::string & to_trim, char c) { + return trim_p(to_trim, [&c](auto && a) { return (a == c); }); } -} // akantu +} // namespace akantu #include namespace akantu { /* -------------------------------------------------------------------------- */ template std::string printMemorySize(UInt size) { Real real_size = size * sizeof(T); UInt mult = 0; if (real_size != 0) mult = (std::log(real_size) / std::log(2)) / 10; std::stringstream sstr; real_size /= Real(1 << (10 * mult)); sstr << std::setprecision(2) << std::fixed << real_size; std::string size_prefix; switch (mult) { case 0: sstr << ""; break; case 1: sstr << "Ki"; break; case 2: sstr << "Mi"; break; case 3: sstr << "Gi"; break; // I started on this type of machines // (32bit computers) (Nicolas) case 4: sstr << "Ti"; break; case 5: sstr << "Pi"; break; case 6: sstr << "Ei"; break; // theoritical limit of RAM of the current // computers in 2014 (64bit computers) (Nicolas) case 7: sstr << "Zi"; break; case 8: sstr << "Yi"; break; default: AKANTU_ERROR( "The programmer in 2014 didn't thought so far (even wikipedia does not " "go further)." << " You have at least 1024 times more than a yobibit of RAM!!!" << " Just add the prefix corresponding in this switch case."); } sstr << "Byte"; return sstr.str(); } -} // akantu +} // namespace akantu #endif /* __AKANTU_AKA_COMMON_INLINE_IMPL_CC__ */ diff --git a/src/io/mesh_io/mesh_io_msh.cc b/src/io/mesh_io/mesh_io_msh.cc index 728d32b15..793428718 100644 --- a/src/io/mesh_io/mesh_io_msh.cc +++ b/src/io/mesh_io/mesh_io_msh.cc @@ -1,1064 +1,1066 @@ /** * @file mesh_io_msh.cc * * @author Dana Christen * @author Mauro Corrado * @author David Simon Kammer * @author Nicolas Richart * * @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 . * */ /* ----------------------------------------------------------------------------- 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 /* -------------------------------------------------------------------------- */ #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 // // #include // #include // #include // #include // #include // #include // #include // #include // #include /* -------------------------------------------------------------------------- */ 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::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 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 tags; // std::vector connectivity; // ElementType type; // }; // /* ------------------------------------------------------------------------ // */ // struct lazy_get_nb_nodes_ { // template struct result { typedef int type; }; // template bool operator()(elem_type et) { // return MeshIOMSH::_msh_nodes_per_elem[et]; // } // }; // /* ------------------------------------------------------------------------ // */ // struct lazy_element_ { // template // struct result { // typedef _Element type; // }; // template // _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 struct result { typedef void type; }; // template 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 // struct result { // typedef bool type; // }; // template // bool operator()(Mesh & mesh, const ID & id, const V & pos, size max, // Map & nodes_mapping) const { // Vector 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 // struct result { // typedef bool type; // }; // template // bool operator()(Mesh & mesh, const EL & element, size max, // const NodeMap & nodes_mapping, // ElemMap & elements_mapping) const { // Vector 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 & connectivity = mesh.getConnectivity(element.type); // UInt i = 0; // for (std::vector::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 * data = // mesh.template getDataPointer(tag_sstr.str(), element.type, // _not_ghost); // data->push_back(element.tags[i]); // } // return (mesh.getNbElement() < UInt(max)); // } // }; // /* ------------------------------------------------------------------------ // */ // template // struct MshMeshGrammar : qi::grammar { // public: // MshMeshGrammar(Mesh & mesh) // : MshMeshGrammar::base_type(start, "msh_mesh_reader"), mesh(mesh) { // phx::function lazy_element; // phx::function lazy_get_nb_nodes; // phx::function lazy_check_value; // phx::function lazy_node_read; // phx::function 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 msh_element_type; // qi::symbols *> sections; // qi::rule start; // qi::rule > // unknown_section; // qi::rule *> > known_section; // qi::rule mesh_format, nodes, elements, // physical_names, ignore_section, // interpolation_scheme, element_data, node_data, any_line; // qi::rule > nodes_; // qi::rule, // vector > > elements_; // qi::rule(), Skipper> version; // qi::rule > // element; // qi::rule(), Skipper, qi::locals > tags; // qi::rule(int), Skipper> connectivity; // qi::rule(), Skipper> position; // qi::real_parser > 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 msh_nodes_to_akantu; // /// correspondance between the element number in the abaqus file and the // /// Element in the akantu mesh // std::map 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(infile), // std::istream_iterator(), // std::back_inserter(storage)); // typedef std::string::const_iterator iterator_t; // typedef ascii::space_type skipper; // typedef _parser::MshMeshGrammar 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::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> 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 & 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 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 * 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 & data = mesh.getDataPointer( 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 * data = &mesh.getDataPointer(tag_name, akantu_type, _not_ghost); data->push_back(tag); sstr_elem >> tag; // reg-elem tag_name = "tag_1"; data = &mesh.getDataPointer(tag_name, akantu_type, _not_ghost); data->push_back(tag); sstr_elem >> tag; // number-of-nodes } Vector 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); }; 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()); auto && real_tags[[gnu::unused]] = read_data_tags(std::vector()); auto && int_tags = read_data_tags(std::vector()); auto && data = mesh.registerNodalData(string_tags[0], int_tags[1]); data.resize(mesh.getNbNodes(), 0.); for (auto _[[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; + + data[node - first_node_number] = trim(value, \"); } }; 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("physical_names"); MeshUtils::fillElementToSubElementsData(mesh); } /* -------------------------------------------------------------------------- */ void MeshIOMSH::write(const std::string & filename, const Mesh & mesh) { std::ofstream outfile; const Array & 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::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 & 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("tag_0", type, _not_ghost)) { const auto & data_tag_0 = mesh.getData("tag_0", type, _not_ghost); tag[0] = data_tag_0.storage(); } if (mesh.hasData("tag_1", type, _not_ghost)) { const auto & data_tag_1 = mesh.getData("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(tag); outfile << "$NodeData" << "\n"; outfile << "1" << "\n"; - outfile << tag << "\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