34 class OPENMEEG_EXPORT
Mesh {
48 Mesh(
Geometry* geometry=
nullptr): geom(create_geometry(geometry)) { }
56 Mesh(
const unsigned nv,
const unsigned nt,
Geometry* geometry=
nullptr);
68 Mesh(
const std::string& filename,
const bool verbose,
Geometry* geometry=
nullptr):
Mesh(geometry) {
69 load(filename,verbose);
75 Mesh(
const std::string& filename,
Geometry* geometry=
nullptr):
Mesh(filename,false,geometry) { }
81 std::string&
name() {
return mesh_name; }
82 const std::string&
name()
const {
return mesh_name; }
104 const TriangleIndices t = { indmap.at(inds[0]), indmap.at(inds[1]), indmap.at(inds[2])};
105 return add_triangle(t);
108 void add(
const std::vector<TriangleIndices>& trgs) {
109 for (
const auto& triangle : trgs)
110 add_triangle(triangle);
113 void add(
const std::vector<TriangleIndices>& trgs,
const IndexMap& indmap) {
114 for (
const auto& triangle : trgs)
115 add_triangle(triangle,indmap);
125 void info(
const bool verbose=
false)
const;
130 void update(
const bool topology_changed);
137 std::vector<size_t> indices;
138 for (
const auto& vertex : vertices())
139 indices.push_back(vertex->index());
140 std::sort(indices.begin(),indices.end());
142 for (
auto it=indices.begin(); it!=indices.end();) {
144 for (
auto it2=it1+1; it2!=indices.end() && *it2==*it1+1; it1=it2++);
145 result.push_back(
Range(*it,*it1));
162 std::map<Triangle*,unsigned> mapt;
164 for (
auto& vertex : triangle)
165 for (
const auto& t2 : triangles(*vertex))
167 result.push_back(t2);
174 for (
auto& triangle : triangles())
175 triangle.change_orientation();
192 void smooth(
const double& smoothing_intensity,
const unsigned& niter);
205 void load(
const std::string& filename,
const bool verbose=
true);
210 void save(
const std::string& filename)
const ;
212 #if !defined(SWIGPYTHON) && !defined(_MSC_VER)
218 void reference_vertices(
const IndexMap& indmap);
224 typedef std::map<std::pair<const Vertex*,const Vertex*>,
int> EdgeMap;
231 void add_mesh(
const Mesh& m);
235 const EdgeMap compute_edge_map()
const;
250 void make_adjacencies() {
251 vertex_triangles.clear();
252 for (
auto& triangle : triangles())
253 for (
const auto& vertex : triangle)
254 vertex_triangles[vertex].push_back(&triangle);
257 typedef std::shared_ptr<Geometry> Geom;
259 std::string mesh_name =
"";
260 VertexTriangles vertex_triangles;
264 bool outermost_ =
false;
268 bool current_barrier_ =
false;
269 bool isolated_ =
false;
Geometry contains the electrophysiological model Vertices, meshes and domains are stored in this geom...
Mesh & operator=(const Mesh &)=default
Triangle & add_triangle(const TriangleIndices inds, const IndexMap &indmap)
void laplacian(SymMatrix &A) const
Compute mesh laplacian.
bool operator!=(const Mesh &m) const
TriangleIndices triangle(const Triangle &t) const
bool outermost() const
Returns True if it is an outermost mesh.
bool has_correct_orientation() const
Check local orientation of mesh triangles.
void correct_local_orientation()
Correct the local orientation of the mesh triangles.
Mesh(const unsigned nv, const unsigned nt, Geometry *geometry=nullptr)
Constructor from scratch (vertices/triangles t be added)
const Triangles & triangles() const
bool intersection(const Mesh &) const
Check whether the mesh intersects another mesh.
void smooth(const double &smoothing_intensity, const unsigned &niter)
Smooth Mesh.
Normal normal(const Vertex &v) const
Get normal at vertex.`.
Geometry & geometry() const
void correct_global_orientation()
Correct the global orientation (if there is one).
void merge(const Mesh &, const Mesh &)
Merge two meshes.
void add(const std::vector< TriangleIndices > &trgs, const IndexMap &indmap)
Mesh(const std::string &filename, Geometry *geometry=nullptr)
const VerticesRefs & vertices() const
void add(const std::vector< TriangleIndices > &trgs)
const std::string & name() const
bool has_self_intersection() const
Check whether the mesh self-intersects.
void generate_indices()
Generate indices (if allocate).
TrianglesRefs triangles(const Vertex &V) const
Get the triangles adjacent to vertex.
std::map< const Vertex *, TrianglesRefs > VertexTriangles
void info(const bool verbose=false) const
Print info Print to std::cout some info about the mesh.
bool operator==(const Mesh &m) const
void save(const std::string &filename) const
Save mesh to file.
Mesh(const std::string &filename, const bool verbose, Geometry *geometry=nullptr)
Constructors.
void change_orientation()
Change mesh orientation.
void update(const bool topology_changed)
Recompute triangles normals, area, and vertex triangles.
Mesh(Geometry *geometry=nullptr)
Default constructor or constructor using a provided geometry.
bool current_barrier() const
Triangle & add_triangle(const TriangleIndices inds)
Add a triangle specified by its indices in the geometry.
TrianglesRefs adjacent_triangles(const Triangle &triangle) const
Get the triangles adjacent to.
VerticesRefs & vertices()
void gradient_norm2(SymMatrix &A) const
Compute the square norm of the surfacic gradient.
Ranges vertices_ranges() const
Get the ranges of the specific mesh in the global matrix.
Mesh(const Mesh &)=default
double solid_angle(const Vect3 &p) const
Given a point p, computes the solid angle of the mesh seen from.
void load(const std::string &filename, const bool verbose=true)
Read mesh from file.
Range triangles_range() const
std::vector< Triangle > Triangles
Vect3 crossprod(const Vect3 &V1, const Vect3 &V2)
std::vector< Mesh > Meshes
std::vector< Triangle * > TrianglesRefs
double det(const Vect3 &V1, const Vect3 &V2, const Vect3 &V3)
std::vector< Vertex * > VerticesRefs
std::map< unsigned, unsigned > IndexMap
double sqr(const double x)
double dotprod(const Vect3 &V1, const Vect3 &V2)