diff --git a/src/common/aka_grid_dynamic.hh b/src/common/aka_grid_dynamic.hh index d717b12c8..9d7550fc4 100644 --- a/src/common/aka_grid_dynamic.hh +++ b/src/common/aka_grid_dynamic.hh @@ -1,512 +1,513 @@ /** * @file aka_grid_dynamic.hh * * @author Aurelia Isabel Cuba Ramos * @author Nicolas Richart * * @date creation: Thu Feb 21 2013 * @date last modification: Wed Nov 08 2017 * * @brief Grid that is auto balanced * * @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 "aka_array.hh" #include "aka_common.hh" #include "aka_types.hh" #include /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_AKA_GRID_DYNAMIC_HH__ #define __AKANTU_AKA_GRID_DYNAMIC_HH__ namespace akantu { class Mesh; template class SpatialGrid { public: explicit SpatialGrid(UInt dimension) : dimension(dimension), spacing(dimension), center(dimension), lower(dimension), upper(dimension), empty_cell() {} SpatialGrid(UInt dimension, const Vector & spacing, const Vector & center) : dimension(dimension), spacing(spacing), center(center), lower(dimension), upper(dimension), empty_cell() { for (UInt i = 0; i < dimension; ++i) { lower(i) = std::numeric_limits::max(); upper(i) = -std::numeric_limits::max(); } } virtual ~SpatialGrid() = default; class neighbor_cells_iterator; class cells_iterator; class CellID { public: CellID() : ids() {} explicit CellID(UInt dimention) : ids(dimention) {} void setID(UInt dir, Int id) { ids(dir) = id; } Int getID(UInt dir) const { return ids(dir); } bool operator<(const CellID & id) const { return std::lexicographical_compare( ids.storage(), ids.storage() + ids.size(), id.ids.storage(), id.ids.storage() + id.ids.size()); } bool operator==(const CellID & id) const { return std::equal(ids.storage(), ids.storage() + ids.size(), id.ids.storage()); } bool operator!=(const CellID & id) const { return !(operator==(id)); } class neighbor_cells_iterator : private std::iterator { public: neighbor_cells_iterator(const CellID & cell_id, bool end) : cell_id(cell_id), position(cell_id.ids.size(), end ? 1 : -1) { this->updateIt(); if (end) this->it++; } neighbor_cells_iterator & operator++() { UInt i = 0; for (; i < position.size() && position(i) == 1; ++i) ; if (i == position.size()) { ++it; return *this; } for (UInt j = 0; j < i; ++j) position(j) = -1; position(i)++; updateIt(); return *this; } neighbor_cells_iterator operator++(int) { neighbor_cells_iterator tmp(*this); operator++(); return tmp; }; bool operator==(const neighbor_cells_iterator & rhs) const { return cell_id == rhs.cell_id && it == rhs.it; }; bool operator!=(const neighbor_cells_iterator & rhs) const { return !operator==(rhs); }; CellID operator*() const { CellID cur_cell_id(cell_id); cur_cell_id.ids += position; return cur_cell_id; }; private: void updateIt() { it = 0; for (UInt i = 0; i < position.size(); ++i) it = it * 3 + (position(i) + 1); } private: /// central cell id const CellID & cell_id; // number representing the current neighbor in base 3; UInt it; // current cell shift Vector position; }; class Neighbors { public: explicit Neighbors(const CellID & cell_id) : cell_id(cell_id) {} decltype(auto) begin() { return neighbor_cells_iterator(cell_id, false); } decltype(auto) end() { return neighbor_cells_iterator(cell_id, true); } private: const CellID & cell_id; }; decltype(auto) neighbors() { return Neighbors(*this); } private: friend class cells_iterator; Vector ids; }; /* ------------------------------------------------------------------------ */ class Cell { public: using iterator = typename std::vector::iterator; using const_iterator = typename std::vector::const_iterator; Cell() : id(), data() {} explicit Cell(const CellID & cell_id) : id(cell_id), data() {} bool operator==(const Cell & cell) const { return id == cell.id; } bool operator!=(const Cell & cell) const { return id != cell.id; } Cell & add(const T & d) { data.push_back(d); return *this; } iterator begin() { return data.begin(); } const_iterator begin() const { return data.begin(); } iterator end() { return data.end(); } const_iterator end() const { return data.end(); } private: CellID id; std::vector data; }; private: using cells_container = std::map; public: const Cell & getCell(const CellID & cell_id) const { auto it = cells.find(cell_id); if (it != cells.end()) return it->second; return empty_cell; } decltype(auto) beginCell(const CellID & cell_id) { auto it = cells.find(cell_id); if (it != cells.end()) return it->second.begin(); return empty_cell.begin(); } decltype(auto) endCell(const CellID & cell_id) { auto it = cells.find(cell_id); if (it != cells.end()) return it->second.end(); return empty_cell.end(); } decltype(auto) beginCell(const CellID & cell_id) const { auto it = cells.find(cell_id); if (it != cells.end()) return it->second.begin(); return empty_cell.begin(); } decltype(auto) endCell(const CellID & cell_id) const { auto it = cells.find(cell_id); if (it != cells.end()) return it->second.end(); return empty_cell.end(); } /* ------------------------------------------------------------------------ */ class cells_iterator : private std::iterator { public: explicit cells_iterator(typename std::map::const_iterator it) : it(it) {} cells_iterator & operator++() { this->it++; return *this; } cells_iterator operator++(int) { cells_iterator tmp(*this); operator++(); return tmp; }; bool operator==(const cells_iterator & rhs) const { return it == rhs.it; }; bool operator!=(const cells_iterator & rhs) const { return !operator==(rhs); }; CellID operator*() const { CellID cur_cell_id(this->it->first); return cur_cell_id; }; private: /// map iterator typename std::map::const_iterator it; }; public: template Cell & insert(const T & d, const vector_type & position) { auto && cell_id = getCellID(position); auto && it = cells.find(cell_id); if (it == cells.end()) { Cell cell(cell_id); auto & tmp = (cells[cell_id] = cell).add(d); for (UInt i = 0; i < dimension; ++i) { Real posl = center(i) + cell_id.getID(i) * spacing(i); Real posu = posl + spacing(i); if (posl < lower(i)) lower(i) = posl; if (posu > upper(i)) upper(i) = posu; } return tmp; } else { return it->second.add(d); } } /* ------------------------------------------------------------------------ */ inline decltype(auto) begin() const { auto begin = this->cells.begin(); return cells_iterator(begin); } inline decltype(auto) end() const { auto end = this->cells.end(); return cells_iterator(end); } template CellID getCellID(const vector_type & position) const { CellID cell_id(dimension); for (UInt i = 0; i < dimension; ++i) { cell_id.setID(i, getCellID(position(i), i)); } return cell_id; } void printself(std::ostream & stream, int indent = 0) const { std::string space; for (Int i = 0; i < indent; i++, space += AKANTU_INDENT) ; std::streamsize prec = stream.precision(); std::ios_base::fmtflags ff = stream.flags(); stream.setf(std::ios_base::showbase); stream.precision(5); stream << space << "SpatialGrid<" << debug::demangle(typeid(T).name()) << "> [" << std::endl; stream << space << " + dimension : " << this->dimension << std::endl; stream << space << " + lower bounds : {"; for (UInt i = 0; i < lower.size(); ++i) { if (i != 0) stream << ", "; stream << lower(i); }; stream << "}" << std::endl; stream << space << " + upper bounds : {"; for (UInt i = 0; i < upper.size(); ++i) { if (i != 0) stream << ", "; stream << upper(i); }; stream << "}" << std::endl; stream << space << " + spacing : {"; for (UInt i = 0; i < spacing.size(); ++i) { if (i != 0) stream << ", "; stream << spacing(i); }; stream << "}" << std::endl; stream << space << " + center : {"; for (UInt i = 0; i < center.size(); ++i) { if (i != 0) stream << ", "; stream << center(i); }; stream << "}" << std::endl; stream << space << " + nb_cells : " << this->cells.size() << "/"; Vector dist(this->dimension); dist = upper; dist -= lower; for (UInt i = 0; i < this->dimension; ++i) { dist(i) /= spacing(i); } UInt nb_cells = std::ceil(dist(0)); for (UInt i = 1; i < this->dimension; ++i) { nb_cells *= std::ceil(dist(i)); } stream << nb_cells << std::endl; stream << space << "]" << std::endl; stream.precision(prec); stream.flags(ff); } void saveAsMesh(Mesh & mesh) const; private: /* -------------------------------------------------------------------------- */ inline UInt getCellID(Real position, UInt direction) const { AKANTU_DEBUG_ASSERT(direction < center.size(), "The direction asked (" << direction << ") is out of range " << center.size()); Real dist_center = position - center(direction); Int id = std::floor(dist_center / spacing(direction)); // if(dist_center < 0) id--; return id; } friend class GridSynchronizer; public: AKANTU_GET_MACRO(LowerBounds, lower, const Vector &); AKANTU_GET_MACRO(UpperBounds, upper, const Vector &); AKANTU_GET_MACRO(Spacing, spacing, const Vector &); protected: UInt dimension; cells_container cells; Vector spacing; Vector center; Vector lower; Vector upper; Cell empty_cell; }; /// standard output stream operator template inline std::ostream & operator<<(std::ostream & stream, const SpatialGrid & _this) { _this.printself(stream); return stream; } } // namespace akantu #include "mesh.hh" namespace akantu { /* -------------------------------------------------------------------------- */ template void SpatialGrid::saveAsMesh(Mesh & mesh) const { auto & nodes = const_cast &>(mesh.getNodes()); - ElementType type; + ElementType type = _not_defined; switch (dimension) { case 1: type = _segment_2; break; case 2: type = _quadrangle_4; break; case 3: type = _hexahedron_8; break; + } mesh.addConnectivityType(type); auto & connectivity = const_cast &>(mesh.getConnectivity(type)); auto & uint_data = mesh.getDataPointer("tag_1", type); Vector pos(dimension); UInt global_id = 0; for (auto & cell_pair : cells) { UInt cur_node = nodes.size(); UInt cur_elem = connectivity.size(); const CellID & cell_id = cell_pair.first; for (UInt i = 0; i < dimension; ++i) pos(i) = center(i) + cell_id.getID(i) * spacing(i); nodes.push_back(pos); for (UInt i = 0; i < dimension; ++i) pos(i) += spacing(i); nodes.push_back(pos); connectivity.push_back(cur_node); switch (dimension) { case 1: connectivity(cur_elem, 1) = cur_node + 1; break; case 2: pos(0) -= spacing(0); nodes.push_back(pos); pos(0) += spacing(0); pos(1) -= spacing(1); nodes.push_back(pos); connectivity(cur_elem, 1) = cur_node + 3; connectivity(cur_elem, 2) = cur_node + 1; connectivity(cur_elem, 3) = cur_node + 2; break; case 3: pos(1) -= spacing(1); pos(2) -= spacing(2); nodes.push_back(pos); pos(1) += spacing(1); nodes.push_back(pos); pos(0) -= spacing(0); nodes.push_back(pos); pos(1) -= spacing(1); pos(2) += spacing(2); nodes.push_back(pos); pos(0) += spacing(0); nodes.push_back(pos); pos(0) -= spacing(0); pos(1) += spacing(1); nodes.push_back(pos); connectivity(cur_elem, 1) = cur_node + 2; connectivity(cur_elem, 2) = cur_node + 3; connectivity(cur_elem, 3) = cur_node + 4; connectivity(cur_elem, 4) = cur_node + 5; connectivity(cur_elem, 5) = cur_node + 6; connectivity(cur_elem, 6) = cur_node + 1; connectivity(cur_elem, 7) = cur_node + 7; break; } uint_data.push_back(global_id); ++global_id; } } } // namespace akantu #endif /* __AKANTU_AKA_GRID_DYNAMIC_HH__ */ diff --git a/src/fe_engine/element_class_structural.hh b/src/fe_engine/element_class_structural.hh index d719c6e11..2dc673bdf 100644 --- a/src/fe_engine/element_class_structural.hh +++ b/src/fe_engine/element_class_structural.hh @@ -1,254 +1,254 @@ /** * @file element_class_structural.hh * * @author Fabian Barras * @author Lucas Frerot * @author Nicolas Richart * @author Damien Spielmann * * @date creation: Thu Feb 21 2013 * @date last modification: Tue Feb 20 2018 * * @brief Specialization of the element classes for structural elements * * @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 "element_class.hh" /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_ELEMENT_CLASS_STRUCTURAL_HH__ #define __AKANTU_ELEMENT_CLASS_STRUCTURAL_HH__ namespace akantu { /// Macro to generate the InterpolationProperty structures for different /// interpolation types #define AKANTU_DEFINE_STRUCTURAL_INTERPOLATION_TYPE_PROPERTY( \ itp_type, itp_geom_type, ndof, nb_stress, nb_dnds_cols) \ template <> struct InterpolationProperty { \ static const InterpolationKind kind{_itk_structural}; \ static const UInt nb_nodes_per_element{ \ InterpolationProperty::nb_nodes_per_element}; \ static const InterpolationType itp_geometry_type{itp_geom_type}; \ static const UInt natural_space_dimension{ \ InterpolationProperty::natural_space_dimension}; \ static const UInt nb_degree_of_freedom{ndof}; \ static const UInt nb_stress_components{nb_stress}; \ static const UInt dnds_columns{nb_dnds_cols}; \ } /* -------------------------------------------------------------------------- */ template class InterpolationElement { public: using interpolation_property = InterpolationProperty; /// compute the shape values for a given set of points in natural coordinates static inline void computeShapes(const Matrix & natural_coord, const Matrix & real_coord, Tensor3 & N) { for (UInt i = 0; i < natural_coord.cols(); ++i) { Matrix n_t = N(i); computeShapes(natural_coord(i), real_coord, n_t); } } /// compute the shape values for a given point in natural coordinates static inline void computeShapes(const Vector & natural_coord, const Matrix & real_coord, Matrix & N); /// compute shape derivatives (input is dxds) for a set of points static inline void computeShapeDerivatives(const Tensor3 & Js, const Tensor3 & DNDSs, const Matrix & R, Tensor3 & Bs) { for (UInt i = 0; i < Js.size(2); ++i) { Matrix J = Js(i); Matrix DNDS = DNDSs(i); Matrix DNDX(DNDS.rows(), DNDS.cols()); auto inv_J = J.inverse(); DNDX.mul(inv_J, DNDS); Matrix B_R = Bs(i); Matrix B(B_R.rows(), B_R.cols()); arrangeInVoigt(DNDX, B); B_R.mul(B, R); } } /** * compute @f$ B_{ij} = \frac{\partial N_j}{\partial S_i} @f$ the variation of * shape functions along with variation of natural coordinates on a given set * of points in natural coordinates */ static inline void computeDNDS(const Matrix & natural_coord, const Matrix & real_coord, Tensor3 & dnds) { for (UInt i = 0; i < natural_coord.cols(); ++i) { Matrix dnds_t = dnds(i); computeDNDS(natural_coord(i), real_coord, dnds_t); } } /** * compute @f$ B_{ij} = \frac{\partial N_j}{\partial S_i} @f$ the variation of * shape functions along with * variation of natural coordinates on a given point in natural * coordinates */ static inline void computeDNDS(const Vector & natural_coord, const Matrix & real_coord, Matrix & dnds); /** * arrange B in Voigt notation from DNDS */ static inline void arrangeInVoigt(const Matrix & dnds, Matrix & B) { // Default implementation assumes dnds is already in Voigt notation B.deepCopy(dnds); } public: static inline constexpr auto getNbNodesPerInterpolationElement() { return interpolation_property::nb_nodes_per_element; } static inline constexpr auto getShapeSize() { return interpolation_property::nb_nodes_per_element * interpolation_property::nb_degree_of_freedom * interpolation_property::nb_degree_of_freedom; } static inline constexpr auto getShapeDerivativesSize() { return interpolation_property::nb_nodes_per_element * interpolation_property::nb_degree_of_freedom * interpolation_property::nb_stress_components; } static inline constexpr auto getNaturalSpaceDimension() { return interpolation_property::natural_space_dimension; } static inline constexpr auto getNbDegreeOfFreedom() { return interpolation_property::nb_degree_of_freedom; } static inline constexpr auto getNbStressComponents() { return interpolation_property::nb_stress_components; } }; /// Macro to generate the element class structures for different structural /// element types /* -------------------------------------------------------------------------- */ #define AKANTU_DEFINE_STRUCTURAL_ELEMENT_CLASS_PROPERTY( \ elem_type, geom_type, interp_type, parent_el_type, elem_kind, sp, \ gauss_int_type, min_int_order) \ template <> struct ElementClassProperty { \ static const GeometricalType geometrical_type{geom_type}; \ static const InterpolationType interpolation_type{interp_type}; \ static const ElementType parent_element_type{parent_el_type}; \ static const ElementKind element_kind{elem_kind}; \ static const UInt spatial_dimension{sp}; \ static const GaussIntegrationType gauss_integration_type{gauss_int_type}; \ static const UInt polynomial_degree{min_int_order}; \ } /* -------------------------------------------------------------------------- */ /* ElementClass for structural elements */ /* -------------------------------------------------------------------------- */ template class ElementClass : public GeometricalElement< ElementClassProperty::geometrical_type>, public InterpolationElement< ElementClassProperty::interpolation_type> { protected: using geometrical_element = GeometricalElement::geometrical_type>; using interpolation_element = InterpolationElement< ElementClassProperty::interpolation_type>; using parent_element = ElementClass::parent_element_type>; public: static inline void computeRotationMatrix(Matrix & /*R*/, const Matrix & /*X*/, const Vector & /*extra_normal*/) { AKANTU_TO_IMPLEMENT(); } /// compute jacobian (or integration variable change factor) for a given point static inline void computeJMat(const Vector & natural_coords, const Matrix & Xs, Matrix & J) { Matrix dnds(Xs.rows(), Xs.cols()); parent_element::computeDNDS(natural_coords, dnds); J.mul(dnds, Xs); } static inline void computeJMat(const Matrix & natural_coords, const Matrix & Xs, Tensor3 & Js) { for (UInt i = 0; i < natural_coords.cols(); ++i) { // because non-const l-value reference does not bind to r-value Matrix J = Js(i); computeJMat(Vector(natural_coords(i)), Xs, J); } } static inline void computeJacobian(const Matrix & natural_coords, const Matrix & node_coords, Vector & jacobians) { using itp = typename interpolation_element::interpolation_property; Tensor3 Js(itp::natural_space_dimension, itp::natural_space_dimension, natural_coords.cols()); computeJMat(natural_coords, node_coords, Js); for (UInt i = 0; i < natural_coords.cols(); ++i) { Matrix J = Js(i); jacobians(i) = J.det(); } } static inline void computeRotation(const Matrix & node_coords, Matrix & rotation); public: static AKANTU_GET_MACRO_NOT_CONST(Kind, _ek_structural, ElementKind); static AKANTU_GET_MACRO_NOT_CONST(P1ElementType, _not_defined, ElementType); static AKANTU_GET_MACRO_NOT_CONST(FacetType, _not_defined, ElementType); static constexpr auto getFacetType(__attribute__((unused)) UInt t = 0) { return _not_defined; } static constexpr AKANTU_GET_MACRO_NOT_CONST( SpatialDimension, ElementClassProperty::spatial_dimension, UInt); static constexpr auto getFacetTypes() { return ElementClass<_not_defined>::getFacetTypes(); } }; } // namespace akantu /* -------------------------------------------------------------------------- */ -#include "element_classes/element_class_hermite_inline_impl.cc" +#include "element_class_hermite_inline_impl.cc" /* keep order */ -#include "element_classes/element_class_bernoulli_beam_inline_impl.cc" -#include "element_classes/element_class_kirchhoff_shell_inline_impl.cc" +#include "element_class_bernoulli_beam_inline_impl.cc" +#include "element_class_kirchhoff_shell_inline_impl.cc" /* -------------------------------------------------------------------------- */ #endif /* __AKANTU_ELEMENT_CLASS_STRUCTURAL_HH__ */ diff --git a/src/mesh_utils/global_ids_updater.cc b/src/mesh_utils/global_ids_updater.cc index 0f906b957..8b64d3235 100644 --- a/src/mesh_utils/global_ids_updater.cc +++ b/src/mesh_utils/global_ids_updater.cc @@ -1,131 +1,136 @@ /** * @file global_ids_updater.cc * * @author Marco Vocialta * * @date creation: Fri Apr 13 2012 * @date last modification: Fri Dec 08 2017 * * @brief Functions of the GlobalIdsUpdater * * @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 . * */ /* -------------------------------------------------------------------------- */ #include "global_ids_updater.hh" #include "element_synchronizer.hh" #include "mesh_accessor.hh" #include "mesh_utils.hh" /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ namespace akantu { UInt GlobalIdsUpdater::updateGlobalIDs(UInt local_nb_new_nodes) { UInt total_nb_new_nodes = this->updateGlobalIDsLocally(local_nb_new_nodes); if (mesh.isDistributed()) { this->synchronizeGlobalIDs(); } return total_nb_new_nodes; } UInt GlobalIdsUpdater::updateGlobalIDsLocally(UInt local_nb_new_nodes) { const auto & comm = mesh.getCommunicator(); Int rank = comm.whoAmI(); Int nb_proc = comm.getNbProc(); if (nb_proc == 1) return local_nb_new_nodes; /// resize global ids array Array & nodes_global_ids = mesh.getGlobalNodesIds(); UInt old_nb_nodes = mesh.getNbNodes() - local_nb_new_nodes; nodes_global_ids.resize(mesh.getNbNodes(), -1); /// compute the number of global nodes based on the number of old nodes Matrix local_master_nodes(2, nb_proc, 0); for (UInt n = 0; n < old_nb_nodes; ++n) if (mesh.isLocalOrMasterNode(n)) ++local_master_nodes(0, rank); /// compute amount of local or master doubled nodes for (UInt n = old_nb_nodes; n < mesh.getNbNodes(); ++n) if (mesh.isLocalOrMasterNode(n)) ++local_master_nodes(1, rank); comm.allGather(local_master_nodes); local_master_nodes = local_master_nodes.transpose(); UInt old_global_nodes = std::accumulate(local_master_nodes(0).storage(), local_master_nodes(0).storage() + nb_proc, 0); /// update global number of nodes UInt total_nb_new_nodes = std::accumulate(local_master_nodes(1).storage(), local_master_nodes(1).storage() + nb_proc, 0); if (total_nb_new_nodes == 0) return 0; /// set global ids of local and master nodes UInt starting_index = std::accumulate(local_master_nodes(1).storage(), local_master_nodes(1).storage() + rank, old_global_nodes); for (UInt n = old_nb_nodes; n < mesh.getNbNodes(); ++n) { if (mesh.isLocalOrMasterNode(n)) { nodes_global_ids(n) = starting_index; ++starting_index; } else { nodes_global_ids(n) = -1; } } MeshAccessor mesh_accessor(mesh); mesh_accessor.setNbGlobalNodes(old_global_nodes + total_nb_new_nodes); return total_nb_new_nodes; } void GlobalIdsUpdater::synchronizeGlobalIDs() { this->reduce = true; this->synchronizer.slaveReductionOnce(*this, _gst_giu_global_conn); +#ifndef AKANTU_NDEBUG for (auto node : nodes_types) { auto node_type = mesh.getNodeType(node.first); if (node_type != _nt_pure_ghost) continue; auto n = 0u; + for (auto & pair : node.second) { if (std::get<1>(pair) == _nt_pure_ghost) ++n; } - if (n == node.second.size()) + + if (n == node.second.size()) { AKANTU_DEBUG_WARNING( "The node " << n << "is ghost on all the neighboring processors"); + } } +#endif this->reduce = false; this->synchronizer.synchronizeOnce(*this, _gst_giu_global_conn); } } // akantu diff --git a/src/mesh_utils/global_ids_updater_inline_impl.cc b/src/mesh_utils/global_ids_updater_inline_impl.cc index 390e9cfbb..14f427431 100644 --- a/src/mesh_utils/global_ids_updater_inline_impl.cc +++ b/src/mesh_utils/global_ids_updater_inline_impl.cc @@ -1,120 +1,128 @@ /** * @file global_ids_updater_inline_impl.cc * * @author Marco Vocialta * * @date creation: Fri Oct 02 2015 * @date last modification: Sun Aug 13 2017 * * @brief Implementation of the inline functions of GlobalIdsUpdater * * @section LICENSE * * Copyright (©) 2015-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 "communicator.hh" #include "global_ids_updater.hh" #include "mesh.hh" /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_GLOBAL_IDS_UPDATER_INLINE_IMPL_CC__ #define __AKANTU_GLOBAL_IDS_UPDATER_INLINE_IMPL_CC__ namespace akantu { /* -------------------------------------------------------------------------- */ inline UInt GlobalIdsUpdater::getNbData(const Array & elements, const SynchronizationTag & tag) const { UInt size = 0; - if (tag == _gst_giu_global_conn) + if (tag == _gst_giu_global_conn){ size += Mesh::getNbNodesPerElementList(elements) * - (sizeof(UInt) + sizeof(NodeType)) + - sizeof(int); - + sizeof(UInt); +#ifndef AKANTU_NDEBUG + size += sizeof(NodeType) * Mesh::getNbNodesPerElementList(elements) + sizeof(int); +#endif + } return size; } /* -------------------------------------------------------------------------- */ inline void GlobalIdsUpdater::packData(CommunicationBuffer & buffer, const Array & elements, const SynchronizationTag & tag) const { if (tag != _gst_giu_global_conn) return; auto & global_nodes_ids = mesh.getGlobalNodesIds(); +#ifndef AKANTU_NDEBUG buffer << int(mesh.getCommunicator().whoAmI()); - +#endif for (auto & element : elements) { /// get element connectivity const Vector current_conn = const_cast(mesh).getConnectivity(element); /// loop on all connectivity nodes for (auto node : current_conn) { UInt index = -1; if ((this->reduce and mesh.isLocalOrMasterNode(node)) or (not this->reduce and not mesh.isPureGhostNode(node))) { index = global_nodes_ids(node); } buffer << index; +#ifndef AKANTU_NDEBUG buffer << mesh.getNodeType(node); +#endif } } } /* -------------------------------------------------------------------------- */ inline void GlobalIdsUpdater::unpackData(CommunicationBuffer & buffer, const Array & elements, const SynchronizationTag & tag) { if (tag != _gst_giu_global_conn) return; auto & global_nodes_ids = mesh.getGlobalNodesIds(); +#ifndef AKANTU_NDEBUG int proc; buffer >> proc; +#endif for (auto & element : elements) { /// get element connectivity Vector current_conn = const_cast(mesh).getConnectivity(element); /// loop on all connectivity nodes for (auto node : current_conn) { UInt index; buffer >> index; +#ifndef AKANTU_NDEBUG NodeType node_type; buffer >> node_type; - if (reduce) nodes_types[node].push_back(std::make_pair(proc, node_type)); +#endif if (index == UInt(-1)) continue; if (mesh.isSlaveNode(node)) global_nodes_ids(node) = index; } } } } // akantu #endif /* __AKANTU_GLOBAL_IDS_UPDATER_INLINE_IMPL_CC__ */