DOLFIN-X
DOLFIN-X C++ interface
CoordinateElement.h
1 // Copyright (C) 2018 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 "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:
38  int geometric_dimension, const std::string& signature,
40  std::function<void(double*, int, const double*, const double*)>
41  compute_physical_coordinates,
42  std::function<void(double*, double*, double*, double*, int, const double*,
43  const double*)>
45 
47  virtual ~CoordinateElement() = default;
48 
51  std::string signature() const;
52 
55  mesh::CellType cell_shape() const;
56 
58  int topological_dimension() const;
59 
61  int geometric_dimension() const;
62 
64  const ElementDofLayout& dof_layout() const;
65 
71  void push_forward(
72  Eigen::Ref<
73  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>>
74  x,
75  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
76  Eigen::Dynamic, Eigen::RowMajor>>& X,
77  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
78  Eigen::Dynamic, Eigen::RowMajor>>&
79  cell_geometry) const;
80 
84  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>& X,
85  Eigen::Tensor<double, 3, Eigen::RowMajor>& J,
86  Eigen::Ref<Eigen::Array<double, Eigen::Dynamic, 1>> detJ,
87  Eigen::Tensor<double, 3, Eigen::RowMajor>& K,
88  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
89  Eigen::Dynamic, Eigen::RowMajor>>& x,
90  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
91  Eigen::Dynamic, Eigen::RowMajor>>&
92  cell_geometry) const;
93 
94 private:
95  int _tdim, _gdim;
96 
97  mesh::CellType _cell;
98 
99  std::string _signature;
100 
101  ElementDofLayout _dof_layout;
102 
103  std::function<void(double*, int, const double*, const double*)>
104  _compute_physical_coordinates;
105 
106  std::function<void(double*, double*, double*, double*, int, const double*,
107  const double*)>
108  _compute_reference_geometry;
109 };
110 } // namespace dolfinx::fem
dolfinx::mesh::CellType
CellType
Cell type identifier.
Definition: cell_types.h:22
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:43
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:38
dolfinx::fem::CoordinateElement::geometric_dimension
int geometric_dimension() const
Return the geometric dimension of the cell shape.
Definition: CoordinateElement.cpp:36
dolfinx::fem::CoordinateElement::CoordinateElement
CoordinateElement(mesh::CellType cell_type, int topological_dimension, int geometric_dimension, const std::string &signature, const ElementDofLayout &dof_layout, std::function< void(double *, int, const double *, const double *)> compute_physical_coordinates, std::function< void(double *, double *, double *, double *, int, const double *, const double *)> compute_reference_geometry)
Create a coordinate element.
Definition: CoordinateElement.cpp:14
dolfinx::fem::ElementDofLayout
The class represents the degree-of-freedom (dofs) for an element. Dofs are associated with a mesh ent...
Definition: ElementDofLayout.h:36
dolfinx::fem
Finite element method functionality.
Definition: assemble_matrix_impl.h:34
dolfinx::fem::CoordinateElement::signature
std::string signature() const
String identifying the finite element.
Definition: CoordinateElement.cpp:30
dolfinx::fem::CoordinateElement::topological_dimension
int topological_dimension() const
Return the topological dimension of the cell shape.
Definition: CoordinateElement.cpp:34
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:60
dolfinx::fem::CoordinateElement::cell_shape
mesh::CellType cell_shape() const
Cell shape.
Definition: CoordinateElement.cpp:32
dolfinx::fem::CoordinateElement
This class manages coordinate mappings for isoparametric cells.
Definition: CoordinateElement.h:23