DOLFIN-X
DOLFIN-X C++ interface
utils.h
1 // Copyright (C) 2019-2020 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 <Eigen/Dense>
10 #include <dolfinx/graph/AdjacencyList.h>
11 
12 namespace dolfinx
13 {
14 namespace fem
15 {
16 class ElementDofLayout;
17 }
18 
19 namespace mesh
20 {
21 enum class CellType;
22 class Mesh;
23 
33 graph::AdjacencyList<std::int64_t>
34 extract_topology(const CellType& cell_type, const fem::ElementDofLayout& layout,
35  const graph::AdjacencyList<std::int64_t>& cells);
36 
38 Eigen::ArrayXd volume_entities(const Mesh& mesh,
39  const Eigen::Ref<const Eigen::ArrayXi>& entities,
40  int dim);
41 
43 Eigen::ArrayXd circumradius(const Mesh& mesh,
44  const Eigen::Ref<const Eigen::ArrayXi>& entities,
45  int dim);
46 
48 Eigen::ArrayXd h(const Mesh& mesh,
49  const Eigen::Ref<const Eigen::ArrayXi>& entities, int dim);
50 
52 Eigen::ArrayXd inradius(const Mesh& mesh,
53  const Eigen::Ref<const Eigen::ArrayXi>& entities);
54 
56 Eigen::ArrayXd radius_ratio(const Mesh& mesh,
57  const Eigen::Ref<const Eigen::ArrayXi>& entities);
58 
60 Eigen::Array<double, Eigen::Dynamic, 3, Eigen::RowMajor>
61 cell_normals(const Mesh& mesh, int dim);
62 
64 Eigen::Array<double, Eigen::Dynamic, 3, Eigen::RowMajor> midpoints(
65  const mesh::Mesh& mesh, int dim,
66  const Eigen::Ref<const Eigen::Array<int, Eigen::Dynamic, 1>>& entities);
67 
78 Eigen::Array<std::int32_t, Eigen::Dynamic, 1> locate_entities(
79  const mesh::Mesh& mesh, const int dim,
80  const std::function<Eigen::Array<bool, Eigen::Dynamic, 1>(
81  const Eigen::Ref<const Eigen::Array<double, 3, Eigen::Dynamic,
82  Eigen::RowMajor>>&)>& marker);
83 
105 Eigen::Array<std::int32_t, Eigen::Dynamic, 1> locate_entities_boundary(
106  const mesh::Mesh& mesh, const int dim,
107  const std::function<Eigen::Array<bool, Eigen::Dynamic, 1>(
108  const Eigen::Ref<const Eigen::Array<double, 3, Eigen::Dynamic,
109  Eigen::RowMajor>>&)>& marker);
110 
111 } // namespace mesh
112 } // namespace dolfinx
dolfinx::mesh::extract_topology
graph::AdjacencyList< std::int64_t > extract_topology(const CellType &cell_type, const fem::ElementDofLayout &layout, const graph::AdjacencyList< std::int64_t > &cells)
Extract topology from cell data, i.e. extract cell vertices.
Definition: utils.cpp:321
dolfinx::mesh::CellType
CellType
Cell type identifier.
Definition: cell_types.h:22
dolfinx::mesh::circumradius
Eigen::ArrayXd circumradius(const Mesh &mesh, const Eigen::Ref< const Eigen::ArrayXi > &entities, int dim)
Compute circumradius of mesh entities.
Definition: utils.cpp:393
dolfinx::mesh::cell_normals
Eigen::Array< double, Eigen::Dynamic, 3, Eigen::RowMajor > cell_normals(const Mesh &mesh, int dim)
Compute normal to given cell (viewed as embedded in 3D)
Definition: utils.cpp:462
dolfinx::mesh::volume_entities
Eigen::ArrayXd volume_entities(const Mesh &mesh, const Eigen::Ref< const Eigen::ArrayXi > &entities, int dim)
Compute (generalized) volume of mesh entities of given dimension.
Definition: utils.cpp:350
dolfinx::mesh::h
Eigen::ArrayXd h(const Mesh &mesh, const Eigen::Ref< const Eigen::ArrayXi > &entities, int dim)
Compute greatest distance between any two vertices.
Definition: utils.cpp:356
dolfinx::mesh::midpoints
Eigen::Array< double, Eigen::Dynamic, 3, Eigen::RowMajor > midpoints(const mesh::Mesh &mesh, int dim, const Eigen::Ref< const Eigen::Array< int, Eigen::Dynamic, 1 >> &entities)
Compute midpoints or mesh entities of a given dimension.
Definition: utils.cpp:538
dolfinx::mesh::locate_entities_boundary
Eigen::Array< std::int32_t, Eigen::Dynamic, 1 > locate_entities_boundary(const mesh::Mesh &mesh, const int dim, const std::function< Eigen::Array< bool, Eigen::Dynamic, 1 >(const Eigen::Ref< const Eigen::Array< double, 3, Eigen::Dynamic, Eigen::RowMajor >> &)> &marker)
Compute indicies of all mesh entities that are attached to an owned boundary facet and evaluate to tr...
Definition: utils.cpp:670
dolfinx::mesh::radius_ratio
Eigen::ArrayXd radius_ratio(const Mesh &mesh, const Eigen::Ref< const Eigen::ArrayXi > &entities)
Compute dim*inradius/circumradius for given cells.
Definition: utils.cpp:451
dolfinx::mesh::inradius
Eigen::ArrayXd inradius(const Mesh &mesh, const Eigen::Ref< const Eigen::ArrayXi > &entities)
Compute inradius of cells.
Definition: utils.cpp:399
dolfinx::mesh::locate_entities
Eigen::Array< std::int32_t, Eigen::Dynamic, 1 > locate_entities(const mesh::Mesh &mesh, const int dim, const std::function< Eigen::Array< bool, Eigen::Dynamic, 1 >(const Eigen::Ref< const Eigen::Array< double, 3, Eigen::Dynamic, Eigen::RowMajor >> &)> &marker)
Compute indicies of all mesh entities that evaluate to true for the provided geometric marking functi...
Definition: utils.cpp:601