DOLFIN-X
DOLFIN-X C++ interface
CoordinateElement.h
1 // Copyright (C) 2018-2020 Garth N. Wells and Chris N. Richardson
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 "ElementDofLayout.h"
10 #include <Eigen/Core>
11 #include <cstdint>
12 #include <dolfinx/common/span.hpp>
13 #include <dolfinx/mesh/cell_types.h>
14 #include <functional>
15 #include <memory>
16 #include <string>
17 #include <unsupported/Eigen/CXX11/Tensor>
18 
19 namespace dolfinx::fem
20 {
21 
22 // FIXME: A dof layout on a reference cell needs to be defined.
24 
26 {
27 public:
37  CoordinateElement(int basix_element_handle, int geometric_dimension,
38  const std::string& signature,
41  std::function<int(int*, const uint32_t)> permute_dofs,
42  std::function<int(int*, const uint32_t)> unpermute_dofs);
43 
45  virtual ~CoordinateElement() = default;
46 
49  std::string signature() const;
50 
53  mesh::CellType cell_shape() const;
54 
56  int topological_dimension() const;
57 
59  int geometric_dimension() const;
60 
62  const ElementDofLayout& dof_layout() const;
63 
65  double non_affine_atol = 1.0e-8;
66 
69 
75  void push_forward(
76  Eigen::Ref<
77  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>>
78  x,
79  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
80  Eigen::Dynamic, Eigen::RowMajor>>& X,
81  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
82  Eigen::Dynamic, Eigen::RowMajor>>&
83  cell_geometry) const;
84 
88  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>& X,
89  Eigen::Tensor<double, 3, Eigen::RowMajor>& J, tcb::span<double> detJ,
90  Eigen::Tensor<double, 3, Eigen::RowMajor>& K,
91  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
92  Eigen::Dynamic, Eigen::RowMajor>>& x,
93  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
94  Eigen::Dynamic, Eigen::RowMajor>>&
95  cell_geometry) const;
96 
98  void permute_dofs(int* dofs, const uint32_t cell_perm) const;
99 
101  void unpermute_dofs(int* dofs, const uint32_t cell_perm) const;
102 
105  bool needs_permutation_data() const;
106 
107 private:
108  // Geometric dimensions
109  int _gdim;
110 
111  // Signature, usually from UFC
112  std::string _signature;
113 
114  // Layout of dofs on element
115  ElementDofLayout _dof_layout;
116 
117  // Flag denoting affine map
118  bool _is_affine;
119 
120  // Basix element
121  int _basix_element_handle;
122 
123  // Does the element need permutation data
124  bool _needs_permutation_data;
125 
126  // Dof permutation maker
127  std::function<int(int*, const uint32_t)> _permute_dofs;
128 
129  // Dof permutation maker
130  std::function<int(int*, const uint32_t)> _unpermute_dofs;
131 };
132 } // namespace dolfinx::fem
This class manages coordinate mappings for isoparametric cells.
Definition: CoordinateElement.h:26
virtual ~CoordinateElement()=default
Destructor.
int topological_dimension() const
Return the topological dimension of the cell shape.
Definition: CoordinateElement.cpp:54
double non_affine_atol
Absolute increment stopping criterium for non-affine Newton solver.
Definition: CoordinateElement.h:65
CoordinateElement(int basix_element_handle, int geometric_dimension, const std::string &signature, const ElementDofLayout &dof_layout, bool needs_permutation_data, std::function< int(int *, const uint32_t)> permute_dofs, std::function< int(int *, const uint32_t)> unpermute_dofs)
Create a coordinate element.
Definition: CoordinateElement.cpp:15
void push_forward(Eigen::Ref< Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> x, const Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &X, const Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &cell_geometry) const
Compute physical coordinates x for points X in the reference configuration.
Definition: CoordinateElement.cpp:67
int non_affine_max_its
Maximum number of iterations for non-affine Newton solver.
Definition: CoordinateElement.h:68
std::string signature() const
String identifying the finite element.
Definition: CoordinateElement.cpp:34
int geometric_dimension() const
Return the geometric dimension of the cell shape.
Definition: CoordinateElement.cpp:60
bool needs_permutation_data() const
Indicates whether the coordinate map needs permutation data passing in (for higher order geometries)
Definition: CoordinateElement.cpp:242
void compute_reference_geometry(Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &X, Eigen::Tensor< double, 3, Eigen::RowMajor > &J, tcb::span< double > detJ, Eigen::Tensor< double, 3, Eigen::RowMajor > &K, const Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &x, const Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &cell_geometry) const
Compute reference coordinates X, and J, detJ and K for physical coordinates x.
Definition: CoordinateElement.cpp:85
void permute_dofs(int *dofs, const uint32_t cell_perm) const
Permutes a list of DOF numbers on a cell.
Definition: CoordinateElement.cpp:231
const ElementDofLayout & dof_layout() const
Return the dof layout.
Definition: CoordinateElement.cpp:62
void unpermute_dofs(int *dofs, const uint32_t cell_perm) const
Reverses a DOF permutation.
Definition: CoordinateElement.cpp:236
mesh::CellType cell_shape() const
Cell shape.
Definition: CoordinateElement.cpp:36
The class represents the degree-of-freedom (dofs) for an element. Dofs are associated with a mesh ent...
Definition: ElementDofLayout.h:37
Finite element method functionality.
Definition: assemble_matrix_impl.h:24
CellType
Cell type identifier.
Definition: cell_types.h:21