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 
18 namespace dolfinx::fem
19 {
20 
21 // FIXME: A dof layout on a reference cell needs to be defined.
23 
25 {
26 public:
36  CoordinateElement(int basix_element_handle, int geometric_dimension,
37  const std::string& signature,
40  std::function<int(int*, const uint32_t)> permute_dofs,
41  std::function<int(int*, const uint32_t)> unpermute_dofs);
42 
44  virtual ~CoordinateElement() = default;
45 
48  std::string signature() const;
49 
52  mesh::CellType cell_shape() const;
53 
55  int topological_dimension() const;
56 
58  int geometric_dimension() const;
59 
61  const ElementDofLayout& dof_layout() const;
62 
64  double non_affine_atol = 1.0e-8;
65 
68 
75  const array2d<double>& cell_geometry) const;
76 
79  void compute_reference_geometry(array2d<double>& X, std::vector<double>& J,
80  tcb::span<double> detJ,
81  std::vector<double>& K,
82  const array2d<double>& x,
83  const array2d<double>& cell_geometry) const;
84 
86  void permute_dofs(int* dofs, const uint32_t cell_perm) const;
87 
89  void unpermute_dofs(int* dofs, const uint32_t cell_perm) const;
90 
93  bool needs_permutation_data() const;
94 
95 private:
96  // Geometric dimensions
97  int _gdim;
98 
99  // Signature, usually from UFC
100  std::string _signature;
101 
102  // Layout of dofs on element
103  ElementDofLayout _dof_layout;
104 
105  // Flag denoting affine map
106  bool _is_affine;
107 
108  // Basix element
109  int _basix_element_handle;
110 
111  // Does the element need permutation data
112  bool _needs_permutation_data;
113 
114  // Dof permutation maker
115  std::function<int(int*, const uint32_t)> _permute_dofs;
116 
117  // Dof permutation maker
118  std::function<int(int*, const uint32_t)> _unpermute_dofs;
119 };
120 } // 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:25
virtual ~CoordinateElement()=default
Destructor.
void push_forward(array2d< double > &x, const array2d< double > &X, const array2d< double > &cell_geometry) const
Compute physical coordinates x for points X in the reference configuration.
Definition: CoordinateElement.cpp:67
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:64
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 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:89
int non_affine_max_its
Maximum number of iterations for non-affine Newton solver.
Definition: CoordinateElement.h:67
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:247
void permute_dofs(int *dofs, const uint32_t cell_perm) const
Permutes a list of DOF numbers on a cell.
Definition: CoordinateElement.cpp:236
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:241
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:36
Finite element method functionality.
Definition: assemble_matrix_impl.h:23
CellType
Cell type identifier.
Definition: cell_types.h:22