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/Dense>
11 #include <cstdint>
12 #include <dolfinx/mesh/cell_types.h>
13 #include <functional>
14 #include <string>
15 #include <unsupported/Eigen/CXX11/Tensor>
16 
17 namespace dolfinx::fem
18 {
19 
20 // FIXME: A dof layout on a reference cell needs to be defined.
22 
24 {
25 public:
35  int geometric_dimension, const std::string& signature,
36  const ElementDofLayout& dof_layout, bool is_affine,
37  std::function<int(double*, int, int, const double*)>
38  evaluate_basis_derivatives);
39 
41  virtual ~CoordinateElement() = default;
42 
45  std::string signature() const;
46 
49  mesh::CellType cell_shape() const;
50 
52  int topological_dimension() const;
53 
55  int geometric_dimension() const;
56 
58  const ElementDofLayout& dof_layout() const;
59 
65  void push_forward(
66  Eigen::Ref<
67  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>>
68  x,
69  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
70  Eigen::Dynamic, Eigen::RowMajor>>& X,
71  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
72  Eigen::Dynamic, Eigen::RowMajor>>&
73  cell_geometry) const;
74 
78  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>& X,
79  Eigen::Tensor<double, 3, Eigen::RowMajor>& J,
80  Eigen::Ref<Eigen::Array<double, Eigen::Dynamic, 1>> detJ,
81  Eigen::Tensor<double, 3, Eigen::RowMajor>& K,
82  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
83  Eigen::Dynamic, Eigen::RowMajor>>& x,
84  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
85  Eigen::Dynamic, Eigen::RowMajor>>&
86  cell_geometry) const;
87 
88 private:
89  // Topological and geometric dimensions
90  int _tdim, _gdim;
91 
92  // Cell type
93  mesh::CellType _cell;
94 
95  // Signature, usually from UFC
96  std::string _signature;
97 
98  // Layout of dofs on element
99  ElementDofLayout _dof_layout;
100 
101  // Flag denoting affine map
102  bool _is_affine;
103 
104  // Function to evaluate the basis on the underlying element
105  // @param basis_values Returned values
106  // @param order
107  // @param num_points
108  // @param reference points
109  std::function<int(double*, int, int, const double*)>
110  _evaluate_basis_derivatives;
111 };
112 } // namespace dolfinx::fem
dolfinx::mesh::CellType
CellType
Cell type identifier.
Definition: cell_types.h:23
dolfinx::fem::CoordinateElement::push_forward
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:40
dolfinx::fem::CoordinateElement::CoordinateElement
CoordinateElement(mesh::CellType cell_type, int topological_dimension, int geometric_dimension, const std::string &signature, const ElementDofLayout &dof_layout, bool is_affine, std::function< int(double *, int, int, const double *)> evaluate_basis_derivatives)
Create a coordinate element.
Definition: CoordinateElement.cpp:14
dolfinx::fem::CoordinateElement::~CoordinateElement
virtual ~CoordinateElement()=default
Destructor.
dolfinx::fem::CoordinateElement::dof_layout
const ElementDofLayout & dof_layout() const
Return the dof layout.
Definition: CoordinateElement.cpp:35
dolfinx::fem::CoordinateElement::geometric_dimension
int geometric_dimension() const
Return the geometric dimension of the cell shape.
Definition: CoordinateElement.cpp:33
dolfinx::fem::ElementDofLayout
The class represents the degree-of-freedom (dofs) for an element. Dofs are associated with a mesh ent...
Definition: ElementDofLayout.h:37
dolfinx::fem
Finite element method functionality.
Definition: assemble_matrix_impl.h:23
dolfinx::fem::CoordinateElement::signature
std::string signature() const
String identifying the finite element.
Definition: CoordinateElement.cpp:27
dolfinx::fem::CoordinateElement::topological_dimension
int topological_dimension() const
Return the topological dimension of the cell shape.
Definition: CoordinateElement.cpp:31
dolfinx::fem::CoordinateElement::compute_reference_geometry
void compute_reference_geometry(Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &X, Eigen::Tensor< double, 3, Eigen::RowMajor > &J, Eigen::Ref< Eigen::Array< double, Eigen::Dynamic, 1 >> 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:61
dolfinx::fem::CoordinateElement::cell_shape
mesh::CellType cell_shape() const
Cell shape.
Definition: CoordinateElement.cpp:29
dolfinx::fem::CoordinateElement
This class manages coordinate mappings for isoparametric cells.
Definition: CoordinateElement.h:24