DOLFIN
DOLFIN C++ interface
Mesh.h
1 // Copyright (C) 2006-2016 Anders Logg
2 //
3 // This file is part of DOLFIN.
4 //
5 // DOLFIN is free software: you can redistribute it and/or modify
6 // it under the terms of the GNU Lesser General Public License as published by
7 // the Free Software Foundation, either version 3 of the License, or
8 // (at your option) any later version.
9 //
10 // DOLFIN is distributed in the hope that it will be useful,
11 // but WITHOUT ANY WARRANTY; without even the implied warranty of
12 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 // GNU Lesser General Public License for more details.
14 //
15 // You should have received a copy of the GNU Lesser General Public License
16 // along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
17 //
18 // Modified by Johan Hoffman 2007
19 // Modified by Magnus Vikstrøm 2007
20 // Modified by Garth N. Wells 2007-2015
21 // Modified by Niclas Jansson 2008
22 // Modified by Kristoffer Selim 2008
23 // Modified by Andre Massing 2009-2010
24 // Modified by Marie E. Rognes 2012
25 // Modified by Mikael Mortensen 2012
26 // Modified by Jan Blechta 2013
27 // Modified by Cecile Daversin-Catty 2018
28 
29 #ifndef __MESH_H
30 #define __MESH_H
31 
32 #include <memory>
33 #include <string>
34 #include <utility>
35 
36 #include <dolfin/common/Hierarchical.h>
37 #include <dolfin/common/MPI.h>
38 #include <dolfin/common/Variable.h>
39 #include "MeshData.h"
40 #include "MeshDomains.h"
41 #include "MeshGeometry.h"
42 #include "MeshConnectivity.h"
43 #include "MeshTopology.h"
44 
45 namespace dolfin
46 {
47  class BoundaryMesh;
48  class CellType;
49  class Expression;
50  class GenericFunction;
51  class LocalMeshData;
52  class MeshEntity;
53  class Point;
54  class SubDomain;
55  class BoundingBoxTree;
56 
82 
83  class Mesh : public Variable, public Hierarchical<Mesh>
84  {
85  public:
86 
88  Mesh();
89 
91  Mesh(MPI_Comm comm);
92 
97  Mesh(const Mesh& mesh);
98 
103  explicit Mesh(std::string filename);
104 
111  Mesh(MPI_Comm comm, std::string filename);
112 
119  Mesh(MPI_Comm comm, LocalMeshData& local_mesh_data);
120 
122  ~Mesh();
123 
128  const Mesh& operator=(const Mesh& mesh);
129 
135  std::size_t num_vertices() const
136  { return _topology.size(0); }
137 
143  std::size_t num_edges() const
144  { return _topology.size(1); }
145 
151  std::size_t num_faces() const
152  { return _topology.size(2); }
153 
159  std::size_t num_facets() const
160  { return _topology.size(_topology.dim() - 1); }
161 
167  std::size_t num_cells() const
168  { return _topology.size(_topology.dim()); }
169 
178  std::size_t num_entities(std::size_t d) const
179  { return _topology.size(d); }
180 
186  std::vector<double>& coordinates()
187  { return _geometry.x(); }
188 
194  const std::vector<double>& coordinates() const
195  { return _geometry.x(); }
196 
202  const std::vector<unsigned int>& cells() const
203  { return _topology(_topology.dim(), 0)(); }
204 
213  std::size_t num_entities_global(std::size_t dim) const
214  { return _topology.size_global(dim); }
215 
221  { return _topology; }
222 
227  const MeshTopology& topology() const
228  { return _topology; }
229 
235  { return _geometry; }
236 
241  const MeshGeometry& geometry() const
242  { return _geometry; }
243 
249  { return _domains; }
250 
255  const MeshDomains& domains() const
256  { return _domains; }
257 
267  std::shared_ptr<BoundingBoxTree> bounding_box_tree() const;
268 
273  MeshData& data();
274 
278  const MeshData& data() const;
279 
285  { dolfin_assert(_cell_type); return *_cell_type; }
286 
288  const CellType& type() const
289  { dolfin_assert(_cell_type); return *_cell_type; }
290 
298  std::size_t init(std::size_t dim) const;
299 
307  void init(std::size_t d0, std::size_t d1) const;
308 
310  void init() const;
311 
313  void init_global(std::size_t dim) const;
314 
318  void clean();
319 
323  void order();
324 
329  bool ordered() const;
330 
338  Mesh renumber_by_color() const;
339 
345  void scale(double factor);
346 
351  void translate(const Point& point);
352 
360  void rotate(double angle, std::size_t axis=2);
361 
370  void rotate(double angle, std::size_t axis, const Point& point);
371 
377  void smooth(std::size_t num_iterations=1);
378 
388  void smooth_boundary(std::size_t num_iterations=1,
389  bool harmonic_smoothing=true);
390 
399  void snap_boundary(const SubDomain& sub_domain,
400  bool harmonic_smoothing=true);
401 
414  const std::vector<std::size_t>& color(std::string coloring_type) const;
415 
427  const std::vector<std::size_t>&
428  color(std::vector<std::size_t> coloring_type) const;
429 
437  double hmin() const;
438 
446  double hmax() const;
447 
453  double rmin() const;
454 
460  double rmax() const;
461 
468  std::size_t hash() const;
469 
478  std::string str(bool verbose) const;
479 
486  const std::vector<int>& cell_orientations() const;
487 
494  void init_cell_orientations(const Expression& global_normal);
495 
498  MPI_Comm mpi_comm() const
499  { return _mpi_comm.comm(); }
500 
506  std::string ghost_mode() const;
507 
508  // Friend in fem_utils.h
510 
513  void build_mapping(std::shared_ptr<const Mesh> other) const;
514 
515  private:
516 
517  // Friends
518  friend class MeshEditor;
519  friend class TopologyComputation;
520  friend class MeshPartitioning;
521 
522  // Mesh topology
523  mutable MeshTopology _topology;
524 
525  // Mesh geometry
526  MeshGeometry _geometry;
527 
528  // Mesh domains
529  MeshDomains _domains;
530 
531  // Auxiliary mesh data
532  MeshData _data;
533 
534  // Bounding box tree used to compute collisions between the mesh
535  // and other objects. The tree is initialized to a zero pointer
536  // and is allocated and built when bounding_box_tree() is called.
537  mutable std::shared_ptr<BoundingBoxTree> _tree;
538 
539  // Cell type
540  std::unique_ptr<CellType> _cell_type;
541 
542  // True if mesh has been ordered
543  mutable bool _ordered;
544 
545  // Orientation of cells relative to a global direction
546  std::vector<int> _cell_orientations;
547 
548  // MPI communicator
549  dolfin::MPI::Comm _mpi_comm;
550 
551  // Ghost mode used for partitioning
552  std::string _ghost_mode;
553 
554  };
555 }
556 
557 #endif
dolfin::Mesh::cell_orientations
const std::vector< int > & cell_orientations() const
Definition: Mesh.cpp:437
dolfin::Mesh::init_cell_orientations
void init_cell_orientations(const Expression &global_normal)
Definition: Mesh.cpp:442
dolfin::Mesh::hash
std::size_t hash() const
Definition: Mesh.cpp:395
dolfin::Mesh::domains
MeshDomains & domains()
Definition: Mesh.h:248
dolfin::Mesh::bounding_box_tree
std::shared_ptr< BoundingBoxTree > bounding_box_tree() const
Definition: Mesh.cpp:347
dolfin::MeshGeometry::x
double x(std::size_t n, std::size_t i) const
Return value of coordinate with local index n in direction i.
Definition: MeshGeometry.h:99
dolfin::TopologyComputation
Definition: TopologyComputation.h:37
dolfin::Mesh::num_entities
std::size_t num_entities(std::size_t d) const
Definition: Mesh.h:178
dolfin::Mesh::init_global
void init_global(std::size_t dim) const
Compute global indices for entity dimension dim.
Definition: Mesh.cpp:233
dolfin::Mesh::hmin
double hmin() const
Definition: Mesh.cpp:359
dolfin::Mesh::renumber_by_color
Mesh renumber_by_color() const
Definition: Mesh.cpp:274
dolfin::Mesh::type
const CellType & type() const
Get mesh cell type (const version).
Definition: Mesh.h:288
dolfin::MeshTopology::size
std::size_t size(std::size_t dim) const
Return number of entities for given dimension.
Definition: MeshTopology.cpp:79
dolfin::Mesh::num_facets
std::size_t num_facets() const
Definition: Mesh.h:159
dolfin::Expression
Definition: Expression.h:50
dolfin::Mesh::init
void init() const
Compute all entities and connectivity.
Definition: Mesh.cpp:221
dolfin::Point
Definition: Point.h:41
dolfin::Mesh::smooth
void smooth(std::size_t num_iterations=1)
Definition: Mesh.cpp:301
dolfin::MPI::Comm::comm
MPI_Comm comm() const
Return the underlying MPI_Comm object.
Definition: MPI.cpp:117
dolfin::Mesh::hmax
double hmax() const
Definition: Mesh.cpp:368
dolfin::Variable
Common base class for DOLFIN variables.
Definition: Variable.h:36
dolfin::Mesh::mpi_comm
MPI_Comm mpi_comm() const
Definition: Mesh.h:498
dolfin::MeshTopology::dim
std::size_t dim() const
Return topological dimension.
Definition: MeshTopology.cpp:74
dolfin::Hierarchical
Definition: Hierarchical.h:44
dolfin::Mesh::Mesh
Mesh()
Create empty mesh.
Definition: Mesh.cpp:59
dolfin::Mesh::num_vertices
std::size_t num_vertices() const
Definition: Mesh.h:135
dolfin::Mesh::~Mesh
~Mesh()
Destructor.
Definition: Mesh.cpp:100
dolfin::Mesh::num_edges
std::size_t num_edges() const
Definition: Mesh.h:143
dolfin::Mesh::operator=
const Mesh & operator=(const Mesh &mesh)
Definition: Mesh.cpp:105
dolfin::Mesh::num_faces
std::size_t num_faces() const
Definition: Mesh.h:151
dolfin::Mesh
Definition: Mesh.h:84
dolfin::MeshDomains
Definition: MeshDomains.h:42
dolfin::Mesh::rmin
double rmin() const
Definition: Mesh.cpp:377
dolfin::Mesh::build_mapping
void build_mapping(std::shared_ptr< const Mesh > other) const
Definition: Mesh.cpp:492
dolfin::Mesh::topology
MeshTopology & topology()
Definition: Mesh.h:220
dolfin::MeshTopology
Definition: MeshTopology.h:53
dolfin::Mesh::type
CellType & type()
Definition: Mesh.h:284
dolfin::Mesh::coordinates
const std::vector< double > & coordinates() const
Definition: Mesh.h:194
dolfin::Mesh::rotate
void rotate(double angle, std::size_t axis=2)
Definition: Mesh.cpp:291
dolfin::Mesh::color
const std::vector< std::size_t > & color(std::string coloring_type) const
Definition: Mesh.cpp:316
dolfin::Mesh::order
void order()
Definition: Mesh.cpp:252
dolfin::Mesh::geometry
const MeshGeometry & geometry() const
Definition: Mesh.h:241
dolfin::Mesh::rmax
double rmax() const
Definition: Mesh.cpp:386
dolfin::Mesh::ghost_mode
std::string ghost_mode() const
Definition: Mesh.cpp:483
dolfin::Function
Definition: Function.h:66
dolfin::Mesh::cells
const std::vector< unsigned int > & cells() const
Definition: Mesh.h:202
dolfin::Mesh::topology
const MeshTopology & topology() const
Definition: Mesh.h:227
dolfin::LocalMeshData
This class stores mesh data on a local processor corresponding to a portion of a (larger) global mesh...
Definition: LocalMeshData.h:59
dolfin::Mesh::snap_boundary
void snap_boundary(const SubDomain &sub_domain, bool harmonic_smoothing=true)
Definition: Mesh.cpp:311
dolfin::Mesh::create_mesh
friend Mesh create_mesh(Function &coordinates)
dolfin::CellType
Definition: CellType.h:47
dolfin::SubDomain
Definition: SubDomain.h:43
dolfin::MeshEditor
Definition: MeshEditor.h:37
dolfin::MeshData
Definition: MeshData.h:59
dolfin::Mesh::domains
const MeshDomains & domains() const
Definition: Mesh.h:255
dolfin::Mesh::smooth_boundary
void smooth_boundary(std::size_t num_iterations=1, bool harmonic_smoothing=true)
Definition: Mesh.cpp:306
dolfin::MeshGeometry
MeshGeometry stores the geometry imposed on a mesh.
Definition: MeshGeometry.h:40
dolfin::MeshPartitioning
Definition: MeshPartitioning.h:63
dolfin::Mesh::translate
void translate(const Point &point)
Definition: Mesh.cpp:286
dolfin::Mesh::str
std::string str(bool verbose) const
Definition: Mesh.cpp:409
dolfin::MeshTopology::size_global
std::size_t size_global(std::size_t dim) const
Return global number of entities for given dimension.
Definition: MeshTopology.cpp:88
dolfin::Mesh::coordinates
std::vector< double > & coordinates()
Definition: Mesh.h:186
dolfin::Mesh::data
MeshData & data()
Definition: Mesh.cpp:129
dolfin::MPI::Comm
Definition: MPI.h:77
dolfin::Mesh::num_cells
std::size_t num_cells() const
Definition: Mesh.h:167
dolfin::Mesh::ordered
bool ordered() const
Definition: Mesh.cpp:264
dolfin::Mesh::geometry
MeshGeometry & geometry()
Definition: Mesh.h:234
dolfin::Mesh::num_entities_global
std::size_t num_entities_global(std::size_t dim) const
Definition: Mesh.h:213
dolfin
Definition: adapt.h:30
dolfin::Mesh::scale
void scale(double factor)
Definition: Mesh.cpp:281
dolfin::Mesh::clean
void clean()
Definition: Mesh.cpp:239