8#include "element-families.h"
11#include "precompute.h"
12#include "sobolev-spaces.h"
30 = std::experimental::mdspan<double,
31 std::experimental::dextents<std::size_t, 2>>;
33 = std::experimental::mdspan<double,
34 std::experimental::dextents<std::size_t, 4>>;
36 = std::experimental::mdspan<
const double,
37 std::experimental::dextents<std::size_t, 2>>;
39 = std::experimental::mdspan<
const double,
40 std::experimental::dextents<std::size_t, 3>>;
42 = std::experimental::mdspan<
const double,
43 std::experimental::dextents<std::size_t, 4>>;
46 = std::experimental::mdarray<double,
47 std::experimental::dextents<std::size_t, 2>>;
49 = std::experimental::mdarray<double,
50 std::experimental::dextents<std::size_t, 4>>;
54inline std::array<std::vector<cmdspan2_t>, 4>
55to_mdspan(std::array<std::vector<mdarray2_t>, 4>& x)
57 std::array<std::vector<cmdspan2_t>, 4> x1;
58 for (std::size_t i = 0; i < x.size(); ++i)
59 for (std::size_t j = 0; j < x[i].size(); ++j)
60 x1[i].emplace_back(x[i][j].data(), x[i][j].extents());
67inline std::array<std::vector<cmdspan4_t>, 4>
68to_mdspan(std::array<std::vector<mdarray4_t>, 4>& M)
70 std::array<std::vector<cmdspan4_t>, 4> M1;
71 for (std::size_t i = 0; i < M.size(); ++i)
72 for (std::size_t j = 0; j < M[i].size(); ++j)
73 M1[i].emplace_back(M[i][j].data(), M[i][j].extents());
80inline std::array<std::vector<cmdspan2_t>, 4>
81to_mdspan(
const std::array<std::vector<std::vector<double>>, 4>& x,
82 const std::array<std::vector<std::array<std::size_t, 2>>, 4>& shape)
84 std::array<std::vector<cmdspan2_t>, 4> x1;
85 for (std::size_t i = 0; i < x.size(); ++i)
86 for (std::size_t j = 0; j < x[i].size(); ++j)
87 x1[i].push_back(cmdspan2_t(x[i][j].data(), shape[i][j]));
94inline std::array<std::vector<cmdspan4_t>, 4>
95to_mdspan(
const std::array<std::vector<std::vector<double>>, 4>& M,
96 const std::array<std::vector<std::array<std::size_t, 4>>, 4>& shape)
98 std::array<std::vector<cmdspan4_t>, 4> M1;
99 for (std::size_t i = 0; i < M.size(); ++i)
100 for (std::size_t j = 0; j < M[i].size(); ++j)
101 M1[i].push_back(cmdspan4_t(M[i][j].data(), shape[i][j]));
130std::tuple<std::array<std::vector<std::vector<double>>, 4>,
131 std::array<std::vector<std::array<std::size_t, 2>>, 4>,
132 std::array<std::vector<std::vector<double>>, 4>,
133 std::array<std::vector<std::array<std::size_t, 4>>, 4>>
135 const std::array<std::vector<cmdspan4_t>, 4>& M,
136 std::size_t tdim, std::size_t value_size);
148 = std::experimental::mdspan<
const double,
149 std::experimental::dextents<std::size_t, 2>>;
151 = std::experimental::mdspan<
const double,
152 std::experimental::dextents<std::size_t, 4>>;
328 const std::array<std::vector<cmdspan2_t>, 4>&
x,
329 const std::array<std::vector<cmdspan4_t>, 4>&
M,
334 std::vector<std::tuple<std::vector<FiniteElement>, std::vector<int>>>
342 const std::array<std::vector<cmdspan2_t>, 4>&
x,
343 const std::array<std::vector<cmdspan4_t>, 4>&
M,
348 std::vector<std::tuple<std::vector<FiniteElement>, std::vector<int>>>
356 const std::array<std::vector<cmdspan2_t>, 4>&
x,
357 const std::array<std::vector<cmdspan4_t>, 4>&
M,
362 std::vector<std::tuple<std::vector<FiniteElement>, std::vector<int>>>
370 const std::array<std::vector<cmdspan2_t>, 4>&
x,
371 const std::array<std::vector<cmdspan4_t>, 4>&
M,
375 std::vector<std::tuple<std::vector<FiniteElement>, std::vector<int>>>
410 std::size_t num_points)
const;
433 std::pair<std::vector<double>, std::array<std::size_t, 4>>
434 tabulate(
int nd, impl::cmdspan2_t
x)
const;
459 std::pair<std::vector<double>, std::array<std::size_t, 4>>
460 tabulate(
int nd,
const std::span<const double>&
x,
461 std::array<std::size_t, 2> shape)
const;
488 void tabulate(
int nd, impl::cmdspan2_t
x, impl::mdspan4_t basis)
const;
515 void tabulate(
int nd,
const std::span<const double>&
x,
516 std::array<std::size_t, 2> xshape,
517 const std::span<double>& basis)
const;
541 const std::vector<std::size_t>&
value_shape()
const;
594 std::pair<std::vector<double>, std::array<std::size_t, 3>>
596 std::span<const double> detJ, impl::cmdspan3_t K)
const;
605 std::pair<std::vector<double>, std::array<std::size_t, 3>>
606 pull_back(impl::cmdspan3_t u, impl::cmdspan3_t J,
607 std::span<const double> detJ, impl::cmdspan3_t K)
const;
640 template <
typename O,
typename P,
typename Q,
typename R>
641 std::function<void(O&,
const P&,
const Q&,
double,
const R&)>
map_fn()
const
645 case maps::type::identity:
646 return [](O& u,
const P& U,
const Q&, double,
const R&)
648 assert(U.extent(0) == u.extent(0));
649 assert(U.extent(1) == u.extent(1));
650 for (std::size_t i = 0; i < U.extent(0); ++i)
651 for (std::size_t j = 0; j < U.extent(1); ++j)
654 case maps::type::covariantPiola:
655 return [](O& u,
const P& U,
const Q& J,
double detJ,
const R& K)
657 case maps::type::contravariantPiola:
658 return [](O& u,
const P& U,
const Q& J,
double detJ,
const R& K)
660 case maps::type::doubleCovariantPiola:
661 return [](O& u,
const P& U,
const Q& J,
double detJ,
const R& K)
663 case maps::type::doubleContravariantPiola:
664 return [](O& u,
const P& U,
const Q& J,
double detJ,
const R& K)
667 throw std::runtime_error(
"Map not implemented");
677 const std::vector<std::vector<std::vector<int>>>&
entity_dofs()
const;
769 std::pair<std::vector<double>, std::array<std::size_t, 3>>
777 std::pair<std::vector<double>, std::array<std::size_t, 3>>>
788 std::uint32_t cell_info)
const;
798 std::uint32_t cell_info)
const;
808 template <
typename T>
810 std::uint32_t cell_info)
const;
820 template <
typename T>
823 std::uint32_t cell_info)
const;
833 template <
typename T>
835 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const;
845 template <
typename T>
848 std::uint32_t cell_info)
const;
858 template <
typename T>
861 std::uint32_t cell_info)
const;
871 template <
typename T>
873 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const;
884 template <
typename T>
886 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const;
896 template <
typename T>
898 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const;
904 const std::pair<std::vector<double>, std::array<std::size_t, 2>>&
959 const std::pair<std::vector<double>, std::array<std::size_t, 2>>&
967 const std::pair<std::vector<double>, std::array<std::size_t, 2>>&
1004 const std::pair<std::vector<double>, std::array<std::size_t, 2>>&
1012 std::vector<std::pair<std::vector<double>, std::array<std::size_t, 2>>>,
1053 std::vector<std::pair<std::vector<double>, std::array<std::size_t, 4>>>,
1062 const std::pair<std::vector<double>, std::array<std::size_t, 2>>&
1085 std::vector<std::tuple<std::vector<FiniteElement>, std::vector<int>>>
1100 std::size_t _cell_tdim;
1103 std::vector<std::vector<cell::type>> _cell_subentity_types;
1118 int _interpolation_nderivs;
1121 int _highest_degree;
1124 int _highest_complete_degree;
1127 std::vector<std::size_t> _value_shape;
1140 std::pair<std::vector<double>, std::array<std::size_t, 2>> _coeffs;
1143 std::vector<std::vector<std::vector<int>>> _edofs;
1146 std::vector<std::vector<std::vector<int>>> _e_closure_dofs;
1148 using array2_t = std::pair<std::vector<double>, std::array<std::size_t, 2>>;
1149 using array3_t = std::pair<std::vector<double>, std::array<std::size_t, 3>>;
1152 std::map<cell::type, array3_t> _entity_transformations;
1159 std::pair<std::vector<double>, std::array<std::size_t, 2>> _points;
1164 std::vector<std::pair<std::vector<double>, std::array<std::size_t, 2>>>,
1169 std::pair<std::vector<double>, std::array<std::size_t, 2>> _matM;
1173 bool _dof_transformations_are_permutations;
1176 bool _dof_transformations_are_identity;
1181 std::map<cell::type, std::vector<std::vector<std::size_t>>> _eperm;
1186 std::map<cell::type, std::vector<std::vector<std::size_t>>> _eperm_rev;
1189 = std::vector<std::pair<std::vector<std::size_t>, array2_t>>;
1192 std::map<cell::type, trans_data_t> _etrans;
1195 std::map<cell::type, trans_data_t> _etransT;
1198 std::map<cell::type, trans_data_t> _etrans_inv;
1201 std::map<cell::type, trans_data_t> _etrans_invT;
1205 bool _discontinuous;
1208 std::pair<std::vector<double>, std::array<std::size_t, 2>> _dual_matrix;
1215 std::vector<std::tuple<std::vector<FiniteElement>, std::vector<int>>>
1219 bool _interpolation_is_identity;
1223 std::pair<std::vector<double>, std::array<std::size_t, 2>> _wcoeffs;
1227 = std::vector<std::pair<std::vector<double>, std::array<std::size_t, 4>>>;
1228 std::array<array4_t, 4> _M;
1260 const std::vector<std::size_t>& value_shape,
1261 const impl::cmdspan2_t& wcoeffs,
1262 const std::array<std::vector<impl::cmdspan2_t>, 4>& x,
1263 const std::array<std::vector<impl::cmdspan4_t>, 4>& M,
1264 int interpolation_nderivs,
maps::type map_type,
1266 int highest_complete_degree,
int highest_degree);
1279 bool discontinuous);
1306 bool discontinuous);
1319 int degree,
bool discontinuous);
1368template <
typename T>
1371 std::uint32_t cell_info)
const
1373 if (_dof_transformations_are_identity)
1376 if (_cell_tdim >= 2)
1380 int face_start = _cell_tdim == 3 ? 3 * _edofs[2].size() : 0;
1382 for (
auto& edofs0 : _edofs[0])
1383 dofstart += edofs0.size();
1387 auto& [v_size_t, matrix] = _etrans.at(cell::type::interval)[0];
1388 for (std::size_t e = 0; e < _edofs[1].size(); ++e)
1391 if (cell_info >> (face_start + e) & 1)
1394 std::span(v_size_t),
1395 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1398 dofstart += _edofs[1][e].size();
1402 if (_cell_tdim == 3)
1405 for (std::size_t f = 0; f < _edofs[2].size(); ++f)
1407 auto& trans = _etrans.at(_cell_subentity_types[2][f]);
1410 if (cell_info >> (3 * f) & 1)
1412 const auto& m = trans[1];
1413 const auto& v_size_t = std::get<0>(m);
1414 const auto& matrix = std::get<1>(m);
1416 std::span(v_size_t),
1417 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1422 for (std::uint32_t r = 0; r < (cell_info >> (3 * f + 1) & 3); ++r)
1424 const auto& m = trans[0];
1425 const auto& v_size_t = std::get<0>(m);
1426 const auto& matrix = std::get<1>(m);
1428 std::span(v_size_t),
1429 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1433 dofstart += _edofs[2][f].size();
1439template <
typename T>
1441 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const
1443 if (_dof_transformations_are_identity)
1446 if (_cell_tdim >= 2)
1450 int face_start = _cell_tdim == 3 ? 3 * _edofs[2].size() : 0;
1452 for (
auto& edofs0 : _edofs[0])
1453 dofstart += edofs0.size();
1457 auto& [v_size_t, matrix] = _etransT.at(cell::type::interval)[0];
1458 for (std::size_t e = 0; e < _edofs[1].size(); ++e)
1461 if (cell_info >> (face_start + e) & 1)
1464 std::span(v_size_t),
1465 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1468 dofstart += _edofs[1][e].size();
1472 if (_cell_tdim == 3)
1475 for (std::size_t f = 0; f < _edofs[2].size(); ++f)
1477 auto& trans = _etransT.at(_cell_subentity_types[2][f]);
1480 for (std::uint32_t r = 0; r < (cell_info >> (3 * f + 1) & 3); ++r)
1483 auto& v_size_t = std::get<0>(m);
1484 auto& matrix = std::get<1>(m);
1486 std::span(v_size_t),
1487 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1492 if (cell_info >> (3 * f) & 1)
1495 auto& v_size_t = std::get<0>(m);
1496 auto& matrix = std::get<1>(m);
1498 std::span(v_size_t),
1499 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1502 dofstart += _edofs[2][f].size();
1508template <
typename T>
1510 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const
1512 if (_dof_transformations_are_identity)
1515 if (_cell_tdim >= 2)
1519 int face_start = _cell_tdim == 3 ? 3 * _edofs[2].size() : 0;
1521 for (
auto& edofs0 : _edofs[0])
1522 dofstart += edofs0.size();
1525 for (std::size_t e = 0; e < _edofs[1].size(); ++e)
1528 if (cell_info >> (face_start + e) & 1)
1530 auto& [v_size_t, matrix] = _etrans_invT.at(cell::type::interval)[0];
1532 cmdspan2_t(matrix.first.data(), matrix.second),
1533 data, dofstart, block_size);
1535 dofstart += _edofs[1][e].size();
1538 if (_cell_tdim == 3)
1541 for (std::size_t f = 0; f < _edofs[2].size(); ++f)
1543 auto& trans = _etrans_invT.at(_cell_subentity_types[2][f]);
1546 if (cell_info >> (3 * f) & 1)
1549 auto& v_size_t = std::get<0>(m);
1550 auto& matrix = std::get<1>(m);
1552 std::span(v_size_t),
1553 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1558 for (std::uint32_t r = 0; r < (cell_info >> (3 * f + 1) & 3); ++r)
1560 const auto& m = trans[0];
1561 const auto& v_size_t = std::get<0>(m);
1562 const auto& matrix = std::get<1>(m);
1564 std::span(v_size_t),
1565 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1569 dofstart += _edofs[2][f].size();
1575template <
typename T>
1577 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const
1579 if (_dof_transformations_are_identity)
1582 if (_cell_tdim >= 2)
1586 int face_start = _cell_tdim == 3 ? 3 * _edofs[2].size() : 0;
1588 for (
auto& edofs0 : _edofs[0])
1589 dofstart += edofs0.size();
1593 auto& [v_size_t, matrix] = _etrans_inv.at(cell::type::interval)[0];
1594 for (std::size_t e = 0; e < _edofs[1].size(); ++e)
1597 if (cell_info >> (face_start + e) & 1)
1600 std::span(v_size_t),
1601 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1604 dofstart += _edofs[1][e].size();
1608 if (_cell_tdim == 3)
1611 for (std::size_t f = 0; f < _edofs[2].size(); ++f)
1613 auto& trans = _etrans_inv.at(_cell_subentity_types[2][f]);
1616 for (std::uint32_t r = 0; r < (cell_info >> (3 * f + 1) & 3); ++r)
1619 auto& v_size_t = std::get<0>(m);
1620 auto& matrix = std::get<1>(m);
1622 std::span(v_size_t),
1623 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1628 if (cell_info >> (3 * f) & 1)
1631 auto& v_size_t = std::get<0>(m);
1632 auto& matrix = std::get<1>(m);
1634 std::span(v_size_t),
1635 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1639 dofstart += _edofs[2][f].size();
1645template <
typename T>
1647 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const
1649 if (_dof_transformations_are_identity)
1652 if (_cell_tdim >= 2)
1656 int face_start = _cell_tdim == 3 ? 3 * _edofs[2].size() : 0;
1658 for (
auto& edofs0 : _edofs[0])
1659 dofstart += edofs0.size();
1663 auto& [v_size_t, matrix] = _etrans.at(cell::type::interval)[0];
1664 for (std::size_t e = 0; e < _edofs[1].size(); ++e)
1667 if (cell_info >> (face_start + e) & 1)
1670 std::span(v_size_t),
1671 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1675 dofstart += _edofs[1][e].size();
1679 if (_cell_tdim == 3)
1682 for (std::size_t f = 0; f < _edofs[2].size(); ++f)
1684 auto& trans = _etrans.at(_cell_subentity_types[2][f]);
1687 if (cell_info >> (3 * f) & 1)
1690 auto& v_size_t = std::get<0>(m);
1691 auto& matrix = std::get<1>(m);
1693 std::span(v_size_t),
1694 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1699 for (std::uint32_t r = 0; r < (cell_info >> (3 * f + 1) & 3); ++r)
1702 auto& v_size_t = std::get<0>(m);
1703 auto& matrix = std::get<1>(m);
1705 std::span(v_size_t),
1706 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1709 dofstart += _edofs[2][f].size();
1715template <
typename T>
1717 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const
1719 if (_dof_transformations_are_identity)
1722 if (_cell_tdim >= 2)
1726 int face_start = _cell_tdim == 3 ? 3 * _edofs[2].size() : 0;
1728 for (
auto& edofs0 : _edofs[0])
1729 dofstart += edofs0.size();
1733 auto& [v_size_t, matrix] = _etrans_invT.at(cell::type::interval)[0];
1734 for (std::size_t e = 0; e < _edofs[1].size(); ++e)
1737 if (cell_info >> (face_start + e) & 1)
1740 std::span(v_size_t),
1741 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1744 dofstart += _edofs[1][e].size();
1748 if (_cell_tdim == 3)
1751 for (std::size_t f = 0; f < _edofs[2].size(); ++f)
1753 auto& trans = _etrans_invT.at(_cell_subentity_types[2][f]);
1756 if (cell_info >> (3 * f) & 1)
1759 auto& v_size_t = std::get<0>(m);
1760 auto& matrix = std::get<1>(m);
1762 std::span(v_size_t),
1763 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1768 for (std::uint32_t r = 0; r < (cell_info >> (3 * f + 1) & 3); ++r)
1771 auto& v_size_t = std::get<0>(m);
1772 auto& matrix = std::get<1>(m);
1774 std::span(v_size_t),
1775 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1778 dofstart += _edofs[2][f].size();
1784template <
typename T>
1786 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const
1788 if (_dof_transformations_are_identity)
1791 if (_cell_tdim >= 2)
1795 int face_start = _cell_tdim == 3 ? 3 * _edofs[2].size() : 0;
1797 for (
auto& edofs0 : _edofs[0])
1798 dofstart += edofs0.size();
1802 auto& [v_size_t, matrix] = _etransT.at(cell::type::interval)[0];
1803 for (std::size_t e = 0; e < _edofs[1].size(); ++e)
1806 if (cell_info >> (face_start + e) & 1)
1809 std::span(v_size_t),
1810 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1813 dofstart += _edofs[1][e].size();
1817 if (_cell_tdim == 3)
1820 for (std::size_t f = 0; f < _edofs[2].size(); ++f)
1822 auto& trans = _etransT.at(_cell_subentity_types[2][f]);
1825 for (std::uint32_t r = 0; r < (cell_info >> (3 * f + 1) & 3); ++r)
1828 auto& v_size_t = std::get<0>(m);
1829 auto& matrix = std::get<1>(m);
1831 std::span(v_size_t),
1832 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1837 if (cell_info >> (3 * f) & 1)
1840 auto& v_size_t = std::get<0>(m);
1841 auto& matrix = std::get<1>(m);
1843 std::span(v_size_t),
1844 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1847 dofstart += _edofs[2][f].size();
1853template <
typename T>
1855 const std::span<T>& data,
int block_size, std::uint32_t cell_info)
const
1857 if (_dof_transformations_are_identity)
1860 if (_cell_tdim >= 2)
1864 int face_start = _cell_tdim == 3 ? 3 * _edofs[2].size() : 0;
1866 for (
auto& edofs0 : _edofs[0])
1867 dofstart += edofs0.size();
1871 auto& [v_size_t, matrix] = _etrans_inv.at(cell::type::interval)[0];
1872 for (std::size_t e = 0; e < _edofs[1].size(); ++e)
1875 if (cell_info >> (face_start + e) & 1)
1878 std::span(v_size_t),
1879 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1882 dofstart += _edofs[1][e].size();
1886 if (_cell_tdim == 3)
1889 for (std::size_t f = 0; f < _edofs[2].size(); ++f)
1891 auto& trans = _etrans_inv.at(_cell_subentity_types[2][f]);
1894 for (std::uint32_t r = 0; r < (cell_info >> (3 * f + 1) & 3); ++r)
1897 auto& v_size_t = std::get<0>(m);
1898 auto& matrix = std::get<1>(m);
1900 std::span(v_size_t),
1901 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1906 if (cell_info >> (3 * f) & 1)
1909 auto& v_size_t = std::get<0>(m);
1910 auto& matrix = std::get<1>(m);
1912 std::span(v_size_t),
1913 cmdspan2_t(matrix.first.data(), matrix.second), data, dofstart,
1917 dofstart += _edofs[2][f].size();
A finite element.
Definition: finite-element.h:146
std::map< cell::type, std::pair< std::vector< double >, std::array< std::size_t, 3 > > > entity_transformations() const
Definition: finite-element.cpp:1427
std::pair< std::vector< double >, std::array< std::size_t, 4 > > tabulate(int nd, impl::cmdspan2_t x) const
Compute basis values and derivatives at set of points.
Definition: finite-element.cpp:1049
void apply_inverse_dof_transformation_to_transpose(const std::span< T > &data, int block_size, std::uint32_t cell_info) const
Definition: finite-element.h:1854
cell::type cell_type() const
Definition: finite-element.cpp:1123
sobolev::space sobolev_space() const
Definition: finite-element.cpp:1145
int highest_complete_degree() const
Definition: finite-element.cpp:1129
FiniteElement & operator=(const FiniteElement &element)=default
Assignment operator.
element::lagrange_variant lagrange_variant() const
Definition: finite-element.cpp:1469
FiniteElement(FiniteElement &&element)=default
Move constructor.
bool operator==(const FiniteElement &e) const
Definition: finite-element.cpp:998
const std::array< std::vector< std::pair< std::vector< double >, std::array< std::size_t, 2 > > >, 4 > & x() const
Definition: finite-element.cpp:1446
void apply_dof_transformation(const std::span< T > &data, int block_size, std::uint32_t cell_info) const
Definition: finite-element.h:1369
bool dof_transformations_are_permutations() const
Definition: finite-element.cpp:1149
maps::type map_type() const
Definition: finite-element.cpp:1143
std::vector< std::tuple< std::vector< FiniteElement >, std::vector< int > > > get_tensor_product_representation() const
Definition: finite-element.cpp:1487
const std::vector< std::vector< std::vector< int > > > & entity_dofs() const
Definition: finite-element.cpp:1166
const std::pair< std::vector< double >, std::array< std::size_t, 2 > > & dual_matrix() const
Definition: finite-element.cpp:1433
void apply_inverse_transpose_dof_transformation_to_transpose(const std::span< T > &data, int block_size, std::uint32_t cell_info) const
Apply inverse transpose DOF transformations to some transposed data.
Definition: finite-element.h:1716
int interpolation_nderivs() const
The number of derivatives needed when interpolating.
Definition: finite-element.cpp:1481
FiniteElement & operator=(FiniteElement &&element)=default
Move assignment operator.
const std::vector< std::vector< std::vector< int > > > & entity_closure_dofs() const
Definition: finite-element.cpp:1172
const std::pair< std::vector< double >, std::array< std::size_t, 2 > > & points() const
Definition: finite-element.cpp:1246
int degree() const
Definition: finite-element.cpp:1125
void apply_transpose_dof_transformation(const std::span< T > &data, int block_size, std::uint32_t cell_info) const
Definition: finite-element.h:1440
const std::vector< std::size_t > & value_shape() const
Definition: finite-element.cpp:1134
const std::pair< std::vector< double >, std::array< std::size_t, 2 > > & coefficient_matrix() const
Definition: finite-element.cpp:1459
bool has_tensor_product_factorisation() const
Definition: finite-element.cpp:1464
~FiniteElement()=default
Destructor.
const std::pair< std::vector< double >, std::array< std::size_t, 2 > > & wcoeffs() const
Definition: finite-element.cpp:1439
int dim() const
Definition: finite-element.cpp:1139
int highest_degree() const
Definition: finite-element.cpp:1127
std::pair< std::vector< double >, std::array< std::size_t, 3 > > push_forward(impl::cmdspan3_t U, impl::cmdspan3_t J, std::span< const double > detJ, impl::cmdspan3_t K) const
Definition: finite-element.cpp:1252
std::array< std::size_t, 4 > tabulate_shape(std::size_t nd, std::size_t num_points) const
Definition: finite-element.cpp:1035
element::dpc_variant dpc_variant() const
Definition: finite-element.cpp:1474
void apply_transpose_dof_transformation_to_transpose(const std::span< T > &data, int block_size, std::uint32_t cell_info) const
Definition: finite-element.h:1785
void apply_dof_transformation_to_transpose(const std::span< T > &data, int block_size, std::uint32_t cell_info) const
Definition: finite-element.h:1646
element::family family() const
Definition: finite-element.cpp:1141
void permute_dofs(const std::span< std::int32_t > &dofs, std::uint32_t cell_info) const
Definition: finite-element.cpp:1319
bool dof_transformations_are_identity() const
Definition: finite-element.cpp:1154
void apply_inverse_dof_transformation(const std::span< T > &data, int block_size, std::uint32_t cell_info) const
Definition: finite-element.h:1576
std::pair< std::vector< double >, std::array< std::size_t, 3 > > pull_back(impl::cmdspan3_t u, impl::cmdspan3_t J, std::span< const double > detJ, impl::cmdspan3_t K) const
Definition: finite-element.cpp:1287
FiniteElement(const FiniteElement &element)=default
Copy constructor.
const std::pair< std::vector< double >, std::array< std::size_t, 2 > > & interpolation_matrix() const
Return a matrix of weights interpolation,.
Definition: finite-element.cpp:1160
FiniteElement(element::family family, cell::type cell_type, int degree, const std::vector< std::size_t > &value_shape, const cmdspan2_t &wcoeffs, const std::array< std::vector< cmdspan2_t >, 4 > &x, const std::array< std::vector< cmdspan4_t >, 4 > &M, int interpolation_nderivs, maps::type map_type, sobolev::space sobolev_space, bool discontinuous, int highest_complete_degree, int highest_degree, element::lagrange_variant lvariant, element::dpc_variant dvariant, std::vector< std::tuple< std::vector< FiniteElement >, std::vector< int > > > tensor_factors={})
Construct a finite element.
Definition: finite-element.cpp:606
bool interpolation_is_identity() const
Definition: finite-element.cpp:1476
std::pair< std::vector< double >, std::array< std::size_t, 3 > > base_transformations() const
Get the base transformations.
Definition: finite-element.cpp:1178
const std::array< std::vector< std::pair< std::vector< double >, std::array< std::size_t, 4 > > >, 4 > & M() const
Definition: finite-element.cpp:1453
std::function< void(O &, const P &, const Q &, double, const R &)> map_fn() const
Definition: finite-element.h:641
bool discontinuous() const
Definition: finite-element.cpp:1147
void apply_inverse_transpose_dof_transformation(const std::span< T > &data, int block_size, std::uint32_t cell_info) const
Definition: finite-element.h:1509
void unpermute_dofs(const std::span< std::int32_t > &dofs, std::uint32_t cell_info) const
Definition: finite-element.cpp:1373
type
Cell type.
Definition: cell.h:20
lagrange_variant
Variants of a Lagrange space that can be created.
Definition: element-families.h:12
dpc_variant
Definition: element-families.h:33
impl::cmdspan4_t cmdspan4_t
Typedef for mdspan.
Definition: finite-element.h:113
impl::cmdspan2_t cmdspan2_t
Typedef for mdspan.
Definition: finite-element.h:111
std::tuple< std::array< std::vector< std::vector< double > >, 4 >, std::array< std::vector< std::array< std::size_t, 2 > >, 4 >, std::array< std::vector< std::vector< double > >, 4 >, std::array< std::vector< std::array< std::size_t, 4 > >, 4 > > make_discontinuous(const std::array< std::vector< cmdspan2_t >, 4 > &x, const std::array< std::vector< cmdspan4_t >, 4 > &M, std::size_t tdim, std::size_t value_size)
Definition: finite-element.cpp:394
family
Available element families.
Definition: element-families.h:46
void covariant_piola(O &&r, const P &U, const Q &, double, const R &K)
Covariant Piola map.
Definition: maps.h:39
void contravariant_piola(O &&r, const P &U, const Q &J, double detJ, const R &)
Contravariant Piola map.
Definition: maps.h:58
void double_contravariant_piola(O &&r, const P &U, const Q &J, double detJ, const R &)
Double contravariant Piola map.
Definition: maps.h:107
void double_covariant_piola(O &&r, const P &U, const Q &J, double, const R &K)
Double covariant Piola map.
Definition: maps.h:79
type
Map type.
Definition: maps.h:17
void apply_matrix_to_transpose(const std::span< const std::size_t > &v_size_t, const std::experimental::mdspan< const T, std::experimental::dextents< std::size_t, 2 > > &M, const std::span< E > &data, std::size_t offset=0, std::size_t block_size=1)
Apply a (precomputed) matrix to some transposed data.
Definition: precompute.h:244
void apply_matrix(const std::span< const std::size_t > &v_size_t, const std::experimental::mdspan< const T, std::experimental::dextents< std::size_t, 2 > > &M, const std::span< E > &data, std::size_t offset=0, std::size_t block_size=1)
Apply a (precomputed) matrix.
Definition: precompute.h:207
space
Sobolev space type.
Definition: sobolev-spaces.h:13
Basix: FEniCS runtime basis evaluation library.
Definition: cell.h:16
FiniteElement create_element(element::family family, cell::type cell, int degree, element::lagrange_variant lvariant, bool discontinuous)
Definition: finite-element.cpp:345
std::string version()
Definition: finite-element.cpp:1494
FiniteElement create_custom_element(cell::type cell_type, const std::vector< std::size_t > &value_shape, const impl::cmdspan2_t &wcoeffs, const std::array< std::vector< impl::cmdspan2_t >, 4 > &x, const std::array< std::vector< impl::cmdspan4_t >, 4 > &M, int interpolation_nderivs, maps::type map_type, sobolev::space sobolev_space, bool discontinuous, int highest_complete_degree, int highest_degree)
Definition: finite-element.cpp:464