diff --git a/examples/cohesive_element/cohesive_extrinsic_ig_tg/CMakeLists.txt b/examples/cohesive_element/cohesive_extrinsic_ig_tg/CMakeLists.txt index d48b8babb..44f6e456c 100644 --- a/examples/cohesive_element/cohesive_extrinsic_ig_tg/CMakeLists.txt +++ b/examples/cohesive_element/cohesive_extrinsic_ig_tg/CMakeLists.txt @@ -1,38 +1,37 @@ #=============================================================================== # @file CMakeLists.txt # # @author Seyedeh Mohadeseh Taheri Mousavi # # @date creation: Mon Jan 18 2016 # # @brief Configuration for cohesive_extrinsic_IG_TG # # @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 . # # @section DESCRIPTION # #=============================================================================== - +add_mesh(cohesive_extrinsic_ig_tg_mesh triangle.geo DIM 2 ORDER 1 OUTPUT square.msh) register_example(cohesive_extrinsic_ig_tg SOURCES cohesive_extrinsic_ig_tg.cc - FILES_TO_COPY material.dat square.msh + DEPENDS cohesive_extrinsic_ig_tg_mesh + FILES_TO_COPY material.dat ) - - diff --git a/examples/cohesive_element/cohesive_extrinsic_ig_tg/cohesive_extrinsic_ig_tg.cc b/examples/cohesive_element/cohesive_extrinsic_ig_tg/cohesive_extrinsic_ig_tg.cc index c609d4fe4..b49ff224a 100644 --- a/examples/cohesive_element/cohesive_extrinsic_ig_tg/cohesive_extrinsic_ig_tg.cc +++ b/examples/cohesive_element/cohesive_extrinsic_ig_tg/cohesive_extrinsic_ig_tg.cc @@ -1,195 +1,150 @@ /** * @file cohesive_extrinsic_ig_tg.cc * * @author Seyedeh Mohadeseh Taheri Mousavi * @author Marco Vocialta * * @date creation: Mon Jan 18 2016 * * @brief Test for considering different cohesive properties for intergranular * (IG) and * transgranular (TG) fractures in extrinsic cohesive elements * * @section LICENSE * * Copyright (©) 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 . * */ /* -------------------------------------------------------------------------- */ #include "solid_mechanics_model_cohesive.hh" /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ using namespace akantu; -class MultiGrainMaterialSelector : public DefaultMaterialCohesiveSelector { +/* -------------------------------------------------------------------------- */ +class Velocity : public BC::Dirichlet::DirichletFunctor { public: - MultiGrainMaterialSelector(const SolidMechanicsModelCohesive & model, - const ID & transgranular_id, - const ID & intergranular_id) - : DefaultMaterialCohesiveSelector(model), - transgranular_id(transgranular_id), intergranular_id(intergranular_id), - model(model), mesh(model.getMesh()), mesh_facets(model.getMeshFacets()), - spatial_dimension(model.getSpatialDimension()), nb_IG(0), nb_TG(0) {} - - UInt operator()(const Element & element) { - if (mesh_facets.getSpatialDimension(element.type) == - (spatial_dimension - 1)) { - const std::vector & element_to_subelement = - mesh_facets.getElementToSubelement(element.type, element.ghost_type)( - element.element); - - const Element & el1 = element_to_subelement[0]; - const Element & el2 = element_to_subelement[1]; - - UInt grain_id1 = - mesh.getData("tag_0", el1.type, el1.ghost_type)(el1.element); - if (el2 != ElementNull) { - UInt grain_id2 = - mesh.getData("tag_0", el2.type, el2.ghost_type)(el2.element); - if (grain_id1 == grain_id2) { - // transgranular = 0 indicator - nb_TG++; - return model.getMaterialIndex(transgranular_id); - } else { - // intergranular = 1 indicator - nb_IG++; - return model.getMaterialIndex(intergranular_id); - } - } else { - // transgranular = 0 indicator - nb_TG++; - return model.getMaterialIndex(transgranular_id); - } - } else { - return DefaultMaterialCohesiveSelector::operator()(element); - } + explicit Velocity(SolidMechanicsModel & model, Real vel, BC::Axis ax = _x) + : DirichletFunctor(ax), model(model), vel(vel) { + disp = vel * model.getTimeStep(); + } + +public: + inline void operator()(UInt node, Vector & /*flags*/, + Vector & disp, + const Vector & coord) const { + Real sign = std::signbit(coord(axis)) ? -1. : 1.; + disp(axis) += sign * this->disp; + model.getVelocity()(node, axis) = sign * vel; } private: - ID transgranular_id, intergranular_id; - const SolidMechanicsModelCohesive & model; - const Mesh & mesh; - const Mesh & mesh_facets; - UInt spatial_dimension; - - UInt nb_IG; - UInt nb_TG; + SolidMechanicsModel & model; + Real vel, disp; }; /* -------------------------------------------------------------------------- */ int main(int argc, char * argv[]) { initialize("material.dat", argc, argv); const UInt spatial_dimension = 2; const UInt max_steps = 1000; Mesh mesh(spatial_dimension); mesh.read("square.msh"); SolidMechanicsModelCohesive model(mesh); + MaterialCohesiveRules rules{{{"top", "bottom"}, "tg_cohesive"}, + {{"top", "top"}, "ig_cohesive"}, + {{"bottom", "bottom"}, "ig_cohesive"}}; /// model initialization - auto material_selector = std::make_shared( - model, "tg_cohesive", "ig_cohesive"); - + auto material_selector = + std::make_shared(model, rules); + auto && current_selector = model.getMaterialSelector(); + material_selector->setFallback(current_selector); + current_selector.setFallback( + std::make_shared>("physical_names", + model)); model.setMaterialSelector(material_selector); + model.initFull(_analysis_method = _explicit_lumped_mass, _is_extrinsic = true); Real time_step = model.getStableTimeStep() * 0.05; model.setTimeStep(time_step); std::cout << "Time step: " << time_step << std::endl; model.assembleMassLumped(); - Array & position = mesh.getNodes(); - Array & velocity = model.getVelocity(); - Array & boundary = model.getBlockedDOFs(); - Array & displacement = model.getDisplacement(); - - UInt nb_nodes = mesh.getNbNodes(); + auto & position = mesh.getNodes(); + auto & velocity = model.getVelocity(); - /// boundary conditions - for (UInt n = 0; n < nb_nodes; ++n) { - if (position(n, 1) > 0.99 || position(n, 1) < -0.99) - boundary(n, 1) = true; + model.applyBC(BC::Dirichlet::FlagOnly(_y), "top"); + model.applyBC(BC::Dirichlet::FlagOnly(_y), "bottom"); - if (position(n, 0) > 0.99 || position(n, 0) < -0.99) - boundary(n, 0) = true; - } + model.applyBC(BC::Dirichlet::FlagOnly(_x), "left"); + model.applyBC(BC::Dirichlet::FlagOnly(_x), "right"); model.setBaseName("extrinsic"); model.addDumpFieldVector("displacement"); model.addDumpField("velocity"); model.addDumpField("acceleration"); model.addDumpField("internal_force"); model.addDumpField("stress"); model.addDumpField("grad_u"); + model.addDumpField("material_index"); model.dump(); /// initial conditions Real loading_rate = 0.1; // bar_height = 2 Real VI = loading_rate * 2 * 0.5; - for (UInt n = 0; n < nb_nodes; ++n) { - velocity(n, 1) = loading_rate * position(n, 1); - velocity(n, 0) = loading_rate * position(n, 0); + for (auto && data : zip(make_view(position, spatial_dimension), + make_view(velocity, spatial_dimension))) { + std::get<1>(data) = loading_rate * std::get<0>(data); } model.dump(); - Real dispy = 0; + Velocity vely(model, VI, _y); + Velocity velx(model, VI, _x); /// Main loop for (UInt s = 1; s <= max_steps; ++s) { - dispy += VI * time_step; - /// update displacement on extreme nodes - for (UInt n = 0; n < mesh.getNbNodes(); ++n) { - if (position(n, 1) > 0.99) { - displacement(n, 1) = dispy; - velocity(n, 1) = VI; - } - if (position(n, 1) < -0.99) { - displacement(n, 1) = -dispy; - velocity(n, 1) = -VI; - } - if (position(n, 0) > 0.99) { - displacement(n, 0) = dispy; - velocity(n, 0) = VI; - } - if (position(n, 0) < -0.99) { - displacement(n, 0) = -dispy; - velocity(n, 0) = -VI; - } - } + + model.applyBC(vely, "top"); + model.applyBC(vely, "bottom"); + + model.applyBC(velx, "left"); + model.applyBC(velx, "right"); model.checkCohesiveStress(); model.solveStep(); if (s % 10 == 0) { model.dump(); std::cout << "passing step " << s << "/" << max_steps << std::endl; } } - finalize(); - return EXIT_SUCCESS; + return 0; } diff --git a/examples/cohesive_element/cohesive_extrinsic_ig_tg/material.dat b/examples/cohesive_element/cohesive_extrinsic_ig_tg/material.dat index 9220b0096..4373bc445 100644 --- a/examples/cohesive_element/cohesive_extrinsic_ig_tg/material.dat +++ b/examples/cohesive_element/cohesive_extrinsic_ig_tg/material.dat @@ -1,20 +1,33 @@ -material elastic [ - name = steel - rho = 1e3 # density - E = 1e3 # young's modulus - nu = 0.0 # poisson's ratio -] +model solid_mechanics_model_cohesive [ + cohesive_inserter [ + cohesive_zones = [btop, bbottom] + ] -material cohesive_linear [ - name = tg_cohesive - beta = 0 - G_c = 10 - sigma_c = 100 -] + material elastic [ + name = btop + rho = 1e3 # density + E = 1e3 # young's modulus + nu = 0.0 # poisson's ratio + ] + + material elastic [ + name = bbottom + rho = 1e3 # density + E = 1e3 # young's modulus + nu = 0.0 # poisson's ratio + ] + + material cohesive_linear [ + name = tg_cohesive + beta = 0 + G_c = 10 + sigma_c = 100 + ] -material cohesive_linear [ - name = ig_cohesive - beta = 0 - G_c = 10 - sigma_c = 20 + material cohesive_linear [ + name = ig_cohesive + beta = 0 + G_c = 10 + sigma_c = 20 + ] ] diff --git a/examples/cohesive_element/cohesive_extrinsic_ig_tg/square.msh b/examples/cohesive_element/cohesive_extrinsic_ig_tg/square.msh deleted file mode 100644 index e34221c8e..000000000 --- a/examples/cohesive_element/cohesive_extrinsic_ig_tg/square.msh +++ /dev/null @@ -1,66 +0,0 @@ -$MeshFormat -2.2 0 8 -$EndMeshFormat -$Nodes -41 -1 1 1 0 -2 -1 1 0 -3 -1 -1 0 -4 1 -1 0 -5 2.081668171172169e-12 1 0 -6 0.5000000000009475 1 0 -7 -0.4999999999989593 1 0 -8 -1 2.081668171172169e-12 0 -9 -1 0.5000000000009475 0 -10 -1 -0.4999999999989593 0 -11 -2.081668171172169e-12 -1 0 -12 -0.5000000000009475 -1 0 -13 0.4999999999989593 -1 0 -14 1 -2.081668171172169e-12 0 -15 1 -0.5000000000009475 0 -16 1 0.4999999999989593 0 -17 -2.888892828660043e-14 2.890280607440824e-14 0 -18 0.4999999999993061 -0.500000000000694 0 -19 0.500000000000694 0.4999999999993061 0 -20 -0.4999999999993061 0.500000000000694 0 -21 -0.5000000000006987 -0.4999999999993013 0 -22 0.7499999999996531 -0.7500000000003471 0 -23 0.2499999999986122 -0.7500000000003471 0 -24 0.7500000000003471 0.7499999999996531 0 -25 0.7500000000003471 0.2499999999986122 0 -26 -0.7499999999996531 0.7500000000003471 0 -27 -0.2499999999986122 0.7500000000003471 0 -28 -0.7500000000003493 -0.7499999999996506 0 -29 -0.7500000000003493 -0.2499999999986098 0 -30 1.026389621442784e-12 0.5000000000000144 0 -31 0.2500000000003326 0.2499999999996675 0 -32 0.2500000000013878 0.7499999999996531 0 -33 0.4999999999999856 -1.02638268254888e-12 0 -34 0.2499999999996386 -0.2500000000003325 0 -35 0.7499999999996531 -0.2500000000013878 0 -36 -0.2499999999996675 0.2500000000003614 0 -37 -0.5000000000000144 1.055285488623288e-12 0 -38 -0.7499999999996531 0.2500000000013878 0 -39 -1.055278549729385e-12 -0.4999999999999855 0 -40 -0.2500000000003638 -0.2499999999996362 0 -41 -0.2500000000013902 -0.7499999999996507 0 -$EndNodes -$Elements -16 -1 9 2 1 1 4 18 11 22 23 13 -2 9 2 2 2 1 19 14 24 25 16 -3 9 2 2 2 2 20 5 26 27 7 -4 9 2 1 1 3 21 8 28 29 10 -5 9 2 2 2 5 17 19 30 31 32 -6 9 2 1 1 14 17 18 33 34 35 -7 9 2 2 2 14 19 17 25 31 33 -8 9 2 2 2 5 20 17 27 36 30 -9 9 2 2 2 8 17 20 37 36 38 -10 9 2 1 1 11 18 17 23 34 39 -11 9 2 1 1 8 21 17 29 40 37 -12 9 2 1 1 11 17 21 39 40 41 -13 9 2 1 1 4 14 18 15 35 22 -14 9 2 2 2 2 8 20 9 38 26 -15 9 2 1 1 3 11 21 12 41 28 -16 9 2 2 2 1 5 19 6 32 24 -$EndElements diff --git a/examples/cohesive_element/cohesive_extrinsic_ig_tg/triangle.geo b/examples/cohesive_element/cohesive_extrinsic_ig_tg/triangle.geo new file mode 100644 index 000000000..374b8b4d9 --- /dev/null +++ b/examples/cohesive_element/cohesive_extrinsic_ig_tg/triangle.geo @@ -0,0 +1,31 @@ +h = 1.; + +Point(1) = { 1, 1, 0, h}; +Point(2) = {-1, 1, 0, h}; +Point(3) = {-1, -1, 0, h}; +Point(4) = { 1, -1, 0, h}; +Point(5) = {-1, 0, 0, h}; +Point(6) = { 1, 0, 0, h}; + +Line(1) = {3, 4}; +Line(2) = {4, 6}; +Line(3) = {6, 5}; +Line(4) = {5, 3}; +Line(5) = {6, 1}; +Line(6) = {1, 2}; +Line(7) = {2, 5}; + +Curve Loop(1) = {1, 2, 3, 4}; +Plane Surface(1) = {1}; + +Curve Loop(2) = {-3, 5, 6, 7}; +Plane Surface(2) = {2}; + +Physical Curve("top") = {6}; +Physical Curve("bottom") = {1}; + +Physical Curve("left") = {4, 7}; +Physical Curve("right") = {2, 5}; + +Physical Surface("btop") = {2}; +Physical Surface("bbottom") = {1}; diff --git a/src/io/parser/parameter_registry_tmpl.hh b/src/io/parser/parameter_registry_tmpl.hh index 24ebc5144..8a879f95d 100644 --- a/src/io/parser/parameter_registry_tmpl.hh +++ b/src/io/parser/parameter_registry_tmpl.hh @@ -1,410 +1,456 @@ /** * @file parameter_registry_tmpl.hh * * @author Nicolas Richart * * @date creation: Wed May 04 2016 * @date last modification: Tue Jan 30 2018 * * @brief implementation of the templated part of ParameterRegistry class and * the derivated ones * * @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 . * */ /* -------------------------------------------------------------------------- */ #include "aka_error.hh" #include "aka_iterators.hh" -#include "parameter_registry.hh" +//#include "parameter_registry.hh" #include "parser.hh" /* -------------------------------------------------------------------------- */ #include +#include #include #include /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_PARAMETER_REGISTRY_TMPL_HH__ #define __AKANTU_PARAMETER_REGISTRY_TMPL_HH__ namespace akantu { namespace debug { class ParameterException : public Exception { public: ParameterException(const std::string & name, const std::string & message) : Exception(message), name(name) {} const std::string & name; }; class ParameterUnexistingException : public ParameterException { public: ParameterUnexistingException(const std::string & name, const ParameterRegistry & registery) : ParameterException(name, "Parameter " + name + " does not exists in this scope") { auto && params = registery.listParameters(); this->_info = std::accumulate(params.begin(), params.end(), this->_info + "\n Possible parameters are: ", [](auto && str, auto && param) { static auto first = true; auto ret = str + (first ? " " : ", ") + param; first = false; return ret; }); } }; class ParameterAccessRightException : public ParameterException { public: ParameterAccessRightException(const std::string & name, const std::string & perm) : ParameterException(name, "Parameter " + name + " is not " + perm) {} }; class ParameterWrongTypeException : public ParameterException { public: ParameterWrongTypeException(const std::string & name, const std::type_info & wrong_type, const std::type_info & type) : ParameterException(name, "Parameter " + name + " type error, cannot convert " + debug::demangle(type.name()) + " to " + debug::demangle(wrong_type.name())) {} }; } // namespace debug /* -------------------------------------------------------------------------- */ template const ParameterTyped & Parameter::getParameterTyped() const { try { const auto & tmp = aka::as_type>(*this); return tmp; } catch (std::bad_cast &) { AKANTU_CUSTOM_EXCEPTION( debug::ParameterWrongTypeException(name, typeid(T), this->type())); } } /* -------------------------------------------------------------------------- */ template ParameterTyped & Parameter::getParameterTyped() { try { auto & tmp = aka::as_type>(*this); return tmp; } catch (std::bad_cast &) { AKANTU_CUSTOM_EXCEPTION( debug::ParameterWrongTypeException(name, typeid(T), this->type())); } } /* ------------------------------------------------------------------------ */ template void Parameter::set(const V & value) { if (!(isWritable())) AKANTU_CUSTOM_EXCEPTION( debug::ParameterAccessRightException(name, "writable")); ParameterTyped & typed_param = getParameterTyped(); typed_param.setTyped(value); } /* ------------------------------------------------------------------------ */ inline void Parameter::setAuto(__attribute__((unused)) const ParserParameter & value) { if (!(isParsable())) AKANTU_CUSTOM_EXCEPTION( debug::ParameterAccessRightException(name, "parsable")); } /* -------------------------------------------------------------------------- */ template const T & Parameter::get() const { if (!(isReadable())) AKANTU_CUSTOM_EXCEPTION( debug::ParameterAccessRightException(name, "readable")); const ParameterTyped & typed_param = getParameterTyped(); return typed_param.getTyped(); } /* -------------------------------------------------------------------------- */ template T & Parameter::get() { ParameterTyped & typed_param = getParameterTyped(); if (!(isReadable()) || !(this->isWritable())) AKANTU_CUSTOM_EXCEPTION( debug::ParameterAccessRightException(name, "accessible")); return typed_param.getTyped(); } /* -------------------------------------------------------------------------- */ template inline Parameter::operator T() const { return this->get(); } /* -------------------------------------------------------------------------- */ template ParameterTyped::ParameterTyped(std::string name, std::string description, ParameterAccessType param_type, T & param) : Parameter(name, description, param_type), param(param) {} /* -------------------------------------------------------------------------- */ template template void ParameterTyped::setTyped(const V & value) { param = value; } /* -------------------------------------------------------------------------- */ template inline void ParameterTyped::setAuto(const ParserParameter & value) { Parameter::setAuto(value); param = static_cast(value); } /* -------------------------------------------------------------------------- */ template <> inline void ParameterTyped::setAuto(const ParserParameter & value) { Parameter::setAuto(value); param = value.getValue(); } /* -------------------------------------------------------------------------- */ template <> inline void ParameterTyped>::setAuto(const ParserParameter & in_param) { Parameter::setAuto(in_param); Vector tmp = in_param; if (param.size() == 0) { param = tmp; } else { for (UInt i = 0; i < param.size(); ++i) { param(i) = tmp(i); } } } /* -------------------------------------------------------------------------- */ template <> inline void ParameterTyped>::setAuto(const ParserParameter & in_param) { Parameter::setAuto(in_param); Matrix tmp = in_param; if (param.size() == 0) { param = tmp; } else { for (UInt i = 0; i < param.rows(); ++i) { for (UInt j = 0; j < param.cols(); ++j) { param(i, j) = tmp(i, j); } } } } /* -------------------------------------------------------------------------- */ template const T & ParameterTyped::getTyped() const { return param; } /* -------------------------------------------------------------------------- */ template T & ParameterTyped::getTyped() { return param; } /* -------------------------------------------------------------------------- */ template inline void ParameterTyped::printself(std::ostream & stream) const { Parameter::printself(stream); stream << param << "\n"; } /* -------------------------------------------------------------------------- */ template class ParameterTyped> : public Parameter { public: ParameterTyped(std::string name, std::string description, ParameterAccessType param_type, std::vector & param) : Parameter(name, description, param_type), param(param) {} - /* ------------------------------------------------------------------------ */ - template void setTyped(const V & value) { param = value; } + /* ------------------------------------------------------------------------ + */ + template + void setTyped(const V & value) { + param = value; + } void setAuto(const ParserParameter & value) override { Parameter::setAuto(value); param.clear(); const std::vector & tmp = value; for (auto && z : tmp) { param.emplace_back(z); } } std::vector & getTyped() { return param; } const std::vector & getTyped() const { return param; } void printself(std::ostream & stream) const override { Parameter::printself(stream); stream << "[ "; for (auto && v : param) stream << v << " "; stream << "]\n"; } inline const std::type_info & type() const override { return typeid(std::vector); } private: /// Value of parameter std::vector & param; }; -/* ------o-------------------------------------------------------------------- - */ +/* -------------------------------------------------------------------------- */ +template class ParameterTyped> : public Parameter { +public: + ParameterTyped(std::string name, std::string description, + ParameterAccessType param_type, std::set & param) + : Parameter(name, description, param_type), param(param) {} + + /* ------------------------------------------------------------------------ + */ + template + void setTyped(const V & value) { + param = value; + } + void setAuto(const ParserParameter & value) override { + Parameter::setAuto(value); + param.clear(); + const std::set & tmp = value; + for (auto && z : tmp) { + param.emplace(z); + } + } + + std::set & getTyped() { return param; } + const std::set & getTyped() const { return param; } + + void printself(std::ostream & stream) const override { + Parameter::printself(stream); + stream << "[ "; + for (auto && v : param) + stream << v << " "; + stream << "]\n"; + } + + inline const std::type_info & type() const override { + return typeid(std::set); + } + +private: + /// Value of parameter + std::set & param; +}; + +/* -------------------------------------------------------------------------- */ template <> inline void ParameterTyped::printself(std::ostream & stream) const { Parameter::printself(stream); stream << std::boolalpha << param << "\n"; } /* -------------------------------------------------------------------------- */ template void ParameterRegistry::registerParam(std::string name, T & variable, ParameterAccessType type, const std::string & description) { auto it = params.find(name); if (it != params.end()) AKANTU_CUSTOM_EXCEPTION(debug::ParameterException( name, "Parameter named " + name + " already registered.")); auto * param = new ParameterTyped(name, description, type, variable); params[name] = param; } /* -------------------------------------------------------------------------- */ template void ParameterRegistry::registerParam(std::string name, T & variable, const T & default_value, ParameterAccessType type, const std::string & description) { variable = default_value; registerParam(name, variable, type, description); } /* -------------------------------------------------------------------------- */ template void ParameterRegistry::setMixed(const std::string & name, const V & value) { auto it = params.find(name); if (it == params.end()) { if (consisder_sub) { for (auto it = sub_registries.begin(); it != sub_registries.end(); ++it) { it->second->setMixed(name, value); } } else { AKANTU_CUSTOM_EXCEPTION(debug::ParameterUnexistingException(name, *this)); } } else { Parameter & param = *(it->second); param.set(value); } } /* -------------------------------------------------------------------------- */ template void ParameterRegistry::set(const std::string & name, const T & value) { this->template setMixed(name, value); } /* -------------------------------------------------------------------------- */ template T & ParameterRegistry::get_(const std::string & name) { auto it = params.find(name); if (it == params.end()) { if (consisder_sub) { for (auto it = sub_registries.begin(); it != sub_registries.end(); ++it) { try { return it->second->get_(name); } catch (...) { } } } // nothing was found not even in sub registries AKANTU_CUSTOM_EXCEPTION(debug::ParameterUnexistingException(name, *this)); } Parameter & param = *(it->second); return param.get(); } /* -------------------------------------------------------------------------- */ const Parameter & ParameterRegistry::get(const std::string & name) const { auto it = params.find(name); if (it == params.end()) { if (consisder_sub) { for (auto it = sub_registries.begin(); it != sub_registries.end(); ++it) { try { return it->second->get(name); } catch (...) { } } } // nothing was found not even in sub registries AKANTU_CUSTOM_EXCEPTION(debug::ParameterUnexistingException(name, *this)); } Parameter & param = *(it->second); return param; } /* -------------------------------------------------------------------------- */ Parameter & ParameterRegistry::get(const std::string & name) { auto it = params.find(name); if (it == params.end()) { if (consisder_sub) { for (auto it = sub_registries.begin(); it != sub_registries.end(); ++it) { try { return it->second->get(name); } catch (...) { } } } // nothing was found not even in sub registries AKANTU_CUSTOM_EXCEPTION(debug::ParameterUnexistingException(name, *this)); } Parameter & param = *(it->second); return param; } /* -------------------------------------------------------------------------- */ namespace { namespace details { template struct CastHelper { static R convert(const T &) { throw std::bad_cast(); } }; template struct CastHelper::value>> { static R convert(const T & val) { return val; } }; } // namespace details } // namespace template inline ParameterTyped::operator Real() const { if (not isReadable()) AKANTU_CUSTOM_EXCEPTION( debug::ParameterAccessRightException(name, "accessible")); return details::CastHelper::convert(param); } } // namespace akantu #endif /* __AKANTU_PARAMETER_REGISTRY_TMPL_HH__ */ diff --git a/src/io/parser/parser_tmpl.hh b/src/io/parser/parser_tmpl.hh index 83644f7a5..fb45a2da9 100644 --- a/src/io/parser/parser_tmpl.hh +++ b/src/io/parser/parser_tmpl.hh @@ -1,112 +1,124 @@ /** * @file parser_tmpl.hh * * @author Nicolas Richart * * @date creation: Wed Nov 13 2013 * @date last modification: Mon Dec 18 2017 * * @brief Implementation of the parser templated methods * * @section LICENSE * * Copyright (©) 2014-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 . * */ /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ namespace akantu { /* -------------------------------------------------------------------------- */ template inline ParserParameter::operator T() const { T t; std::stringstream sstr(value); sstr >> t; if (sstr.bad()) AKANTU_EXCEPTION("No known conversion of a ParserParameter \"" << name << "\" to the type " << typeid(T).name()); return t; } /* -------------------------------------------------------------------------- */ template <> inline ParserParameter::operator const char *() const { return value.c_str(); } /* -------------------------------------------------------------------------- */ template <> inline ParserParameter::operator Real() const { return Parser::parseReal(value, *parent_section); } /* --------------------------------------------------------- ----------------- */ template <> inline ParserParameter::operator bool() const { bool b; std::stringstream sstr(value); sstr >> std::boolalpha >> b; if (sstr.fail()) { sstr.clear(); sstr >> std::noboolalpha >> b; } return b; } /* -------------------------------------------------------------------------- */ template <> inline ParserParameter::operator std::vector() const { std::vector tmp; auto string = std::regex_replace(value, std::regex("[[:space:]]|\\[|\\]"), ""); std::smatch sm; while (std::regex_search(string, sm, std::regex("[^,]+"))) { tmp.push_back(sm.str()); string = sm.suffix(); } return tmp; } -/* --------------------------------------------------------- ----------------- - */ +/* -------------------------------------------------------------------------- */ +template <> inline ParserParameter::operator std::set() const { + std::set tmp; + auto string = + std::regex_replace(value, std::regex("[[:space:]]|\\[|\\]"), ""); + std::smatch sm; + while (std::regex_search(string, sm, std::regex("[^,]+"))) { + tmp.emplace(sm.str()); + string = sm.suffix(); + } + return tmp; +} + +/* -------------------------------------------------------------------------- */ template <> inline ParserParameter::operator Vector() const { return Parser::parseVector(value, *parent_section); } /* --------------------------------------------------------- ----------------- */ template <> inline ParserParameter::operator Vector() const { Vector tmp = Parser::parseVector(value, *parent_section); Vector tmp_uint(tmp.size()); for (UInt i = 0; i < tmp.size(); ++i) { tmp_uint(i) = UInt(tmp(i)); } return tmp_uint; } /* --------------------------------------------------------- ----------------- */ template <> inline ParserParameter::operator Matrix() const { return Parser::parseMatrix(value, *parent_section); } /* -------------------------------------------------------------------------- */ template <> inline ParserParameter::operator RandomParameter() const { return Parser::parseRandomParameter(value, *parent_section); } } // namespace akantu diff --git a/src/mesh_utils/cohesive_element_inserter.cc b/src/mesh_utils/cohesive_element_inserter.cc index e09aaf7c5..3c5f4a38e 100644 --- a/src/mesh_utils/cohesive_element_inserter.cc +++ b/src/mesh_utils/cohesive_element_inserter.cc @@ -1,281 +1,308 @@ /** * @file cohesive_element_inserter.cc * * @author Marco Vocialta * * @date creation: Wed Dec 04 2013 * @date last modification: Mon Feb 19 2018 * * @brief Cohesive element inserter functions * * @section LICENSE * * Copyright (©) 2014-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 . * */ /* -------------------------------------------------------------------------- */ #include "cohesive_element_inserter.hh" #include "communicator.hh" #include "element_group.hh" #include "element_synchronizer.hh" #include "global_ids_updater.hh" #include "mesh_accessor.hh" #include "mesh_iterators.hh" /* -------------------------------------------------------------------------- */ #include #include /* -------------------------------------------------------------------------- */ namespace akantu { CohesiveElementInserter::CohesiveElementInserter(Mesh & mesh, const ID & id) : Parsable(ParserType::_cohesive_inserter), id(id), mesh(mesh), mesh_facets(mesh.initMeshFacets()), insertion_facets("insertion_facets", id), insertion_limits(mesh.getSpatialDimension(), 2), check_facets("check_facets", id) { - this->registerParam("cohesive_surfaces", physical_groups, _pat_parsable, + this->registerParam("cohesive_surfaces", physical_surfaces, _pat_parsable, + "List of groups to consider for insertion"); + this->registerParam("cohesive_zones", physical_zones, _pat_parsable, "List of groups to consider for insertion"); this->registerParam("bounding_box", insertion_limits, _pat_parsable, "Global limit for insertion"); UInt spatial_dimension = mesh.getSpatialDimension(); MeshUtils::resetFacetToDouble(mesh_facets); /// init insertion limits for (UInt dim = 0; dim < spatial_dimension; ++dim) { insertion_limits(dim, 0) = std::numeric_limits::max() * Real(-1.); insertion_limits(dim, 1) = std::numeric_limits::max(); } insertion_facets.initialize(mesh_facets, _spatial_dimension = spatial_dimension - 1, _with_nb_element = true, _default_value = false); } /* -------------------------------------------------------------------------- */ CohesiveElementInserter::~CohesiveElementInserter() = default; /* -------------------------------------------------------------------------- */ void CohesiveElementInserter::parseSection(const ParserSection & section) { Parsable::parseSection(section); if (is_extrinsic) limitCheckFacets(this->check_facets); } /* -------------------------------------------------------------------------- */ void CohesiveElementInserter::limitCheckFacets() { limitCheckFacets(this->check_facets); } /* -------------------------------------------------------------------------- */ void CohesiveElementInserter::setLimit(SpatialDirection axis, Real first_limit, Real second_limit) { AKANTU_DEBUG_ASSERT( axis < mesh.getSpatialDimension(), "You are trying to limit insertion in a direction that doesn't exist"); insertion_limits(axis, 0) = std::min(first_limit, second_limit); insertion_limits(axis, 1) = std::max(first_limit, second_limit); } /* -------------------------------------------------------------------------- */ UInt CohesiveElementInserter::insertIntrinsicElements() { limitCheckFacets(insertion_facets); return insertElements(); } /* -------------------------------------------------------------------------- */ void CohesiveElementInserter::limitCheckFacets( ElementTypeMapArray & check_facets) { AKANTU_DEBUG_IN(); UInt spatial_dimension = mesh.getSpatialDimension(); check_facets.initialize(mesh_facets, _spatial_dimension = spatial_dimension - 1, _with_nb_element = true, _default_value = true); check_facets.set(true); // remove the pure ghost elements for_each_element( mesh_facets, [&](auto && facet) { const auto & element_to_facet = mesh_facets.getElementToSubelement( facet.type, facet.ghost_type)(facet.element); auto & left = element_to_facet[0]; auto & right = element_to_facet[1]; if (right == ElementNull || (left.ghost_type == _ghost && right.ghost_type == _ghost)) { check_facets(facet) = false; return; } #ifndef AKANTU_NDEBUG if (left == ElementNull) { AKANTU_DEBUG_WARNING("By convention element should not have " "ElementNull on there first side: " << facet); } #endif if (left.kind() == _ek_cohesive or right.kind() == _ek_cohesive) { check_facets(facet) = false; } }, _spatial_dimension = spatial_dimension - 1); auto tolerance = Math::getTolerance(); Vector bary_facet(spatial_dimension); // set the limits to the bounding box for_each_element( mesh_facets, [&](auto && facet) { auto & need_check = check_facets(facet); if (not need_check) return; mesh_facets.getBarycenter(facet, bary_facet); UInt coord_in_limit = 0; while (coord_in_limit < spatial_dimension and bary_facet(coord_in_limit) > (insertion_limits(coord_in_limit, 0) - tolerance) and bary_facet(coord_in_limit) < (insertion_limits(coord_in_limit, 1) + tolerance)) ++coord_in_limit; if (coord_in_limit != spatial_dimension) need_check = false; }, _spatial_dimension = spatial_dimension - 1); - if (physical_groups.size() == 0) { + // remove the physical zones + if(mesh.hasData("physical_names") and physical_zones.size() > 0) { + auto && physical_names = mesh.getData("physical_names"); + for_each_element( + mesh_facets, + [&](auto && facet) { + const auto & element_to_facet = mesh_facets.getElementToSubelement( + facet.type, facet.ghost_type)(facet.element); + auto count = 0; + for(auto i : arange(2)) { + const auto & element = element_to_facet[i]; + if(element == ElementNull) + continue; + const auto & name = physical_names(element); + count += find(physical_zones.begin(), + physical_zones.end(), name) != physical_zones.end(); + + } + + if(count != 2) + check_facets(facet) = false; + }, + _spatial_dimension = spatial_dimension - 1); + } + + if (physical_surfaces.size() == 0) { AKANTU_DEBUG_OUT(); return; } if (not mesh_facets.hasData("physical_names")) { AKANTU_DEBUG_ASSERT( - physical_groups.size() == 0, + physical_surfaces.size() == 0, "No physical names in the mesh but insertion limited to a group"); AKANTU_DEBUG_OUT(); return; } const auto & physical_ids = mesh_facets.getData("physical_names"); // set the limits to the physical surfaces for_each_element(mesh_facets, [&](auto && facet) { auto & need_check = check_facets(facet); if (not need_check) return; const auto & physical_id = physical_ids(facet); - auto it = find(physical_groups.begin(), - physical_groups.end(), physical_id); + auto it = find(physical_surfaces.begin(), + physical_surfaces.end(), physical_id); - need_check = (it != physical_groups.end()); + need_check = (it != physical_surfaces.end()); }, _spatial_dimension = spatial_dimension - 1); AKANTU_DEBUG_OUT(); } /* -------------------------------------------------------------------------- */ UInt CohesiveElementInserter::insertElements(bool only_double_facets) { CohesiveNewNodesEvent node_event(AKANTU_CURRENT_FUNCTION); NewElementsEvent element_event(AKANTU_CURRENT_FUNCTION); Array new_pairs(0, 2); if (mesh_facets.isDistributed()) { mesh_facets.getElementSynchronizer().synchronizeOnce( *this, SynchronizationTag::_ce_groups); } UInt nb_new_elements = MeshUtils::insertCohesiveElements( mesh, mesh_facets, insertion_facets, new_pairs, element_event.getList(), only_double_facets); UInt nb_new_nodes = new_pairs.size(); node_event.getList().reserve(nb_new_nodes); node_event.getOldNodesList().reserve(nb_new_nodes); for (UInt n = 0; n < nb_new_nodes; ++n) { node_event.getList().push_back(new_pairs(n, 1)); node_event.getOldNodesList().push_back(new_pairs(n, 0)); } if (nb_new_elements > 0) { updateInsertionFacets(); } MeshAccessor mesh_accessor(mesh); std::tie(nb_new_nodes, nb_new_elements) = mesh_accessor.updateGlobalData(node_event, element_event); return nb_new_elements; } /* -------------------------------------------------------------------------- */ void CohesiveElementInserter::updateInsertionFacets() { AKANTU_DEBUG_IN(); UInt spatial_dimension = mesh.getSpatialDimension(); for (auto && facet_gt : ghost_types) { for (auto && facet_type : mesh_facets.elementTypes(spatial_dimension - 1, facet_gt)) { auto & ins_facets = insertion_facets(facet_type, facet_gt); // this is the intrinsic case if (not is_extrinsic) continue; auto & f_check = check_facets(facet_type, facet_gt); for (auto && pair : zip(ins_facets, f_check)) { bool & ins = std::get<0>(pair); bool & check = std::get<1>(pair); if (ins) ins = check = false; } } } // resize for the newly added facets insertion_facets.initialize(mesh_facets, _spatial_dimension = spatial_dimension - 1, _with_nb_element = true, _default_value = false); // resize for the newly added facets if (is_extrinsic) { check_facets.initialize(mesh_facets, _spatial_dimension = spatial_dimension - 1, _with_nb_element = true, _default_value = false); } else { insertion_facets.set(false); } AKANTU_DEBUG_OUT(); } } // namespace akantu diff --git a/src/mesh_utils/cohesive_element_inserter.hh b/src/mesh_utils/cohesive_element_inserter.hh index 9c767dab1..8551193de 100644 --- a/src/mesh_utils/cohesive_element_inserter.hh +++ b/src/mesh_utils/cohesive_element_inserter.hh @@ -1,171 +1,174 @@ /** * @file cohesive_element_inserter.hh * * @author Fabian Barras * @author Nicolas Richart * @author Marco Vocialta * * @date creation: Wed Dec 04 2013 * @date last modification: Sun Feb 04 2018 * * @brief Cohesive element inserter * * @section LICENSE * * Copyright (©) 2014-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 . * */ /* -------------------------------------------------------------------------- */ #include "data_accessor.hh" #include "mesh_utils.hh" #include "parsable.hh" /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_COHESIVE_ELEMENT_INSERTER_HH__ #define __AKANTU_COHESIVE_ELEMENT_INSERTER_HH__ namespace akantu { class GlobalIdsUpdater; class FacetSynchronizer; } // namespace akantu namespace akantu { class CohesiveElementInserter : public DataAccessor, public Parsable { /* ------------------------------------------------------------------------ */ /* Constructors/Destructors */ /* ------------------------------------------------------------------------ */ public: CohesiveElementInserter(Mesh & mesh, const ID & id = "cohesive_element_inserter"); ~CohesiveElementInserter() override; /* ------------------------------------------------------------------------ */ /* Methods */ /* ------------------------------------------------------------------------ */ public: /// set range limitation for intrinsic cohesive element insertion void setLimit(SpatialDirection axis, Real first_limit, Real second_limit); /// insert intrinsic cohesive elements in a predefined range UInt insertIntrinsicElements(); /// insert extrinsic cohesive elements (returns the number of new /// cohesive elements) UInt insertElements(bool only_double_facets = false); /// limit check facets to match given insertion limits void limitCheckFacets(); protected: void parseSection(const ParserSection & section) override; protected: /// internal version of limitCheckFacets void limitCheckFacets(ElementTypeMapArray & check_facets); /// update facet insertion arrays after facets doubling void updateInsertionFacets(); /// functions for parallel communications inline UInt getNbData(const Array & elements, const SynchronizationTag & tag) const override; inline void packData(CommunicationBuffer & buffer, const Array & elements, const SynchronizationTag & tag) const override; inline void unpackData(CommunicationBuffer & buffer, const Array & elements, const SynchronizationTag & tag) override; /* ------------------------------------------------------------------------ */ /* Accessors */ /* ------------------------------------------------------------------------ */ public: AKANTU_GET_MACRO_NOT_CONST(InsertionFacetsByElement, insertion_facets, ElementTypeMapArray &); AKANTU_GET_MACRO(InsertionFacetsByElement, insertion_facets, const ElementTypeMapArray &); AKANTU_GET_MACRO_BY_ELEMENT_TYPE(InsertionFacets, insertion_facets, bool); AKANTU_GET_MACRO(CheckFacets, check_facets, const ElementTypeMapArray &); AKANTU_GET_MACRO_BY_ELEMENT_TYPE(CheckFacets, check_facets, bool); AKANTU_GET_MACRO_BY_ELEMENT_TYPE_CONST(CheckFacets, check_facets, bool); AKANTU_GET_MACRO(MeshFacets, mesh_facets, const Mesh &); AKANTU_GET_MACRO_NOT_CONST(MeshFacets, mesh_facets, Mesh &); AKANTU_SET_MACRO(IsExtrinsic, is_extrinsic, bool); public: friend class SolidMechanicsModelCohesive; /* ------------------------------------------------------------------------ */ /* Class Members */ /* ------------------------------------------------------------------------ */ private: /// object id ID id; /// main mesh where to insert cohesive elements Mesh & mesh; /// mesh containing facets Mesh & mesh_facets; /// list of facets where to insert elements ElementTypeMapArray insertion_facets; /// limits for element insertion Matrix insertion_limits; /// list of groups to consider for insertion, ignored if empty - std::vector physical_groups; + std::set physical_surfaces; + + /// list of groups in between which an inside which element are insterted + std::set physical_zones; /// vector containing facets in which extrinsic cohesive elements can be /// inserted ElementTypeMapArray check_facets; /// global connectivity ids updater std::unique_ptr global_ids_updater; /// is this inserter used in extrinsic bool is_extrinsic{false}; }; class CohesiveNewNodesEvent : public NewNodesEvent { public: CohesiveNewNodesEvent(const std::string & origin) : NewNodesEvent(origin) {} ~CohesiveNewNodesEvent() override = default; AKANTU_GET_MACRO_NOT_CONST(OldNodesList, old_nodes, Array &); AKANTU_GET_MACRO(OldNodesList, old_nodes, const Array &); private: Array old_nodes; }; } // namespace akantu #include "cohesive_element_inserter_inline_impl.cc" #endif /* __AKANTU_COHESIVE_ELEMENT_INSERTER_HH__ */