11 #include "element-families.h"
17 #include <xtensor/xadapt.hpp>
18 #include <xtensor/xcomplex.hpp>
19 #include <xtensor/xtensor.hpp>
20 #include <xtensor/xview.hpp>
183 const std::vector<std::vector<xt::xtensor<double, 3>>>& M,
184 const std::vector<std::vector<xt::xtensor<double, 2>>>& x,
int degree,
185 double kappa_tol = 0.0);
211 const xt::xtensor<double, 3>& coeffs,
212 const std::vector<xt::xtensor<double, 2>>& entity_transformations,
213 const std::array<std::vector<xt::xtensor<double, 2>>, 4>& x,
214 const std::array<std::vector<xt::xtensor<double, 3>>, 4>& M,
249 xt::xtensor<double, 4>
tabulate(
int nd,
const xt::xarray<double>& x)
const;
255 void tabulate(
int nd,
const xt::xarray<double>& x,
double* basis_data)
const;
281 element::family
family()
const;
301 xt::xtensor<double, 3>
303 const xt::xtensor<double, 3>& J,
304 const tcb::span<const double>& detJ,
305 const xt::xtensor<double, 3>& K)
const;
313 template <
typename T,
typename E>
315 const xt::xtensor<double, 3>& J,
316 const tcb::span<const double>& detJ,
317 const xt::xtensor<double, 3>& K, E&& u)
const;
325 xt::xtensor<double, 3>
map_pull_back(
const xt::xtensor<double, 3>& u,
326 const xt::xtensor<double, 3>& J,
327 const tcb::span<const double>& detJ,
328 const xt::xtensor<double, 3>& K)
const;
337 template <
typename T,
typename E>
339 const xt::xtensor<double, 3>& J,
340 const tcb::span<const double>& detJ,
341 const xt::xtensor<double, 3>& K, E&& U)
const;
350 const std::vector<std::vector<int>>&
entity_dofs()
const;
436 std::uint32_t cell_info)
const;
442 std::uint32_t cell_info)
const;
448 const xt::xtensor<double, 2>&
points()
const;
469 std::size_t _cell_tdim;
472 element::family _family;
478 std::vector<int> _value_shape;
488 xt::xtensor<double, 2> _coeffs;
491 std::vector<int> _cell_sub_entity_count;
499 std::vector<std::vector<int>> _entity_dofs;
502 std::vector<xt::xtensor<double, 2>> _entity_transformations;
509 xt::xtensor<double, 2> _points;
513 std::array<std::vector<xt::xtensor<double, 2>>, 4> _x;
516 xt::xtensor<double, 2> _matM;
519 std::array<std::vector<xt::xtensor<double, 3>>, 4> _matM_new;
522 bool _dof_transformations_are_permutations;
525 bool _dof_transformations_are_identity;
529 std::vector<std::vector<int>> _entity_permutations;
533 std::vector<std::vector<int>> _reverse_entity_permutations;
548 template <
typename T,
typename E>
550 const xt::xtensor<double, 3>& J,
551 const tcb::span<const double>& detJ,
552 const xt::xtensor<double, 3>& K,
561 for (std::size_t p = 0; p < U.shape(0); ++p)
563 auto J_p = xt::view(J, p, xt::all(), xt::all());
564 auto K_p = xt::view(K, p, xt::all(), xt::all());
567 for (std::size_t i = 0; i < U.shape(1); ++i)
569 auto U_data = xt::view(U, p, i, xt::all());
570 auto u_data = xt::view(u, p, i, xt::all());
576 template <
typename T,
typename E>
578 const xt::xtensor<double, 3>& J,
579 const tcb::span<const double>& detJ,
580 const xt::xtensor<double, 3>& K,
584 for (std::size_t p = 0; p < u.shape(0); ++p)
586 auto J_p = xt::view(J, p, xt::all(), xt::all());
587 auto K_p = xt::view(K, p, xt::all(), xt::all());
590 for (std::size_t i = 0; i < u.shape(1); ++i)
592 auto u_data = xt::view(u, p, i, xt::all());
593 auto U_data = xt::view(U, p, i, xt::all());
Definition: finite-element.h:191
FiniteElement & operator=(const FiniteElement &element)=default
Assignment operator.
cell::type cell_type() const
Definition: finite-element.cpp:415
FiniteElement(element::family family, cell::type cell_type, int degree, const std::vector< std::size_t > &value_shape, const xt::xtensor< double, 3 > &coeffs, const std::vector< xt::xtensor< double, 2 >> &entity_transformations, const std::array< std::vector< xt::xtensor< double, 2 >>, 4 > &x, const std::array< std::vector< xt::xtensor< double, 3 >>, 4 > &M, maps::type map_type=maps::type::identity)
Definition: finite-element.cpp:250
FiniteElement & operator=(FiniteElement &&element)=default
Move assignment operator.
xt::xtensor< double, 4 > tabulate(int nd, const xt::xarray< double > &x) const
Definition: finite-element.cpp:457
xt::xtensor< double, 3 > map_pull_back(const xt::xtensor< double, 3 > &u, const xt::xtensor< double, 3 > &J, const tcb::span< const double > &detJ, const xt::xtensor< double, 3 > &K) const
Definition: finite-element.cpp:608
int value_size() const
Definition: finite-element.cpp:419
FiniteElement(FiniteElement &&element)=default
Move constructor.
bool dof_transformations_are_permutations() const
Definition: finite-element.cpp:436
void map_pull_back_m(const xt::xtensor< T, 3 > &u, const xt::xtensor< double, 3 > &J, const tcb::span< const double > &detJ, const xt::xtensor< double, 3 > &K, E &&U) const
Definition: finite-element.h:577
xt::xtensor< double, 3 > map_push_forward(const xt::xtensor< double, 3 > &U, const xt::xtensor< double, 3 > &J, const tcb::span< const double > &detJ, const xt::xtensor< double, 3 > &K) const
Definition: finite-element.cpp:596
const xt::xtensor< double, 2 > & points() const
Definition: finite-element.cpp:594
void unpermute_dofs(tcb::span< std::int32_t > &dofs, std::uint32_t cell_info) const
Definition: finite-element.cpp:684
const xt::xtensor< double, 2 > & interpolation_matrix() const
Definition: finite-element.cpp:446
int degree() const
Definition: finite-element.cpp:417
const std::vector< int > & value_shape() const
Definition: finite-element.cpp:425
~FiniteElement()=default
Destructor.
int num_points() const
Return the number of interpolation points.
Definition: finite-element.cpp:592
const std::vector< std::vector< int > > & entity_dofs() const
Definition: finite-element.cpp:451
int dim() const
Definition: finite-element.cpp:430
element::family family() const
Definition: finite-element.cpp:432
maps::type map_type
Element map type.
Definition: finite-element.h:462
void map_push_forward_m(const xt::xtensor< T, 3 > &U, const xt::xtensor< double, 3 > &J, const tcb::span< const double > &detJ, const xt::xtensor< double, 3 > &K, E &&u) const
Definition: finite-element.h:549
bool dof_transformations_are_identity() const
Definition: finite-element.cpp:441
FiniteElement(const FiniteElement &element)=default
Copy constructor.
maps::type mapping_type() const
Definition: finite-element.cpp:434
void permute_dofs(tcb::span< std::int32_t > &dofs, std::uint32_t cell_info) const
Definition: finite-element.cpp:618
xt::xtensor< double, 3 > base_transformations() const
Definition: finite-element.cpp:539
type
Cell type.
Definition: cell.h:21
void apply_map(O &&u, const P &U, const Mat0 &J, double detJ, const Mat1 &K, maps::type map_type)
Definition: maps.h:121
type
Cell type.
Definition: maps.h:24
Placeholder.
Definition: brezzi-douglas-marini.h:10
const char * cell_type(int handle)
Definition: basix.cpp:189
std::string version()
Definition: finite-element.cpp:750
FiniteElement create_element(std::string family, std::string cell, int degree)
Create an element by name.
Definition: finite-element.cpp:75
xt::xtensor< double, 3 > compute_expansion_coefficients(cell::type cell_type, const xt::xtensor< double, 2 > &B, const std::vector< std::vector< xt::xtensor< double, 3 >>> &M, const std::vector< std::vector< xt::xtensor< double, 2 >>> &x, int degree, double kappa_tol=0.0)
Definition: finite-element.cpp:150
int degree(int handle)
Degree.
Definition: basix.cpp:195