DOLFIN-X
DOLFIN-X C++ interface
Mesh.h
1 // Copyright (C) 2006-2020 Anders Logg, Chris Richardson and Garth N. Wells
2 //
3 // This file is part of DOLFINX (https://www.fenicsproject.org)
4 //
5 // SPDX-License-Identifier: LGPL-3.0-or-later
6 
7 #pragma once
8 
9 #include "Geometry.h"
10 #include "Topology.h"
11 #include "cell_types.h"
12 #include <Eigen/Dense>
13 #include <dolfinx/common/MPI.h>
14 #include <dolfinx/common/UniqueIdGenerator.h>
15 #include <string>
16 #include <utility>
17 
18 namespace dolfinx
19 {
20 
21 namespace fem
22 {
23 class CoordinateElement;
24 }
25 
26 namespace graph
27 {
28 template <typename T>
29 class AdjacencyList;
30 }
31 
32 namespace mesh
33 {
34 
36 enum class GhostMode : int
37 {
38  none,
39  shared_facet,
40  shared_vertex
41 };
42 
45 
46 class Mesh
47 {
48 public:
53  template <typename Topology, typename Geometry>
54  Mesh(MPI_Comm comm, Topology&& topology, Geometry&& geometry)
55  : _topology(std::forward<Topology>(topology)),
56  _geometry(std::forward<Geometry>(geometry)), _mpi_comm(comm)
57  {
58  // Do nothing
59  }
60 
80  [[deprecated]] Mesh(
81  MPI_Comm comm, mesh::CellType type,
82  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
83  Eigen::Dynamic, Eigen::RowMajor>>& x,
84  const Eigen::Ref<const Eigen::Array<std::int64_t, Eigen::Dynamic,
85  Eigen::Dynamic, Eigen::RowMajor>>&
86  cells,
87  const fem::CoordinateElement& element,
88  const std::vector<std::int64_t>& global_cell_indices,
89  const GhostMode ghost_mode, std::int32_t num_ghost_cells = 0);
90 
93  Mesh(const Mesh& mesh) = default;
94 
97  Mesh(Mesh&& mesh) = default;
98 
100  ~Mesh() = default;
101 
102  // Assignment operator
103  Mesh& operator=(const Mesh& mesh) = delete;
104 
107  Mesh& operator=(Mesh&& mesh) = default;
108 
109  // TODO: Is there any use for this? In many situations one has to get the
110  // topology of a const Mesh, which is done by Mesh::topology_mutable. Note
111  // that the python interface (calls Mesh::topology()) may still rely on it.
114  Topology& topology();
115 
118  const Topology& topology() const;
119 
122  Topology& topology_mutable() const;
123 
126  Geometry& geometry();
127 
130  const Geometry& geometry() const;
131 
136  double hmin() const;
137 
142  double hmax() const;
143 
146  double rmin() const;
147 
150  double rmax() const;
151 
156  std::size_t hash() const;
157 
160  std::size_t id() const { return _unique_id; }
161 
164  MPI_Comm mpi_comm() const;
165 
167  std::string name = "mesh";
168 
169 private:
170  // Mesh topology:
171  // TODO: This is mutable because of the current memory management within
172  // mesh::Topology. It allows to obtain a non-const Topology from a
173  // const mesh (via Mesh::topology_mutable()).
174  //
175  mutable Topology _topology;
176 
177  // Mesh geometry
178  Geometry _geometry;
179 
180  // MPI communicator
181  dolfinx::MPI::Comm _mpi_comm;
182 
183  // Unique identifier
184  std::size_t _unique_id = common::UniqueIdGenerator::id();
185 };
186 
188 Mesh create(MPI_Comm comm, const graph::AdjacencyList<std::int64_t>& cells,
189  const fem::CoordinateElement& element,
190  const Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic,
191  Eigen::RowMajor>& x,
192  GhostMode ghost_mode);
193 
194 } // namespace mesh
195 } // namespace dolfinx
dolfinx::mesh::Mesh::topology_mutable
Topology & topology_mutable() const
Get mesh topology if one really needs the mutable version.
Definition: Mesh.cpp:141
dolfinx::mesh::CellType
CellType
Cell type identifier.
Definition: cell_types.h:22
dolfinx::mesh::Mesh::id
std::size_t id() const
Get unique identifier for the mesh.
Definition: Mesh.h:160
dolfinx::graph::AdjacencyList
This class provides a static adjacency list data structure. It is commonly used to store directed gra...
Definition: assemble_matrix_impl.h:26
dolfinx::mesh::GhostMode
GhostMode
Enum for different partitioning ghost modes.
Definition: Mesh.h:36
dolfinx::mesh::create
Mesh create(MPI_Comm comm, const graph::AdjacencyList< std::int64_t > &cells, const fem::CoordinateElement &element, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &x, GhostMode ghost_mode)
Create a mesh.
Definition: Mesh.cpp:64
dolfinx::mesh::Mesh::topology
Topology & topology()
Get mesh topology.
Definition: Mesh.cpp:137
dolfinx::mesh::Mesh::Mesh
Mesh(MPI_Comm comm, Topology &&topology, Geometry &&geometry)
Create a mesh.
Definition: Mesh.h:54
dolfinx::mesh::Mesh::rmax
double rmax() const
Compute maximum cell inradius.
Definition: Mesh.cpp:153
dolfinx::mesh::Mesh
A Mesh consists of a set of connected and numbered mesh topological entities, and geometry data.
Definition: Mesh.h:46
dolfinx::mesh::Mesh::~Mesh
~Mesh()=default
Destructor.
dolfinx::MPI::Comm
A duplicate MPI communicator and manage lifetime of the communicator.
Definition: MPI.h:34
dolfinx::mesh::Mesh::name
std::string name
Name.
Definition: Mesh.h:167
dolfinx::mesh::Mesh::rmin
double rmin() const
Compute minimum cell inradius.
Definition: Mesh.cpp:151
dolfinx::mesh::Topology
Topology stores the topology of a mesh, consisting of mesh entities and connectivity (incidence relat...
Definition: Topology.h:58
dolfinx::mesh::Mesh::geometry
Geometry & geometry()
Get mesh geometry.
Definition: Mesh.cpp:143
dolfinx::mesh::Mesh::hash
std::size_t hash() const
Compute hash of mesh, currently based on the has of the mesh geometry and mesh topology.
Definition: Mesh.cpp:155
dolfinx::mesh::Mesh::mpi_comm
MPI_Comm mpi_comm() const
Mesh MPI communicator.
Definition: Mesh.cpp:169
dolfinx::mesh::Mesh::hmin
double hmin() const
Compute minimum cell size in mesh, measured greatest distance between any two vertices of a cell.
Definition: Mesh.cpp:147
dolfinx::mesh::Geometry
Geometry stores the geometry imposed on a mesh.
Definition: Geometry.h:38
dolfinx::common::UniqueIdGenerator::id
static std::size_t id()
Generate a unique ID.
Definition: UniqueIdGenerator.cpp:22
dolfinx::fem::CoordinateElement
This class manages coordinate mappings for isoparametric cells.
Definition: CoordinateElement.h:23
dolfinx::mesh::Mesh::hmax
double hmax() const
Compute maximum cell size in mesh, measured greatest distance between any two vertices of a cell.
Definition: Mesh.cpp:149