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 <cstdint>
11 #include <dolfinx/common/array2d.h>
12 #include <dolfinx/common/span.hpp>
13 #include <dolfinx/mesh/cell_types.h>
14 #include <functional>
15 #include <memory>
16 #include <string>
17 #include <xtensor/xtensor.hpp>
18 
19 namespace basix
20 {
21 class FiniteElement;
22 }
23 
24 namespace dolfinx::fem
25 {
26 
27 // FIXME: A dof layout on a reference cell needs to be defined.
29 
31 {
32 public:
38  CoordinateElement(std::shared_ptr<basix::FiniteElement> element,
39  int geometric_dimension, const std::string& signature,
41 
43  virtual ~CoordinateElement() = default;
44 
47  std::string signature() const;
48 
51  mesh::CellType cell_shape() const;
52 
54  int topological_dimension() const;
55 
57  int geometric_dimension() const;
58 
62  xt::xtensor<double, 4> tabulate(int n, const array2d<double>& X) const;
63 
66  void compute_jacobian_data(const xt::xtensor<double, 4>& tabulated_data,
67  const array2d<double>& X,
68  const array2d<double>& cell_geometry,
69  std::vector<double>& J, tcb::span<double> detJ,
70  std::vector<double>& K) const;
71 
73  const ElementDofLayout& dof_layout() const;
74 
76  double non_affine_atol = 1.0e-8;
77 
80 
86  void push_forward(array2d<double>& x, const array2d<double>& cell_geometry,
87  const xt::xtensor<double, 2>& phi) const;
88 
91  void compute_reference_geometry(array2d<double>& X, std::vector<double>& J,
92  tcb::span<double> detJ,
93  std::vector<double>& K,
94  const array2d<double>& x,
95  const array2d<double>& cell_geometry) const;
96 
98  void permute_dofs(std::int32_t* dofs, const uint32_t cell_perm) const;
99 
101  void unpermute_dofs(std::int32_t* 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  // TODO: This should be removed now as we are transitioning to
121  // basix::FiniteElement
122  // Basix element
123  int _basix_element_handle;
124 
125  // Basix Element (basix::FiniteElement)
126  std::shared_ptr<basix::FiniteElement> _element;
127 };
128 } // namespace dolfinx::fem
This class provides a dynamic 2-dimensional row-wise array data structure.
Definition: array2d.h:21
This class manages coordinate mappings for isoparametric cells.
Definition: CoordinateElement.h:31
virtual ~CoordinateElement()=default
Destructor.
xt::xtensor< double, 4 > tabulate(int n, const array2d< double > &X) const
Tabulate shape functions up to n-th order derivative at points X in the reference geometry Note: Dyna...
Definition: CoordinateElement.cpp:248
int topological_dimension() const
Return the topological dimension of the cell shape.
Definition: CoordinateElement.cpp:55
double non_affine_atol
Absolute increment stopping criterium for non-affine Newton solver.
Definition: CoordinateElement.h:76
void push_forward(array2d< double > &x, const array2d< double > &cell_geometry, const xt::xtensor< double, 2 > &phi) const
Compute physical coordinates x for points X in the reference configuration.
Definition: CoordinateElement.cpp:67
void compute_reference_geometry(array2d< double > &X, std::vector< double > &J, tcb::span< double > detJ, std::vector< double > &K, const array2d< double > &x, const array2d< double > &cell_geometry) const
Compute reference coordinates X, and J, detJ and K for physical coordinates x.
Definition: CoordinateElement.cpp:83
int non_affine_max_its
Maximum number of iterations for non-affine Newton solver.
Definition: CoordinateElement.h:79
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 permute_dofs(std::int32_t *dofs, const uint32_t cell_perm) const
Permutes a list of DOF numbers on a cell.
Definition: CoordinateElement.cpp:230
void unpermute_dofs(std::int32_t *dofs, const uint32_t cell_perm) const
Reverses a DOF permutation.
Definition: CoordinateElement.cpp:236
void compute_jacobian_data(const xt::xtensor< double, 4 > &tabulated_data, const array2d< double > &X, const array2d< double > &cell_geometry, std::vector< double > &J, tcb::span< double > detJ, std::vector< double > &K) const
Compute J, K and detJ for a cell with given geometry, and the basis functions and first order derivat...
Definition: CoordinateElement.cpp:255
const ElementDofLayout & dof_layout() const
Return the dof layout.
Definition: CoordinateElement.cpp:62
CoordinateElement(std::shared_ptr< basix::FiniteElement > element, int geometric_dimension, const std::string &signature, const ElementDofLayout &dof_layout)
Create a coordinate element.
Definition: CoordinateElement.cpp:17
std::string signature() const
String identifying the finite element.
mesh::CellType cell_shape() const
Cell shape.
Definition: CoordinateElement.cpp:37
The class represents the degree-of-freedom (dofs) for an element. Dofs are associated with a mesh ent...
Definition: ElementDofLayout.h:36
Finite element method functionality.
Definition: assemble_matrix_impl.h:23
CellType
Cell type identifier.
Definition: cell_types.h:23