11#include "FunctionSpace.h"
14#include <basix/mdspan.hpp>
15#include <dolfinx/common/IndexMap.h>
16#include <dolfinx/mesh/Geometry.h>
17#include <dolfinx/mesh/Mesh.h>
18#include <dolfinx/mesh/Topology.h>
22namespace dolfinx::fem::impl
25template <dolfinx::scalar T>
26T assemble_cells(mdspan2_t x_dofmap,
27 md::mdspan<
const scalar_value_t<T>,
28 md::extents<std::size_t, md::dynamic_extent, 3>>
30 std::span<const std::int32_t> cells, FEkernel<T>
auto fn,
31 std::span<const T> constants,
32 md::mdspan<
const T, md::dextents<std::size_t, 2>> coeffs)
39 std::vector<scalar_value_t<T>> cdofs(3 * x_dofmap.extent(1));
42 for (std::size_t index = 0; index <
cells.size(); ++index)
44 std::int32_t c =
cells[index];
47 auto x_dofs = md::submdspan(x_dofmap, c, md::full_extent);
48 for (std::size_t i = 0; i < x_dofs.size(); ++i)
50 std::copy_n(&x(x_dofs[i], 0), 3, std::next(cdofs.begin(), 3 * i));
53 fn(&value, &coeffs(index, 0), constants.data(), cdofs.data(),
nullptr,
61template <dolfinx::scalar T>
62T assemble_exterior_facets(
64 md::mdspan<
const scalar_value_t<T>,
65 md::extents<std::size_t, md::dynamic_extent, 3>>
67 md::mdspan<
const std::int32_t,
68 md::extents<std::size_t, md::dynamic_extent, 2>>
70 FEkernel<T>
auto fn, std::span<const T> constants,
71 md::mdspan<
const T, md::dextents<std::size_t, 2>> coeffs,
72 md::mdspan<
const std::uint8_t, md::dextents<std::size_t, 2>> perms)
79 std::vector<scalar_value_t<T>> cdofs(3 * x_dofmap.extent(1));
82 for (std::size_t f = 0; f < facets.extent(0); ++f)
84 std::int32_t
cell = facets(f, 0);
85 std::int32_t local_facet = facets(f, 1);
88 auto x_dofs = md::submdspan(x_dofmap,
cell, md::full_extent);
89 for (std::size_t i = 0; i < x_dofs.size(); ++i)
90 std::copy_n(&x(x_dofs[i], 0), 3, std::next(cdofs.begin(), 3 * i));
93 std::uint8_t perm = perms.empty() ? 0 : perms(
cell, local_facet);
94 fn(&value, &coeffs(f, 0), constants.data(), cdofs.data(), &local_facet,
102template <dolfinx::scalar T>
103T assemble_interior_facets(
105 md::mdspan<
const scalar_value_t<T>,
106 md::extents<std::size_t, md::dynamic_extent, 3>>
108 md::mdspan<
const std::int32_t,
109 md::extents<std::size_t, md::dynamic_extent, 2, 2>>
111 FEkernel<T>
auto fn, std::span<const T> constants,
112 md::mdspan<
const T, md::extents<std::size_t, md::dynamic_extent, 2,
115 md::mdspan<
const std::uint8_t, md::dextents<std::size_t, 2>> perms)
122 using X = scalar_value_t<T>;
123 std::vector<X> cdofs(2 * x_dofmap.extent(1) * 3);
124 std::span<X> cdofs0(cdofs.data(), x_dofmap.extent(1) * 3);
125 std::span<X> cdofs1(cdofs.data() + x_dofmap.extent(1) * 3,
126 x_dofmap.extent(1) * 3);
129 for (std::size_t f = 0; f < facets.extent(0); ++f)
131 std::array
cells = {facets(f, 0, 0), facets(f, 1, 0)};
132 std::array local_facet = {facets(f, 0, 1), facets(f, 1, 1)};
135 auto x_dofs0 = md::submdspan(x_dofmap, cells[0], md::full_extent);
136 for (std::size_t i = 0; i < x_dofs0.size(); ++i)
137 std::copy_n(&x(x_dofs0[i], 0), 3, std::next(cdofs0.begin(), 3 * i));
138 auto x_dofs1 = md::submdspan(x_dofmap, cells[1], md::full_extent);
139 for (std::size_t i = 0; i < x_dofs1.size(); ++i)
140 std::copy_n(&x(x_dofs1[i], 0), 3, std::next(cdofs1.begin(), 3 * i));
142 std::array perm = perms.empty()
143 ? std::array<std::uint8_t, 2>{0, 0}
144 : std::array{perms(cells[0], local_facet[0]),
145 perms(cells[1], local_facet[1])};
146 fn(&value, &coeffs(f, 0, 0), constants.data(), cdofs.data(),
147 local_facet.data(), perm.data(),
nullptr);
154template <dolfinx::scalar T, std::
floating_po
int U>
156 const fem::Form<T, U>& M, mdspan2_t x_dofmap,
157 md::mdspan<
const scalar_value_t<T>,
158 md::extents<std::size_t, md::dynamic_extent, 3>>
160 std::span<const T> constants,
161 const std::map<std::pair<IntegralType, int>,
162 std::pair<std::span<const T>,
int>>& coefficients)
164 std::shared_ptr<const mesh::Mesh<U>> mesh = M.mesh();
174 assert(
cells.size() * cstride == coeffs.size());
175 value += impl::assemble_cells(
176 x_dofmap, x, cells, fn, constants,
177 md::mdspan(coeffs.data(),
cells.size(), cstride));
181 int num_facets_per_cell
183 md::mdspan<const std::uint8_t, md::dextents<std::size_t, 2>> perms;
184 if (M.needs_facet_permutations())
186 mesh->topology_mutable()->create_entity_permutations();
187 const std::vector<std::uint8_t>& p
188 = mesh->topology()->get_facet_permutations();
189 perms = md::mdspan(p.data(), p.size() / num_facets_per_cell,
190 num_facets_per_cell);
197 auto& [coeffs, cstride]
201 assert((facets.size() / 2) * cstride == coeffs.size());
202 value += impl::assemble_exterior_facets(
204 md::mdspan<
const std::int32_t,
205 md::extents<std::size_t, md::dynamic_extent, 2>>(
206 facets.data(), facets.size() / 2, 2),
207 fn, constants, md::mdspan(coeffs.data(), facets.size() / 2, cstride),
215 auto& [coeffs, cstride]
218 assert((facets.size() / 4) * 2 * cstride == coeffs.size());
219 value += impl::assemble_interior_facets(
221 md::mdspan<
const std::int32_t,
222 md::extents<std::size_t, md::dynamic_extent, 2, 2>>(
223 facets.data(), facets.size() / 4, 2, 2),
225 md::mdspan<
const T, md::extents<std::size_t, md::dynamic_extent, 2,
226 md::dynamic_extent>>(
227 coeffs.data(), facets.size() / 4, 2, cstride),
Functions supporting finite element method operations.
void cells(la::SparsityPattern &pattern, std::array< std::span< const std::int32_t >, 2 > cells, std::array< std::reference_wrapper< const DofMap >, 2 > dofmaps)
Iterate over cells and insert entries into sparsity pattern.
Definition sparsitybuild.cpp:16
@ interior_facet
Interior facet.
Definition Form.h:39
@ cell
Cell.
Definition Form.h:37
@ exterior_facet
Exterior facet.
Definition Form.h:38
int cell_num_entities(CellType type, int dim)
Definition cell_types.cpp:139
CellType
Cell type identifier.
Definition cell_types.h:22