diff --git a/src/geometry/aabb_primitives/aabb_primitive.hh b/src/geometry/aabb_primitives/aabb_primitive.hh index 43a51b85d..98f6d85fc 100644 --- a/src/geometry/aabb_primitives/aabb_primitive.hh +++ b/src/geometry/aabb_primitives/aabb_primitive.hh @@ -1,75 +1,77 @@ /** * @file aabb_primitive.hh * * @author Lucas Frérot * * @date creation: Fri Mar 13 2015 * @date last modification: Fri Mar 13 2015 * * @brief Macro classe (primitive) for AABB CGAL algos * * @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 . * */ /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_AABB_PRIMITIVE_HH__ #define __AKANTU_AABB_PRIMITIVE_HH__ #include "aka_common.hh" #include "triangle.hh" #include "tetrahedron.hh" __BEGIN_AKANTU__ typedef CGAL::Cartesian Kernel; #define AKANTU_AABB_CLASS(name) \ class name##_primitive { \ typedef std::list< name >::iterator Iterator; \ \ public: \ typedef UInt Id; \ typedef Kernel::Point_3 Point; \ typedef Kernel::name##_3 Datum; \ \ public: \ name##_primitive() : meshId(0), primitive() {} \ name##_primitive(Iterator it) : meshId(it->id()), primitive(*it) {} \ \ public: \ const Datum & datum() const { return primitive; } \ const Point & reference_point() const { return primitive.vertex(0); } \ const Id & id() const { return meshId; } \ \ protected: \ Id meshId; \ name primitive; \ \ } +// If the primitive is supported by CGAL::intersection() then the +// implementation process is really easy with this macro AKANTU_AABB_CLASS(Triangle); AKANTU_AABB_CLASS(Tetrahedron); #undef AKANTU_AABB_CLASS __END_AKANTU__ #endif // __AKANTU_AABB_PRIMITIVE_HH__ diff --git a/src/geometry/geom_helper_functions.hh b/src/geometry/geom_helper_functions.hh index 45bef4dcd..e8a520336 100644 --- a/src/geometry/geom_helper_functions.hh +++ b/src/geometry/geom_helper_functions.hh @@ -1,71 +1,121 @@ /** * @file geom_helper_functions.hh * * @author Lucas Frérot * * @date creation: Wed Mar 4 2015 * @date last modification: Thu Mar 5 2015 * * @brief Helper functions for the computational geometry algorithms * * @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 . * */ /* -------------------------------------------------------------------------- */ #ifndef _AKANTU_GEOM_HELPER_FUNCTIONS_HH__ #define _AKANTU_GEOM_HELPER_FUNCTIONS_HH__ #include "aka_common.hh" +#include "aka_math.hh" +#include "tree_type_helper.hh" #include __BEGIN_AKANTU__ typedef CGAL::Cartesian K; #define EPS 1e-10 inline bool comparePoints(const K::Point_3 & a, const K::Point_3 & b) { Math::setTolerance(EPS); return Math::are_float_equal(a.x(), b.x()) && Math::are_float_equal(a.y(), b.y()); } inline bool compareSegments(const K::Segment_3 & a, const K::Segment_3 & b) { - Math::setTolerance(EPS); return (comparePoints(a.source(), b.source()) && comparePoints(a.target(), b.target())) || (comparePoints(a.source(), b.target()) && comparePoints(a.target(), b.source())); } -struct CompareSegments { - bool operator()(const K::Segment_3 & a, const K::Segment_3 & b) { - return compareSegments(a, b); +inline bool compareSegmentPairs(const std::pair & a, const std::pair & b) { + return compareSegments(a.first, b.first); +} + +inline bool comparePairElement(const std::pair & a, const std::pair & b) { + return a.second < b.second; +} + +/* -------------------------------------------------------------------------- */ +/* Predicates */ +/* -------------------------------------------------------------------------- */ + +// Predicate used to eliminate faces of mesh not belonging to a specific element +template +class BelongsNotToElement { + +public: + BelongsNotToElement(UInt el): + el(el) + {} + + bool operator()(const typename TreeTypeHelper::primitive_type & primitive) { + return primitive.id() != el; } + +protected: + const UInt el; +}; + +// Predicate used to determine if point is on edge of faces of mesh +template +class HasOnEdge { + +public: + HasOnEdge(const K::Point_3 & point): + point(point) + {} + + bool operator()(const typename TreeTypeHelper::primitive_type & primitive) { + return primitive.has_on(point); + } + +protected: + const K::Point_3 & point; }; -struct CompareSegmentPairs { - bool operator()(const std::pair & a, const std::pair & b) { - return compareSegments(a.first, b.first); +class IsSameSegment { + +public: + IsSameSegment(const K::Segment_3 & segment): + segment(segment) + {} + + bool operator()(const std::pair & test_pair) { + return compareSegments(segment, test_pair.first); } +protected: + const K::Segment_3 segment; }; + __END_AKANTU__ #endif // _AKANTU_GEOM_HELPER_FUNCTIONS_HH__ diff --git a/src/geometry/mesh_geom_container.cc b/src/geometry/mesh_geom_container.cc index 9f3329653..84925c1b5 100644 --- a/src/geometry/mesh_geom_container.cc +++ b/src/geometry/mesh_geom_container.cc @@ -1,153 +1,153 @@ /** * @file mesh_geom_container.cc * * @author Lucas Frérot * * @date creation: Fri Feb 27 2015 * @date last modification: Fri Mar 6 2015 * * @brief Contains the CGAL representation of a mesh * * @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 "aka_common.hh" #include "mesh_geom_container.hh" #include "mesh_geom_factory.hh" #include "mesh.hh" #include /* -------------------------------------------------------------------------- */ #define MESH_GEOM_CASE(d, type) \ MeshGeomFactory * factory = new MeshGeomFactory(mesh); \ factory_map(factory, type); \ factory->constructData(); #define MESH_GEOM_CASE_2D(type) MESH_GEOM_CASE(2, type) #define MESH_GEOM_CASE_3D(type) MESH_GEOM_CASE(3, type) /* -------------------------------------------------------------------------- */ __BEGIN_AKANTU__ typedef CGAL::Cartesian K; MeshGeomContainer::MeshGeomContainer(const Mesh & mesh): MeshGeomAbstract(mesh), factory_map(), interface_mesh(mesh.getSpatialDimension(), "mesh_geom") { interface_mesh.addConnectivityType(_segment_2, _not_ghost); interface_mesh.addConnectivityType(_segment_2, _ghost); interface_mesh.registerData("associated_element").alloc(0, 1, _segment_2); interface_mesh.registerData("material").alloc(0, 1, _segment_2); } MeshGeomContainer::~MeshGeomContainer() {} void MeshGeomContainer::constructData() { AKANTU_DEBUG_IN(); const UInt spatial_dim = mesh.getSpatialDimension(); Mesh::type_iterator it = mesh.firstType(spatial_dim, _not_ghost); Mesh::type_iterator end = mesh.lastType(spatial_dim, _not_ghost); /// Loop over the element types of the mesh and construct the primitive trees for (; it != end ; ++it) { ElementType type = *it; // for AKANTU_BOOST_ELEMENT_SWITCH macro switch(spatial_dim) { case 1: AKANTU_DEBUG_WARNING("Geometry in 1D is undefined"); break; case 2: // Expand the list of elements when they are implemented AKANTU_BOOST_ELEMENT_SWITCH(MESH_GEOM_CASE_2D, (_triangle_3)); break; case 3: - //AKANTU_BOOST_ELEMENT_SWITCH(MESH_GEOM_CASE_3D, (_tetrahedron_4)); + AKANTU_BOOST_ELEMENT_SWITCH(MESH_GEOM_CASE_3D, (_tetrahedron_4)); break; } } AKANTU_DEBUG_OUT(); } UInt MeshGeomContainer::numberOfIntersectionsWithInterface(const K::Segment_3 & interface) const { AKANTU_DEBUG_IN(); UInt total = 0; GeomMap::type_iterator it = factory_map.firstType(); GeomMap::type_iterator end = factory_map.lastType(); for (; it != end ; ++it) { total += factory_map(*it)->numberOfIntersectionsWithInterface(interface); } AKANTU_DEBUG_OUT(); return total; } void MeshGeomContainer::meshOfLinearInterface(const Interface & interface, Mesh & interface_mesh) { AKANTU_DEBUG_IN(); GeomMap::type_iterator it = factory_map.firstType(); GeomMap::type_iterator end = factory_map.lastType(); for (; it != end ; ++it) { factory_map(*it)->meshOfLinearInterface(interface, interface_mesh); } AKANTU_DEBUG_OUT(); } Mesh & MeshGeomContainer::meshOfLinearInterfaces(const std::list & interfaces) { std::list::const_iterator interfaces_it = interfaces.begin(); std::list::const_iterator interfaces_end = interfaces.end(); for (; interfaces_it != interfaces_end ; ++interfaces_it) { meshOfLinearInterface(*interfaces_it, interface_mesh); } return interface_mesh; } const MeshGeomAbstract * MeshGeomContainer::getFactoryForElementType(ElementType el_type) const { return factory_map(el_type); } /* -------------------------------------------------------------------------- */ #undef MESH_GEOM_CASE_2D #undef MESH_GEOM_CASE_3D #undef MESH_GEOM_CASE __END_AKANTU__ diff --git a/src/geometry/mesh_geom_factory.hh b/src/geometry/mesh_geom_factory.hh index 2be6fa7c9..1acc72d9a 100644 --- a/src/geometry/mesh_geom_factory.hh +++ b/src/geometry/mesh_geom_factory.hh @@ -1,87 +1,90 @@ /** * @file mesh_geom_factory.hh * * @author Lucas Frérot * * @date creation: Thu Feb 26 2015 * @date last modification: Fri Mar 6 2015 * * @brief Class for constructing the CGAL primitives of a mesh * * @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 . * */ /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_MESH_GEOM_FACTORY_HH__ #define __AKANTU_MESH_GEOM_FACTORY_HH__ #include "aka_common.hh" #include "mesh.hh" #include "mesh_geom_abstract.hh" #include "tree_type_helper.hh" #include /* -------------------------------------------------------------------------- */ __BEGIN_AKANTU__ typedef CGAL::Cartesian K; template class MeshGeomFactory : public MeshGeomAbstract { public: /// Construct from mesh explicit MeshGeomFactory(const Mesh & mesh); /// Desctructor virtual ~MeshGeomFactory(); public: /// Construct AABB tree for fast intersection computing virtual void constructData(); /// Compute the number of intersections with primitive virtual UInt numberOfIntersectionsWithInterface(const K::Segment_3 & interface) const; /// Create a mesh of the intersection with a linear interface virtual void meshOfLinearInterface(const Interface & interface, Mesh & interface_mesh); /// Construct a primitive and add it to the list void addPrimitive(const Matrix & node_coordinates, UInt id); /// Construct segment list from intersections and remove duplicates - void constructSegments(const std::list::linear_intersection> & intersections, - std::list > & segments); + void constructSegments( + const std::list::linear_intersection> & intersections, + std::list > & segments, + const K::Segment_3 & interface); const typename TreeTypeHelper::tree & getTree() const { return *data_tree; } protected: typename TreeTypeHelper::tree * data_tree; typename TreeTypeHelper::container_type primitive_list; }; + __END_AKANTU__ #include "mesh_geom_factory_tmpl.hh" #endif // __AKANTU_MESH_GEOM_FACTORY_HH__ diff --git a/src/geometry/mesh_geom_factory_tmpl.hh b/src/geometry/mesh_geom_factory_tmpl.hh index 218ad5723..e8411b9e0 100644 --- a/src/geometry/mesh_geom_factory_tmpl.hh +++ b/src/geometry/mesh_geom_factory_tmpl.hh @@ -1,216 +1,265 @@ /** * @file mesh_geom_factory_tmpl.hh * * @author Lucas Frérot * * @date creation: Thu Feb 26 2015 * @date last modification: Fri Mar 6 2015 * * @brief Class for constructing the CGAL primitives of a mesh * * @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 . * */ /* -------------------------------------------------------------------------- */ #ifndef _AKANTU_MESH_GEOM_FACTORY_TMPL_HH__ #define _AKANTU_MESH_GEOM_FACTORY_TMPL_HH__ #include "mesh_geom_factory.hh" #include "aka_common.hh" #include "tree_type_helper.hh" #include "triangle.hh" #include "geom_helper_functions.hh" +#include + #include /* -------------------------------------------------------------------------- */ __BEGIN_AKANTU__ typedef CGAL::Cartesian K; -template -MeshGeomFactory::MeshGeomFactory(const Mesh & mesh) : +template +MeshGeomFactory::MeshGeomFactory(const Mesh & mesh) : MeshGeomAbstract(mesh), data_tree(NULL), primitive_list() {} -template -MeshGeomFactory::~MeshGeomFactory() { +template +MeshGeomFactory::~MeshGeomFactory() { delete data_tree; } -template -void MeshGeomFactory::constructData() { +template +void MeshGeomFactory::constructData() { AKANTU_DEBUG_IN(); const GhostType ghost_type = _not_ghost; primitive_list.clear(); UInt nb_nodes_per_element = mesh.getNbNodesPerElement(type); const Array & connectivity = mesh.getConnectivity(type, ghost_type); const Array & nodes = mesh.getNodes(); Array::const_vector_iterator begin = connectivity.begin(nb_nodes_per_element); Array::const_vector_iterator it = connectivity.begin(nb_nodes_per_element); Array::const_vector_iterator end = connectivity.end(nb_nodes_per_element); /// This loop builds the list of primitives for (; it != end ; ++it) { const Vector & el_connectivity = *it; - Matrix node_coordinates(d, nb_nodes_per_element); + Matrix node_coordinates(dim, nb_nodes_per_element); for (UInt i = 0 ; i < nb_nodes_per_element ; i++) - for (UInt j = 0 ; j < d ; j++) + for (UInt j = 0 ; j < dim ; j++) node_coordinates(j, i) = nodes(el_connectivity(i), j); addPrimitive(node_coordinates, it - begin); } delete data_tree; - data_tree = new typename TreeTypeHelper::tree(primitive_list.begin(), primitive_list.end()); + data_tree = new typename TreeTypeHelper::tree(primitive_list.begin(), primitive_list.end()); AKANTU_DEBUG_OUT(); } // 2D and _triangle_3 implementation template<> void MeshGeomFactory<2, _triangle_3>::addPrimitive(const Matrix & node_coordinates, UInt id) { TreeTypeHelper<2, _triangle_3>::point_type a(node_coordinates(0, 0), node_coordinates(1, 0), 0.); TreeTypeHelper<2, _triangle_3>::point_type b(node_coordinates(0, 1), node_coordinates(1, 1), 0.); TreeTypeHelper<2, _triangle_3>::point_type c(node_coordinates(0, 2), node_coordinates(1, 2), 0.); Triangle t(a, b, c); t.setId(id); primitive_list.push_back(t); } // 3D and _tetrahedron_4 implementation -/* template<> void MeshGeomFactory<3, _tetrahedron_4>::addPrimitive(const Matrix & node_coordinates, UInt id) { - TreeTypeHelper<3, _tetrahedron_4>::point_type a(node_coordinates(0, 0), - node_coordinates(1, 0), - node_coordinates(2, 0)); - TreeTypeHelper<3, _tetrahedron_4>::point_type b(node_coordinates(0, 1), - node_coordinates(1, 1), - node_coordinates(2, 1)); - TreeTypeHelper<3, _tetrahedron_4>::point_type c(node_coordinates(0, 2), - node_coordinates(1, 2), - node_coordinates(2, 2)); - TreeTypeHelper<3, _tetrahedron_4>::point_type d(node_coordinates(0, 3), - node_coordinates(1, 3), - node_coordinates(2, 3)); - - Tetrahedron t(a, b, c, d); - t.setId(id); - primitive_list.push_back(t); + TreeTypeHelper<3, _tetrahedron_4>::point_type + a(node_coordinates(0, 0), node_coordinates(1, 0), node_coordinates(2, 0)), + b(node_coordinates(0, 1), node_coordinates(1, 1), node_coordinates(2, 1)), + c(node_coordinates(0, 2), node_coordinates(1, 2), node_coordinates(2, 2)), + d(node_coordinates(0, 3), node_coordinates(1, 3), node_coordinates(2, 3)); + + Triangle + t1(a, b, c), + t2(b, c, d), + t3(c, d, a), + t4(d, a, b); + + t1.setId(id); + t2.setId(id); + t3.setId(id); + t4.setId(id); + + primitive_list.push_back(t1); + primitive_list.push_back(t2); + primitive_list.push_back(t3); + primitive_list.push_back(t4); } -*/ -template -UInt MeshGeomFactory::numberOfIntersectionsWithInterface(const K::Segment_3 & interface) const { +template +UInt MeshGeomFactory::numberOfIntersectionsWithInterface(const K::Segment_3 & interface) const { return data_tree->number_of_intersected_primitives(interface); } -template -void MeshGeomFactory::meshOfLinearInterface(const Interface & interface_pair, Mesh & interface_mesh) { +template +void MeshGeomFactory::meshOfLinearInterface(const Interface & interface_pair, Mesh & interface_mesh) { AKANTU_DEBUG_IN(); const K::Segment_3 & interface = interface_pair.first; UInt number_of_intersections = this->numberOfIntersectionsWithInterface(interface); if (!number_of_intersections) { return; } - std::list::linear_intersection> list_of_intersections; + std::list::linear_intersection> list_of_intersections; std::list< std::pair > list_of_segments; // Contains no duplicate elements /// Compute all the intersection pairs (segment + element id) and remove duplicates data_tree->all_intersections(interface, std::back_inserter(list_of_intersections)); - this->constructSegments(list_of_intersections, list_of_segments); + this->constructSegments(list_of_intersections, list_of_segments, interface); /// Arrays for storing nodes and connectivity Array & nodes = interface_mesh.getNodes(); Array & connectivity = interface_mesh.getConnectivity(_segment_2); /// Arrays for storing associated element id and type Array & associated_element = interface_mesh.getData("associated_element", _segment_2); Array & associated_material = interface_mesh.getData("material", _segment_2); std::list >::iterator it = list_of_segments.begin(); std::list >::iterator end = list_of_segments.end(); /// Loop over the intersections pairs (segment, id) for (; it != end ; ++it) { Vector segment_connectivity(2); segment_connectivity(0) = interface_mesh.getNbNodes(); segment_connectivity(1) = interface_mesh.getNbNodes() + 1; connectivity.push_back(segment_connectivity); /// Copy nodes - Vector source(d), target(d); - for (UInt j = 0 ; j < d ; j++) { + Vector source(dim), target(dim); + for (UInt j = 0 ; j < dim ; j++) { source(j) = it->first.source()[j]; target(j) = it->first.target()[j]; } nodes.push_back(source); nodes.push_back(target); /// Copy associated element info associated_element.push_back(Element(el_type, it->second)); associated_material.push_back(interface_pair.second); } AKANTU_DEBUG_OUT(); } -template -void MeshGeomFactory::constructSegments(const std::list< typename TreeTypeHelper::linear_intersection > & intersections, - std::list > & segments) { +template +void MeshGeomFactory::constructSegments( + const std::list< typename TreeTypeHelper::linear_intersection > & intersections, + std::list > & segments, + const K::Segment_3 & interface) +{ AKANTU_DEBUG_IN(); - typename std::list::linear_intersection>::const_iterator int_it = intersections.begin(), - int_end = intersections.end(); + typename std::list::linear_intersection>::const_iterator + int_it = intersections.begin(), + int_end = intersections.end(); for (; int_it != int_end ; ++int_it) { + UInt el = (*int_it)->second; + if (const K::Segment_3 * segment = boost::get(&((*int_it)->first))) { - segments.push_back(std::make_pair(*segment, (*int_it)->second)); + if (std::find_if(segments.begin(), segments.end(), IsSameSegment(*segment)) == segments.end()) + segments.push_back(std::make_pair(*segment, (*int_it)->second)); } - } - segments.unique(CompareSegmentPairs()); + else if (const K::Point_3 * point = boost::get(&((*int_it)->first))) { + // We only want to treat points differently if we're in 3D with Tetra4 elements + // This should be optimized by compilator + if (dim == 3 && el_type == _tetrahedron_4) { + UInt nb_facets = Mesh::getNbFacetsPerElement(el_type); + typename TreeTypeHelper::container_type facets(nb_facets); + + // TODO Use mesh facets instead of this + std::remove_copy_if(primitive_list.begin(), + primitive_list.end(), + facets.begin(), + BelongsNotToElement(el)); + + typename TreeTypeHelper::tree * local_tree = + new typename TreeTypeHelper::tree(facets.begin(), facets.end()); + + std::list< typename TreeTypeHelper::linear_intersection> + local_intersections; + + local_tree->all_intersections(interface, std::back_inserter(local_intersections)); + + typename std::list::linear_intersection>::const_iterator + local_it = local_intersections.begin(), + local_end = local_intersections.end(); + + for (; local_it != local_end ; ++local_it) { + if (const K::Point_3 * local_point = boost::get(&((*local_it)->first))) { + if (!comparePoints(*point, *local_point)) { + K::Segment_3 seg(*point, *local_point); + if (std::find_if(segments.begin(), segments.end(), IsSameSegment(seg)) == segments.end()) + segments.push_back(std::make_pair(seg, el)); + } + } + } + + delete local_tree; + } + } + } AKANTU_DEBUG_OUT(); } __END_AKANTU__ #endif // _AKANTU_MESH_GEOM_FACTORY_TMPL_HH__ diff --git a/src/geometry/tree_type_helper.hh b/src/geometry/tree_type_helper.hh index 2f0362155..e79f6f28f 100644 --- a/src/geometry/tree_type_helper.hh +++ b/src/geometry/tree_type_helper.hh @@ -1,90 +1,81 @@ /** * @file tree_type_helper.hh * * @author Lucas Frérot * * @date creation: Fri Feb 27 2015 * @date last modification: Thu Mar 5 2015 * * @brief Converts element types of a mesh to CGAL primitive types * * @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 . * */ /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_TREE_TYPE_HELPER_HH__ #define __AKANTU_TREE_TYPE_HELPER_HH__ #include "aka_common.hh" #include "triangle.hh" #include "tetrahedron.hh" #include "aabb_primitive.hh" #include #include #include __BEGIN_AKANTU__ typedef CGAL::Cartesian Kernel; /* -------------------------------------------------------------------------- */ -template +template struct TreeTypeHelper; -/// 2D and _triangle_3 implementation -template<> -struct TreeTypeHelper<2, _triangle_3> { - typedef Triangle primitive_type; - typedef Kernel::Point_3 point_type; - - typedef std::list container_type; - typedef container_type::iterator iterator; - typedef Triangle_primitive aabb_primitive_type; - typedef CGAL::AABB_traits aabb_traits_type; - typedef CGAL::AABB_tree tree; - - typedef boost::optional < tree::Intersection_and_primitive_id< CGAL::Segment_3 >::Type > linear_intersection; -}; - -/// 3D and _tetrahedron_4 implementation -/* -template<> -struct TreeTypeHelper<3, _tetrahedron_4> { - typedef Tetrahedron primitive_type; - typedef Kernel::Point_3 point_type; - - typedef std::list container_type; - typedef container_type::iterator iterator; - typedef Tetrahedron_primitive aabb_primitive_type; - typedef CGAL::AABB_traits aabb_traits_type; - typedef CGAL::AABB_tree tree; - - typedef boost::optional < tree::Intersection_and_primitive_id< CGAL::Segment_3 >::Type > linear_intersection; -}; -*/ +#define TREE_TYPE_HELPER_MACRO(dim, el_type, my_primitive) \ + template<> \ + struct TreeTypeHelper { \ + typedef my_primitive primitive_type; \ + typedef my_primitive##_primitive aabb_primitive_type; \ + typedef Kernel::Point_3 point_type; \ + \ + typedef std::list container_type; \ + typedef container_type::iterator iterator; \ + typedef CGAL::AABB_traits aabb_traits_type; \ + typedef CGAL::AABB_tree tree; \ + \ + typedef boost::optional< \ + tree::Intersection_and_primitive_id< \ + CGAL::Segment_3 \ + >::Type > linear_intersection; \ + } + +TREE_TYPE_HELPER_MACRO(2, _triangle_3, Triangle); +TREE_TYPE_HELPER_MACRO(3, _tetrahedron_4, Triangle); + +#undef TREE_TYPE_HELPER_MACRO __END_AKANTU__ #endif // __AKANTU_TREE_TYPE_HELPER_HH__ diff --git a/src/model/solid_mechanics/embedded_interface.cc b/src/model/solid_mechanics/embedded_interface.cc index 7065bedad..9a62662eb 100644 --- a/src/model/solid_mechanics/embedded_interface.cc +++ b/src/model/solid_mechanics/embedded_interface.cc @@ -1,69 +1,71 @@ /** * @file embedded_interface.cc * * @author Lucas Frerot * * @date creation: Mon Mar 23 2015 * @date last modification: Mon Mar 23 2015 * * @brief Class used to represent embedded interfaces * * @section LICENSE * * Copyright (©) 2010-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 "aka_common.hh" #include "embedded_interface.hh" #include __BEGIN_AKANTU__ typedef CGAL::Cartesian K; EmbeddedInterface::EmbeddedInterface(UInt dim, const ID & id) : Parsable(_st_embedded_interface, id), points(2, dim) { registerParam("name", name, _pat_parsable | _pat_readable, "ID"); registerParam("material", mat_name, _pat_parsable | _pat_readable, "Material Name"); registerParam< Matrix >("points", points, _pat_parsable | _pat_readable, "Interface Points"); } EmbeddedInterface::~EmbeddedInterface() {} K::Segment_3 EmbeddedInterface::getPrimitive() const { + AKANTU_DEBUG_INFO("Creating interface " << points); + if (points.cols() == 2) { K::Point_3 a(points(0, 0), points(0, 1), 0.); K::Point_3 b(points(1, 0), points(1, 1), 0.); return K::Segment_3(a, b); } else if (points.cols() == 3) { K::Point_3 a(points(0, 0), points(0, 1), points(0, 2)); K::Point_3 b(points(1, 0), points(1, 1), points(1, 2)); return K::Segment_3(a, b); } else { AKANTU_DEBUG_ERROR("Error in reading embedded interface : wrong number of components"); return K::Segment_3(CGAL::ORIGIN, CGAL::ORIGIN); } } __END_AKANTU__ diff --git a/src/model/solid_mechanics/embedded_interface_model.cc b/src/model/solid_mechanics/embedded_interface_model.cc index c967f2eed..9b4242161 100644 --- a/src/model/solid_mechanics/embedded_interface_model.cc +++ b/src/model/solid_mechanics/embedded_interface_model.cc @@ -1,139 +1,144 @@ /** * @file embedded_interface_model.cc * * @author Lucas Frérot * * @date creation: Mon Mar 9 2015 * @date last modification: Mon Mar 9 2015 * * @brief Model of Solid Mechanics with embedded interfaces * * @section LICENSE * * Copyright (©) 2010-2012, 2014 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_common.hh" #include "embedded_interface_model.hh" #include "material_reinforcement.hh" #include /* -------------------------------------------------------------------------- */ __BEGIN_AKANTU__ typedef CGAL::Cartesian K; EmbeddedInterfaceModel::EmbeddedInterfaceModel(Mesh & mesh, UInt spatial_dimension, const ID & id, const MemoryID & memory_id) : SolidMechanicsModel(mesh, spatial_dimension, id, memory_id), interface_mesh(NULL), interface_material_selector(NULL), interface_container(mesh) {} EmbeddedInterfaceModel::~EmbeddedInterfaceModel() { delete interface_material_selector; } void EmbeddedInterfaceModel::initModel() { interface_container.constructData(); instanciateInterfaces(); SolidMechanicsModel::initModel(); } void EmbeddedInterfaceModel::initMaterials() { delete interface_material_selector; interface_material_selector = new EmbeddedInterfaceMaterialSelector("material", *this); this->setMaterialSelector(*interface_material_selector); Element element; //Material ** mat_val = &(materials.at(0)); for (ghost_type_t::iterator gt = ghost_type_t::begin(); gt != ghost_type_t::end(); ++gt) { element.ghost_type = *gt; Mesh::type_iterator it = interface_mesh->firstType(1, *gt); Mesh::type_iterator end = interface_mesh->lastType(1, *gt); for (; it != end ; ++it) { UInt nb_element = interface_mesh->getNbElement(*it, *gt); element.type = *it; Array & mat_indexes = material_index.alloc(nb_element, 1, *it, *gt); for (UInt el = 0 ; el < nb_element ; el++) { element.element = el; UInt mat_index = (*interface_material_selector)(element); AKANTU_DEBUG_ASSERT(mat_index < materials.size(), "The material selector returned an index that does not exist"); mat_indexes(element.element) = mat_index; materials.at(mat_index)->addElement(*it, el, *gt); } } } SolidMechanicsModel::initMaterials(); } void EmbeddedInterfaceModel::instanciateInterfaces() { const std::pair & sub_sect = this->parser->getSubSections(_st_embedded_interface); Parser::const_section_iterator section_it = sub_sect.first; std::list > interface_list; for (; section_it != sub_sect.second ; ++section_it) { EmbeddedInterface interface(this->spatial_dimension); interface.parseSection(*section_it); interface_list.push_back(std::make_pair(interface.getPrimitive(), interface.getMaterialName())); } initInterface(interface_list); } void EmbeddedInterfaceModel::initInterface(const std::list > & interface_list) { AKANTU_DEBUG_IN(); interface_mesh = &(interface_container.meshOfLinearInterfaces(interface_list)); registerFEEngineObject("EmbeddedInterfaceFEEngine", *interface_mesh, 1); FEEngine & engine = getFEEngine("EmbeddedInterfaceFEEngine"); engine.initShapeFunctions(_not_ghost); engine.initShapeFunctions(_ghost); AKANTU_DEBUG_OUT(); } void EmbeddedInterfaceModel::assembleStiffnessMatrix() { SolidMechanicsModel::assembleStiffnessMatrix(); } +void EmbeddedInterfaceModel::dump() { + SolidMechanicsModel::dump(); + //interface_mesh->dump(); +} + __END_AKANTU__ diff --git a/src/model/solid_mechanics/embedded_interface_model.hh b/src/model/solid_mechanics/embedded_interface_model.hh index 25a81a98e..77961b3fc 100644 --- a/src/model/solid_mechanics/embedded_interface_model.hh +++ b/src/model/solid_mechanics/embedded_interface_model.hh @@ -1,155 +1,156 @@ /** * @file embedded_interface_model.hh * * @author Lucas Frérot * * @date creation: Mon Mar 9 2015 * @date last modification: Mon Mar 9 2015 * * @brief Model of Solid Mechanics with embedded interfaces * * @section LICENSE * * Copyright (©) 2010-2012, 2014 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 . * */ /* -------------------------------------------------------------------------- */ #ifndef __AKANTU_EMBEDDED_INTERFACE_MODEL_HH__ #define __AKANTU_EMBEDDED_INTERFACE_MODEL_HH__ #include "aka_common.hh" #include "solid_mechanics_model.hh" #include "mesh.hh" #include "parsable.hh" #include "embedded_interface.hh" #include "mesh_geom_container.hh" #include /* -------------------------------------------------------------------------- */ __BEGIN_AKANTU__ typedef CGAL::Cartesian K; class EmbeddedInterfaceModel : public SolidMechanicsModel { typedef FEEngineTemplate MyFEEngineType; /* ------------------------------------------------------------------------ */ /* Constructors/Destructors */ /* ------------------------------------------------------------------------ */ public: /// Constructor EmbeddedInterfaceModel(Mesh & mesh, UInt spatial_dimension = _all_dimensions, const ID & id = "embedded_interface_model", const MemoryID & memory_id = 0); /// Destructor virtual ~EmbeddedInterfaceModel(); /* ------------------------------------------------------------------------ */ /* Methods */ /* ------------------------------------------------------------------------ */ public: /// Initialise the model virtual void initModel(); /// Initialise the materials virtual void initMaterials(); /// Instanciate the interfaces virtual void instanciateInterfaces(); /// Initialise the interface mesh void initInterface(const std::list > & interface_list); -public: /// Assemble the stiffness matrix of the model virtual void assembleStiffnessMatrix(); + /// Dump + virtual void dump(); /* ------------------------------------------------------------------------ */ /* Accessors */ /* ------------------------------------------------------------------------ */ public: /// get interface mesh AKANTU_GET_MACRO(InterfaceMesh, *interface_mesh, Mesh &); /// get associated elements AKANTU_GET_MACRO_BY_ELEMENT_TYPE(InterfaceAssociatedElements, interface_mesh->getData("associated_element"), Element); /* ------------------------------------------------------------------------ */ /* Class Members */ /* ------------------------------------------------------------------------ */ protected: /// Interface mesh (weak reference) Mesh * interface_mesh; /// Material selector for interface MaterialSelector * interface_material_selector; /// Geom object to build the interface mesh MeshGeomContainer interface_container; }; /// Material selector for EmbeddedInterfaceModel template class EmbeddedInterfaceMaterialSelector : public DefaultMaterialSelector {}; /// Material selector based on associated material name template<> class EmbeddedInterfaceMaterialSelector : public DefaultMaterialSelector { public: EmbeddedInterfaceMaterialSelector(const std::string & name, const EmbeddedInterfaceModel & model): DefaultMaterialSelector(model.getMaterialByElement()), names(model.getInterfaceMesh().getData(name)), model(model) {} UInt operator() (const Element & element) { try { DebugLevel dbl = debug::getDebugLevel(); debug::setDebugLevel(dblError); std::string material_name = names(element.type, element.ghost_type)(element.element); debug::setDebugLevel(dbl); return model.getMaterialIndex(material_name); } catch (...) { /// TODO this returns the last material, which should be a bulk material. Try to find cleaner way return model.getNbMaterials() - 1; } } protected: const ElementTypeMapArray & names; const EmbeddedInterfaceModel & model; }; __END_AKANTU__ #endif // __AKANTU_EMBEDDED_INTERFACE_MODEL_HH__ diff --git a/test/test_geometry/CMakeLists.txt b/test/test_geometry/CMakeLists.txt index a9a0abf8e..307793fab 100644 --- a/test/test_geometry/CMakeLists.txt +++ b/test/test_geometry/CMakeLists.txt @@ -1,44 +1,57 @@ #=============================================================================== # @file CMakeLists.txt # # @author Lucas Frerot # # @date creation: Fri Feb 27 2015 # @date last modification: Fri Feb 27 2015 # # @brief configuration for solver tests # # @section LICENSE # # Copyright (©) 2014 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(test_geometry_intersection_mesh mesh.geo 2 1) +add_mesh(test_geometry_intersection_triangle_3_mesh triangle_3.geo 2 1) +add_mesh(test_geometry_intersection_tetrahedron_4_mesh tetrahedron_4.geo 3 1) + +register_test(test_geometry_predicates + SOURCES test_geometry_predicates.cc + PACKAGE CGAL + ) + register_test(test_geometry_intersection SOURCES test_geometry_intersection.cc - DEPENDENCIES test_geometry_intersection_mesh + DEPENDENCIES test_geometry_intersection_triangle_3_mesh + PACKAGE CGAL + ) + +register_test(test_geometry_intersection_triangle_3 + SOURCES test_geometry_intersection_triangle_3.cc + DEPENDENCIES test_geometry_intersection_triangle_3_mesh PACKAGE CGAL ) -register_test(test_geometry_mesh - SOURCES test_geometry_mesh.cc - DEPENDENCIES test_geometry_intersection_mesh +register_test(test_geometry_intersection_tetrahedron_4 + SOURCES test_geometry_intersection_tetrahedron_4.cc + DEPENDENCIES test_geometry_intersection_tetrahedron_4_mesh PACKAGE CGAL ) diff --git a/test/test_geometry/test_geometry_intersection.cc b/test/test_geometry/test_geometry_intersection.cc index ed939cb60..6112dbd68 100644 --- a/test/test_geometry/test_geometry_intersection.cc +++ b/test/test_geometry/test_geometry_intersection.cc @@ -1,102 +1,102 @@ /** * @file test_geometry_intersection.cc * * @author Lucas Frérot * * @date creation: Fri Feb 27 2015 * @date last modification: Thu Mar 5 2015 * * @brief Tests the intersection module * * @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 "aka_common.hh" #include "mesh_geom_container.hh" #include "mesh_geom_factory.hh" #include "tree_type_helper.hh" #include "geom_helper_functions.hh" #include #include /* -------------------------------------------------------------------------- */ using namespace akantu; typedef CGAL::Cartesian K; typedef TreeTypeHelper<2, _triangle_3>::linear_intersection Line_intersection; typedef TreeTypeHelper<2, _triangle_3>::tree::Primitive_id Primitive_id; /* -------------------------------------------------------------------------- */ int main (int argc, char * argv[]) { debug::setDebugLevel(dblWarning); initialize("", argc, argv); Math::setTolerance(1e-10); Mesh mesh(2); - mesh.read("mesh.msh"); + mesh.read("triangle_3.msh"); MeshGeomContainer container(mesh); container.constructData(); const MeshGeomFactory<2, _triangle_3> * factory = dynamic_cast *>(container.getFactoryForElementType(_triangle_3)); const TreeTypeHelper<2, _triangle_3>::tree & tree = factory->getTree(); K::Point_3 a(0., 0.25, 0.), b(1., 0.25, 0.); K::Segment_3 line(a, b); if (container.numberOfIntersectionsWithInterface(line) != 2) return EXIT_FAILURE; K::Point_3 begin(a), intermediate(0.25, 0.25, 0.), end(0.75, 0.25, 0.); K::Segment_3 result_0(begin, intermediate), result_1(intermediate, end); std::list list_of_intersections; tree.all_intersections(line, std::back_inserter(list_of_intersections)); const Line_intersection & intersection_0 = list_of_intersections.front(); const Line_intersection & intersection_1 = list_of_intersections.back(); if (!intersection_0 || !intersection_1) return EXIT_FAILURE; /// *-> first is the intersection ; *->second is the primitive id if (const K::Segment_3 * segment = boost::get(&(intersection_0->first))) { if (!compareSegments(*segment, result_0)) { return EXIT_FAILURE; } } if (const K::Segment_3 * segment = boost::get(&(intersection_1->first))) { if (!compareSegments(*segment, result_1)) { return EXIT_FAILURE; } } finalize(); return EXIT_SUCCESS; } diff --git a/test/test_geometry/test_geometry_mesh.cc b/test/test_geometry/test_geometry_intersection_tetrahedron_4.cc similarity index 65% copy from test/test_geometry/test_geometry_mesh.cc copy to test/test_geometry/test_geometry_intersection_tetrahedron_4.cc index 45d756ac8..5487390dd 100644 --- a/test/test_geometry/test_geometry_mesh.cc +++ b/test/test_geometry/test_geometry_intersection_tetrahedron_4.cc @@ -1,93 +1,75 @@ /** - * @file test_geometry_mesh.cc + * @file test_geometry_intersection_tetrahedron_4.cc * * @author Lucas Frérot * - * @date creation: Fri Mar 13 2015 - * @date last modification: Fri Mar 13 2015 + * @date creation: Thu Mar 26 2015 + * @date last modification: Thu Mar 26 2015 * - * @brief Tests the interface mesh generation + * @brief Tests the intersection module with _tetrahedron_4 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 "aka_common.hh" - #include "mesh_geom_container.hh" +#include "mesh_geom_factory.hh" +#include "tree_type_helper.hh" #include "geom_helper_functions.hh" #include #include /* -------------------------------------------------------------------------- */ using namespace akantu; typedef CGAL::Cartesian K; +typedef K::Point_3 Point; +typedef K::Segment_3 Segment; /* -------------------------------------------------------------------------- */ int main (int argc, char * argv[]) { debug::setDebugLevel(dblWarning); initialize("", argc, argv); - Math::setTolerance(1e-10); - - Mesh mesh(2); - mesh.read("mesh.msh"); + Mesh mesh(3); + mesh.read("tetrahedron_4.msh"); MeshGeomContainer container(mesh); container.constructData(); - K::Point_3 a(0, 0.25, 0), - b(1, 0.25, 0), - c(0.25, 0, 0), - d(0.25, 1, 0); + Point point(1., 1., 1.); + Segment segment(CGAL::ORIGIN, point); - K::Segment_3 h_interface(a, b), - v_interface(c, d); + std::list > interfaces; + interfaces.push_back(std::make_pair(segment, "mat")); - std::list > interface_list; - interface_list.push_back(std::make_pair(h_interface, "mat")); - interface_list.push_back(std::make_pair(v_interface, "mat")); - - Mesh & interface_mesh = container.meshOfLinearInterfaces(interface_list); + Mesh & interface_mesh = container.meshOfLinearInterfaces(interfaces); if (interface_mesh.getNbElement(_segment_2) != 4) return EXIT_FAILURE; - Vector bary(2); - Element test; - test.element = 1; - test.type = _segment_2; - - interface_mesh.getBarycenter(test, bary); - - if (!Math::are_float_equal(bary(0), 0.5) || !Math::are_float_equal(bary(1), 0.25)) - return EXIT_FAILURE; - - finalize(); return EXIT_SUCCESS; } - - diff --git a/test/test_geometry/test_geometry_mesh.cc b/test/test_geometry/test_geometry_intersection_triangle_3.cc similarity index 98% copy from test/test_geometry/test_geometry_mesh.cc copy to test/test_geometry/test_geometry_intersection_triangle_3.cc index 45d756ac8..8c437e2e8 100644 --- a/test/test_geometry/test_geometry_mesh.cc +++ b/test/test_geometry/test_geometry_intersection_triangle_3.cc @@ -1,93 +1,93 @@ /** * @file test_geometry_mesh.cc * * @author Lucas Frérot * * @date creation: Fri Mar 13 2015 * @date last modification: Fri Mar 13 2015 * * @brief Tests the interface mesh generation * * @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 "aka_common.hh" #include "mesh_geom_container.hh" #include "geom_helper_functions.hh" #include #include /* -------------------------------------------------------------------------- */ using namespace akantu; typedef CGAL::Cartesian K; /* -------------------------------------------------------------------------- */ int main (int argc, char * argv[]) { debug::setDebugLevel(dblWarning); initialize("", argc, argv); Math::setTolerance(1e-10); Mesh mesh(2); - mesh.read("mesh.msh"); + mesh.read("triangle_3.msh"); MeshGeomContainer container(mesh); container.constructData(); K::Point_3 a(0, 0.25, 0), b(1, 0.25, 0), c(0.25, 0, 0), d(0.25, 1, 0); K::Segment_3 h_interface(a, b), v_interface(c, d); std::list > interface_list; interface_list.push_back(std::make_pair(h_interface, "mat")); interface_list.push_back(std::make_pair(v_interface, "mat")); Mesh & interface_mesh = container.meshOfLinearInterfaces(interface_list); if (interface_mesh.getNbElement(_segment_2) != 4) return EXIT_FAILURE; Vector bary(2); Element test; test.element = 1; test.type = _segment_2; interface_mesh.getBarycenter(test, bary); if (!Math::are_float_equal(bary(0), 0.5) || !Math::are_float_equal(bary(1), 0.25)) return EXIT_FAILURE; finalize(); return EXIT_SUCCESS; } diff --git a/test/test_geometry/test_geometry_mesh.cc b/test/test_geometry/test_geometry_predicates.cc similarity index 53% rename from test/test_geometry/test_geometry_mesh.cc rename to test/test_geometry/test_geometry_predicates.cc index 45d756ac8..e9bdbd1cb 100644 --- a/test/test_geometry/test_geometry_mesh.cc +++ b/test/test_geometry/test_geometry_predicates.cc @@ -1,93 +1,63 @@ /** - * @file test_geometry_mesh.cc + * @file test_geometry_predicates.cc * * @author Lucas Frérot * - * @date creation: Fri Mar 13 2015 - * @date last modification: Fri Mar 13 2015 + * @date creation: Thu Mar 26 2015 + * @date last modification: Thu Mar 26 2015 * - * @brief Tests the interface mesh generation + * @brief Tests the geometry predicates * * @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 "aka_common.hh" - -#include "mesh_geom_container.hh" #include "geom_helper_functions.hh" #include #include /* -------------------------------------------------------------------------- */ using namespace akantu; typedef CGAL::Cartesian K; +typedef K::Point_3 Point; +typedef K::Segment_3 Segment; -/* -------------------------------------------------------------------------- */ - -int main (int argc, char * argv[]) { +int main(int argc, char * argv[]) { debug::setDebugLevel(dblWarning); initialize("", argc, argv); - Math::setTolerance(1e-10); - - Mesh mesh(2); - mesh.read("mesh.msh"); - - MeshGeomContainer container(mesh); - container.constructData(); - - K::Point_3 a(0, 0.25, 0), - b(1, 0.25, 0), - c(0.25, 0, 0), - d(0.25, 1, 0); - - K::Segment_3 h_interface(a, b), - v_interface(c, d); + Point a(0, 1, 0); + Point b(0, 1, 1); - std::list > interface_list; - interface_list.push_back(std::make_pair(h_interface, "mat")); - interface_list.push_back(std::make_pair(v_interface, "mat")); + Segment seg1(a, b); + Segment seg2(b, a); - Mesh & interface_mesh = container.meshOfLinearInterfaces(interface_list); - - if (interface_mesh.getNbElement(_segment_2) != 4) - return EXIT_FAILURE; - - Vector bary(2); - Element test; - test.element = 1; - test.type = _segment_2; - - interface_mesh.getBarycenter(test, bary); - - if (!Math::are_float_equal(bary(0), 0.5) || !Math::are_float_equal(bary(1), 0.25)) + if (!compareSegments(seg1, seg2)) return EXIT_FAILURE; finalize(); return EXIT_SUCCESS; } - - diff --git a/test/test_geometry/tetrahedron_4.geo b/test/test_geometry/tetrahedron_4.geo new file mode 100644 index 000000000..dce326158 --- /dev/null +++ b/test/test_geometry/tetrahedron_4.geo @@ -0,0 +1,14 @@ +f = 1.0; +Point(1) = {0, 0, 0, f}; +Point(2) = {1, 0, 0, f}; +Point(3) = {1, 1, 0, f}; +Point(4) = {0, 1, 0, f}; +Line(1) = {4, 3}; +Line(2) = {3, 2}; +Line(3) = {1, 2}; +Line(4) = {4, 1}; +Line Loop(5) = {4, 3, -2, -1}; +Plane Surface(6) = {5}; +Extrude {0, 0, 1} { + Surface{6}; +} diff --git a/test/test_geometry/mesh.geo b/test/test_geometry/triangle_3.geo similarity index 100% rename from test/test_geometry/mesh.geo rename to test/test_geometry/triangle_3.geo