DOLFIN-X
DOLFIN-X C++ interface
assemble_matrix_impl.h
1 // Copyright (C) 2018-2019 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 "DofMap.h"
10 #include "Form.h"
11 #include "utils.h"
12 #include <Eigen/Dense>
13 #include <dolfinx/function/FunctionSpace.h>
14 #include <dolfinx/graph/AdjacencyList.h>
15 #include <dolfinx/la/utils.h>
16 #include <dolfinx/mesh/Geometry.h>
17 #include <dolfinx/mesh/Mesh.h>
18 #include <dolfinx/mesh/Topology.h>
19 #include <functional>
20 #include <vector>
21 
22 namespace dolfinx::fem::impl
23 {
24 
30 
31 template <typename T>
32 void assemble_matrix(
33  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
34  const std::int32_t*, const T*)>& mat_set_values,
35  const Form<T>& a, const std::vector<bool>& bc0,
36  const std::vector<bool>& bc1);
37 
39 template <typename T>
40 void assemble_cells(
41  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
42  const std::int32_t*, const T*)>& mat_set_values,
43  const mesh::Geometry& geometry,
44  const std::vector<std::int32_t>& active_cells,
47  const std::vector<bool>& bc0, const std::vector<bool>& bc1,
48  const std::function<void(T*, const T*, const T*, const double*, const int*,
49  const std::uint8_t*, const std::uint32_t)>& kernel,
50  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>&
51  coeffs,
52  const Eigen::Array<T, Eigen::Dynamic, 1>& constants,
53  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info);
54 
56 template <typename T>
57 void assemble_exterior_facets(
58  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
59  const std::int32_t*, const T*)>& mat_set_values,
60  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_facets,
63  const std::vector<bool>& bc0, const std::vector<bool>& bc1,
64  const std::function<void(T*, const T*, const T*, const double*, const int*,
65  const std::uint8_t*, const std::uint32_t)>& fn,
66  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>&
67  coeffs,
68  const Eigen::Array<T, Eigen::Dynamic, 1> constants,
69  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info,
70  const Eigen::Array<std::uint8_t, Eigen::Dynamic, Eigen::Dynamic>& perms);
71 
73 template <typename T>
74 void assemble_interior_facets(
75  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
76  const std::int32_t*, const T*)>& mat_set_values,
77  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_facets,
78  const DofMap& dofmap0, const DofMap& dofmap1, const std::vector<bool>& bc0,
79  const std::vector<bool>& bc1,
80  const std::function<void(T*, const T*, const T*, const double*, const int*,
81  const std::uint8_t*, const std::uint32_t)>& kernel,
82  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>&
83  coeffs,
84  const std::vector<int>& offsets,
85  const Eigen::Array<T, Eigen::Dynamic, 1>& constants,
86  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info,
87  const Eigen::Array<std::uint8_t, Eigen::Dynamic, Eigen::Dynamic>& perms);
88 
89 //-----------------------------------------------------------------------------
90 template <typename T>
91 void assemble_matrix(
92  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
93  const std::int32_t*, const T*)>& mat_set_values,
94  const Form<T>& a, const std::vector<bool>& bc0,
95  const std::vector<bool>& bc1)
96 {
97  std::shared_ptr<const mesh::Mesh> mesh = a.mesh();
98  assert(mesh);
99  const int tdim = mesh->topology().dim();
100  const std::int32_t num_cells
101  = mesh->topology().connectivity(tdim, 0)->num_nodes();
102 
103  // Get dofmap data
104  std::shared_ptr<const fem::DofMap> dofmap0 = a.function_space(0)->dofmap();
105  std::shared_ptr<const fem::DofMap> dofmap1 = a.function_space(1)->dofmap();
106  assert(dofmap0);
107  assert(dofmap1);
108  const graph::AdjacencyList<std::int32_t>& dofs0 = dofmap0->list();
109  const graph::AdjacencyList<std::int32_t>& dofs1 = dofmap1->list();
110 
111  // Prepare constants
112  if (!a.all_constants_set())
113  throw std::runtime_error("Unset constant in Form");
114  const Eigen::Array<T, Eigen::Dynamic, 1> constants = pack_constants(a);
115 
116  // Prepare coefficients
117  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> coeffs
118  = pack_coefficients(a);
119 
120  const FormIntegrals<T>& integrals = a.integrals();
121  const bool needs_permutation_data = integrals.needs_permutation_data();
122  if (needs_permutation_data)
123  mesh->topology_mutable().create_entity_permutations();
124  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info
125  = needs_permutation_data
126  ? mesh->topology().get_cell_permutation_info()
127  : Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>(num_cells);
128 
129  for (int i = 0; i < integrals.num_integrals(IntegralType::cell); ++i)
130  {
131  const auto& fn = integrals.get_tabulate_tensor(IntegralType::cell, i);
132  const std::vector<std::int32_t>& active_cells
133  = integrals.integral_domains(IntegralType::cell, i);
134  impl::assemble_cells<T>(mat_set_values, mesh->geometry(), active_cells,
135  dofs0, dofs1, bc0, bc1, fn, coeffs, constants,
136  cell_info);
137  }
138 
139  if (integrals.num_integrals(IntegralType::exterior_facet) > 0
140  or integrals.num_integrals(IntegralType::interior_facet) > 0)
141  {
142  // FIXME: cleanup these calls? Some of the happen internally again.
143  mesh->topology_mutable().create_entities(tdim - 1);
144  mesh->topology_mutable().create_connectivity(tdim - 1, tdim);
145 
146  const int facets_per_cell
147  = mesh::cell_num_entities(mesh->topology().cell_type(), tdim - 1);
148  const Eigen::Array<std::uint8_t, Eigen::Dynamic, Eigen::Dynamic>& perms
149  = needs_permutation_data
150  ? mesh->topology().get_facet_permutations()
151  : Eigen::Array<std::uint8_t, Eigen::Dynamic, Eigen::Dynamic>(
152  facets_per_cell, num_cells);
153 
154  for (int i = 0; i < integrals.num_integrals(IntegralType::exterior_facet);
155  ++i)
156  {
157  const auto& fn
158  = integrals.get_tabulate_tensor(IntegralType::exterior_facet, i);
159  const std::vector<std::int32_t>& active_facets
160  = integrals.integral_domains(IntegralType::exterior_facet, i);
161  impl::assemble_exterior_facets<T>(mat_set_values, *mesh, active_facets,
162  dofs0, dofs1, bc0, bc1, fn, coeffs,
163  constants, cell_info, perms);
164  }
165 
166  const std::vector<int> c_offsets = a.coefficients().offsets();
167  for (int i = 0; i < integrals.num_integrals(IntegralType::interior_facet);
168  ++i)
169  {
170  const auto& fn
171  = integrals.get_tabulate_tensor(IntegralType::interior_facet, i);
172  const std::vector<std::int32_t>& active_facets
173  = integrals.integral_domains(IntegralType::interior_facet, i);
174  impl::assemble_interior_facets<T>(
175  mat_set_values, *mesh, active_facets, *dofmap0, *dofmap1, bc0, bc1,
176  fn, coeffs, c_offsets, constants, cell_info, perms);
177  }
178  }
179 }
180 //-----------------------------------------------------------------------------
181 template <typename T>
182 void assemble_cells(
183  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
184  const std::int32_t*, const T*)>& mat_set,
185  const mesh::Geometry& geometry,
186  const std::vector<std::int32_t>& active_cells,
187  const graph::AdjacencyList<std::int32_t>& dofmap0,
188  const graph::AdjacencyList<std::int32_t>& dofmap1,
189  const std::vector<bool>& bc0, const std::vector<bool>& bc1,
190  const std::function<void(T*, const T*, const T*, const double*, const int*,
191  const std::uint8_t*, const std::uint32_t)>& kernel,
192  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>&
193  coeffs,
194  const Eigen::Array<T, Eigen::Dynamic, 1>& constants,
195  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info)
196 {
197  const int gdim = geometry.dim();
198 
199  // Prepare cell geometry
200  const graph::AdjacencyList<std::int32_t>& x_dofmap = geometry.dofmap();
201 
202  // FIXME: Add proper interface for num coordinate dofs
203  const int num_dofs_g = x_dofmap.num_links(0);
204  const Eigen::Array<double, Eigen::Dynamic, 3, Eigen::RowMajor>& x_g
205  = geometry.x();
206 
207  // Data structures used in assembly
208  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
209  coordinate_dofs(num_dofs_g, gdim);
210  const int num_dofs0 = dofmap0.links(0).size();
211  const int num_dofs1 = dofmap1.links(0).size();
212  Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> Ae(
213  num_dofs0, num_dofs1);
214 
215  // Iterate over active cells
216  for (std::int32_t c : active_cells)
217  {
218  // Get cell coordinates/geometry
219  auto x_dofs = x_dofmap.links(c);
220  for (int i = 0; i < x_dofs.rows(); ++i)
221  for (int j = 0; j < gdim; ++j)
222  coordinate_dofs(i, j) = x_g(x_dofs[i], j);
223 
224  // Tabulate tensor
225  std::fill(Ae.data(), Ae.data() + num_dofs0 * num_dofs1, 0);
226  kernel(Ae.data(), coeffs.row(c).data(), constants.data(),
227  coordinate_dofs.data(), nullptr, nullptr, cell_info[c]);
228 
229  // Zero rows/columns for essential bcs
230  auto dofs0 = dofmap0.links(c);
231  auto dofs1 = dofmap1.links(c);
232  if (!bc0.empty())
233  {
234  for (Eigen::Index i = 0; i < Ae.rows(); ++i)
235  {
236  if (bc0[dofs0[i]])
237  Ae.row(i).setZero();
238  }
239  }
240  if (!bc1.empty())
241  {
242  for (Eigen::Index j = 0; j < Ae.cols(); ++j)
243  {
244  if (bc1[dofs1[j]])
245  Ae.col(j).setZero();
246  }
247  }
248 
249  mat_set(dofs0.size(), dofs0.data(), dofs1.size(), dofs1.data(), Ae.data());
250  }
251 }
252 //-----------------------------------------------------------------------------
253 template <typename T>
254 void assemble_exterior_facets(
255  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
256  const std::int32_t*, const T*)>& mat_set_values,
257  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_facets,
258  const graph::AdjacencyList<std::int32_t>& dofmap0,
259  const graph::AdjacencyList<std::int32_t>& dofmap1,
260  const std::vector<bool>& bc0, const std::vector<bool>& bc1,
261  const std::function<void(T*, const T*, const T*, const double*, const int*,
262  const std::uint8_t*, const std::uint32_t)>& kernel,
263  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>&
264  coeffs,
265  const Eigen::Array<T, Eigen::Dynamic, 1> constants,
266  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info,
267  const Eigen::Array<std::uint8_t, Eigen::Dynamic, Eigen::Dynamic>& perms)
268 {
269  const int gdim = mesh.geometry().dim();
270  const int tdim = mesh.topology().dim();
271 
272  // Prepare cell geometry
273  const graph::AdjacencyList<std::int32_t>& x_dofmap = mesh.geometry().dofmap();
274 
275  // FIXME: Add proper interface for num coordinate dofs
276  const int num_dofs_g = x_dofmap.num_links(0);
277  const Eigen::Array<double, Eigen::Dynamic, 3, Eigen::RowMajor>& x_g
278  = mesh.geometry().x();
279 
280  // Data structures used in assembly
281  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
282  coordinate_dofs(num_dofs_g, gdim);
283  const int num_dofs0 = dofmap0.links(0).size();
284  const int num_dofs1 = dofmap1.links(0).size();
285  Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> Ae(
286  num_dofs0, num_dofs1);
287 
288  // Iterate over all facets
289  auto f_to_c = mesh.topology().connectivity(tdim - 1, tdim);
290  assert(f_to_c);
291  auto c_to_f = mesh.topology().connectivity(tdim, tdim - 1);
292  assert(c_to_f);
293  for (std::int32_t f : active_facets)
294  {
295  auto cells = f_to_c->links(f);
296  assert(cells.rows() == 1);
297 
298  // Get local index of facet with respect to the cell
299  auto facets = c_to_f->links(cells[0]);
300  const auto* it = std::find(facets.data(), facets.data() + facets.rows(), f);
301  assert(it != (facets.data() + facets.rows()));
302  const int local_facet = std::distance(facets.data(), it);
303 
304  // Get cell vertex coordinates
305  auto x_dofs = x_dofmap.links(cells[0]);
306  for (int i = 0; i < num_dofs_g; ++i)
307  for (int j = 0; j < gdim; ++j)
308  coordinate_dofs(i, j) = x_g(x_dofs[i], j);
309 
310  // Tabulate tensor
311  std::fill(Ae.data(), Ae.data() + num_dofs0 * num_dofs1, 0);
312  kernel(Ae.data(), coeffs.row(cells[0]).data(), constants.data(),
313  coordinate_dofs.data(), &local_facet, &perms(local_facet, cells[0]),
314  cell_info[cells[0]]);
315 
316  // Zero rows/columns for essential bcs
317  auto dmap0 = dofmap0.links(cells[0]);
318  auto dmap1 = dofmap1.links(cells[0]);
319  if (!bc0.empty())
320  {
321  for (Eigen::Index i = 0; i < Ae.rows(); ++i)
322  {
323  if (bc0[dmap0[i]])
324  Ae.row(i).setZero();
325  }
326  }
327  if (!bc1.empty())
328  {
329  for (Eigen::Index j = 0; j < Ae.cols(); ++j)
330  {
331  if (bc1[dmap1[j]])
332  Ae.col(j).setZero();
333  }
334  }
335 
336  mat_set_values(dmap0.size(), dmap0.data(), dmap1.size(), dmap1.data(),
337  Ae.data());
338  }
339 }
340 //-----------------------------------------------------------------------------
341 template <typename T>
342 void assemble_interior_facets(
343  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
344  const std::int32_t*, const T*)>& mat_set_values,
345  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_facets,
346  const DofMap& dofmap0, const DofMap& dofmap1, const std::vector<bool>& bc0,
347  const std::vector<bool>& bc1,
348  const std::function<void(T*, const T*, const T*, const double*, const int*,
349  const std::uint8_t*, const std::uint32_t)>& fn,
350  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>&
351  coeffs,
352  const std::vector<int>& offsets,
353  const Eigen::Array<T, Eigen::Dynamic, 1>& constants,
354  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info,
355  const Eigen::Array<std::uint8_t, Eigen::Dynamic, Eigen::Dynamic>& perms)
356 {
357  const int gdim = mesh.geometry().dim();
358  const int tdim = mesh.topology().dim();
359 
360  // Prepare cell geometry
361  const graph::AdjacencyList<std::int32_t>& x_dofmap = mesh.geometry().dofmap();
362 
363  // FIXME: Add proper interface for num coordinate dofs
364  const int num_dofs_g = x_dofmap.num_links(0);
365 
366  const Eigen::Array<double, Eigen::Dynamic, 3, Eigen::RowMajor>& x_g
367  = mesh.geometry().x();
368 
369  // Data structures used in assembly
370  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
371  coordinate_dofs(2 * num_dofs_g, gdim);
372  Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> Ae;
373  Eigen::Array<T, Eigen::Dynamic, 1> coeff_array(2 * offsets.back());
374  assert(offsets.back() == coeffs.cols());
375 
376  // Temporaries for joint dofmaps
377  Eigen::Array<std::int32_t, Eigen::Dynamic, 1> dmapjoint0, dmapjoint1;
378 
379  // Iterate over all facets
380  auto c = mesh.topology().connectivity(tdim - 1, tdim);
381  assert(c);
382  auto c_to_f = mesh.topology().connectivity(tdim, tdim - 1);
383  assert(c_to_f);
384  for (std::int32_t facet_index : active_facets)
385  {
386  // Create attached cells
387  auto cells = c->links(facet_index);
388  assert(cells.rows() == 2);
389 
390  // Get local index of facet with respect to the cell
391  auto facets0 = c_to_f->links(cells[0]);
392  const auto* it0 = std::find(facets0.data(), facets0.data() + facets0.rows(),
393  facet_index);
394  assert(it0 != (facets0.data() + facets0.rows()));
395  const int local_facet0 = std::distance(facets0.data(), it0);
396  auto facets1 = c_to_f->links(cells[1]);
397  const auto* it1 = std::find(facets1.data(), facets1.data() + facets1.rows(),
398  facet_index);
399  assert(it1 != (facets1.data() + facets1.rows()));
400  const int local_facet1 = std::distance(facets1.data(), it1);
401 
402  const std::array local_facet{local_facet0, local_facet1};
403 
404  // Get cell geometry
405  auto x_dofs0 = x_dofmap.links(cells[0]);
406  auto x_dofs1 = x_dofmap.links(cells[1]);
407  for (int i = 0; i < num_dofs_g; ++i)
408  {
409  for (int j = 0; j < gdim; ++j)
410  {
411  coordinate_dofs(i, j) = x_g(x_dofs0[i], j);
412  coordinate_dofs(i + num_dofs_g, j) = x_g(x_dofs1[i], j);
413  }
414  }
415 
416  // Get dof maps for cells and pack
417  auto dmap0_cell0 = dofmap0.cell_dofs(cells[0]);
418  auto dmap0_cell1 = dofmap0.cell_dofs(cells[1]);
419  dmapjoint0.resize(dmap0_cell0.size() + dmap0_cell1.size());
420  dmapjoint0.head(dmap0_cell0.size()) = dmap0_cell0;
421  dmapjoint0.tail(dmap0_cell1.size()) = dmap0_cell1;
422 
423  auto dmap1_cell0 = dofmap1.cell_dofs(cells[0]);
424  auto dmap1_cell1 = dofmap1.cell_dofs(cells[1]);
425  dmapjoint1.resize(dmap1_cell0.size() + dmap1_cell1.size());
426  dmapjoint1.head(dmap1_cell0.size()) = dmap1_cell0;
427  dmapjoint1.tail(dmap1_cell1.size()) = dmap1_cell1;
428 
429  // Layout for the restricted coefficients is flattened
430  // w[coefficient][restriction][dof]
431  auto coeff_cell0 = coeffs.row(cells[0]);
432  auto coeff_cell1 = coeffs.row(cells[1]);
433 
434  // Loop over coefficients
435  for (std::size_t i = 0; i < offsets.size() - 1; ++i)
436  {
437  // Loop over entries for coefficient i
438  const int num_entries = offsets[i + 1] - offsets[i];
439  coeff_array.segment(2 * offsets[i], num_entries)
440  = coeff_cell0.segment(offsets[i], num_entries);
441  coeff_array.segment(offsets[i + 1] + offsets[i], num_entries)
442  = coeff_cell1.segment(offsets[i], num_entries);
443  }
444 
445  // Tabulate tensor
446  Ae.setZero(dmapjoint0.size(), dmapjoint1.size());
447  const std::array perm{perms(local_facet[0], cells[0]),
448  perms(local_facet[1], cells[1])};
449  fn(Ae.data(), coeff_array.data(), constants.data(), coordinate_dofs.data(),
450  local_facet.data(), perm.data(), cell_info[cells[0]]);
451 
452  // Zero rows/columns for essential bcs
453  if (!bc0.empty())
454  {
455  for (Eigen::Index i = 0; i < dmapjoint0.size(); ++i)
456  {
457  if (bc0[dmapjoint0[i]])
458  Ae.row(i).setZero();
459  }
460  }
461  if (!bc1.empty())
462  {
463  for (Eigen::Index j = 0; j < dmapjoint1.size(); ++j)
464  {
465  if (bc1[dmapjoint1[j]])
466  Ae.col(j).setZero();
467  }
468  }
469 
470  mat_set_values(dmapjoint0.size(), dmapjoint0.data(), dmapjoint1.size(),
471  dmapjoint1.data(), Ae.data());
472  }
473 }
474 //-----------------------------------------------------------------------------
475 
476 } // namespace dolfinx::fem::impl
dolfinx::fem::FormIntegrals::num_integrals
int num_integrals(IntegralType type) const
Number of integrals of given type.
Definition: FormIntegrals.h:131
dolfinx::graph::AdjacencyList::num_links
int num_links(int node) const
Number of connections for given node.
Definition: AdjacencyList.h:143
dolfinx::fem::Form::all_constants_set
bool all_constants_set() const
Check if all constants associated with the form have been set.
Definition: Form.h:164
dolfinx::mesh::cell_num_entities
int cell_num_entities(mesh::CellType type, int dim)
Number of entities of dimension dim.
Definition: cell_types.cpp:381
dolfinx::fem::FormIntegrals::integral_domains
const std::vector< std::int32_t > & integral_domains(IntegralType type, int i) const
Get the list of active entities for the ith integral of type t. Note, these are not retrieved by ID,...
Definition: FormIntegrals.h:157
dolfinx::mesh::Topology::dim
int dim() const
Return topological dimension.
Definition: Topology.cpp:153
dolfinx::graph::AdjacencyList::links
Eigen::Array< T, Eigen::Dynamic, 1 >::SegmentReturnType links(int node)
Links (edges) for given node.
Definition: AdjacencyList.h:153
dolfinx::mesh::Geometry::dofmap
const graph::AdjacencyList< std::int32_t > & dofmap() const
DOF map.
Definition: Geometry.cpp:22
dolfinx::fem::Form
Class for variational forms.
Definition: assembler.h:28
dolfinx::graph::AdjacencyList
This class provides a static adjacency list data structure. It is commonly used to store directed gra...
Definition: AdjacencyList.h:28
dolfinx::mesh::Mesh::geometry
Geometry & geometry()
Get mesh geometry.
Definition: Mesh.cpp:147
dolfinx::fem::DofMap
Degree-of-freedom map.
Definition: DofMap.h:66
dolfinx::fem::Form::coefficients
FormCoefficients< T > & coefficients()
Access coefficients.
Definition: Form.h:230
dolfinx::mesh::Geometry::dim
int dim() const
Return Euclidean dimension of coordinate system.
Definition: Geometry.cpp:20
dolfinx::fem::FormIntegrals::get_tabulate_tensor
const std::function< void(T *, const T *, const T *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> & get_tabulate_tensor(IntegralType type, int i) const
Get the function for 'tabulate_tensor' for integral i of given type.
Definition: FormIntegrals.h:75
dolfinx::fem::Form::integrals
const FormIntegrals< T > & integrals() const
Access form integrals.
Definition: Form.h:236
dolfinx::mesh::Mesh
A Mesh consists of a set of connected and numbered mesh topological entities, and geometry data.
Definition: Mesh.h:47
dolfinx::fem::pack_constants
Eigen::Array< T, Eigen::Dynamic, 1 > pack_constants(const fem::Form< T > &form)
Pack form constants ready for assembly.
Definition: utils.h:359
dolfinx::mesh::Topology::connectivity
std::shared_ptr< const graph::AdjacencyList< std::int32_t > > connectivity(int d0, int d1) const
Return connectivity from entities of dimension d0 to entities of dimension d1.
Definition: Topology.cpp:255
dolfinx::mesh::Mesh::topology
Topology & topology()
Get mesh topology.
Definition: Mesh.cpp:141
dolfinx::fem::Form::function_space
std::shared_ptr< const function::FunctionSpace > function_space(int i) const
Return function space for given argument.
Definition: Form.h:204
dolfinx::fem::Form::mesh
std::shared_ptr< const mesh::Mesh > mesh() const
Extract common mesh from form.
Definition: Form.h:199
dolfinx::fem::DofMap::cell_dofs
Eigen::Array< std::int32_t, Eigen::Dynamic, 1 >::ConstSegmentReturnType cell_dofs(int cell) const
Local-to-global mapping of dofs on a cell.
Definition: DofMap.h:95
dolfinx::fem::pack_coefficients
Eigen::Array< T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > pack_coefficients(const fem::Form< T > &form)
Pack form coefficients ready for assembly.
Definition: utils.h:314
dolfinx::fem::FormIntegrals::needs_permutation_data
bool needs_permutation_data() const
Get bool indicating whether permutation data needs to be passed into these integrals.
Definition: FormIntegrals.h:357
dolfinx::fem::FormIntegrals
Integrals of a Form, including those defined over cells, interior and exterior facets,...
Definition: FormIntegrals.h:39
dolfinx::mesh::Geometry::x
Eigen::Array< double, Eigen::Dynamic, 3, Eigen::RowMajor > & x()
Geometry degrees-of-freedom.
Definition: Geometry.cpp:32
dolfinx::mesh::Geometry
Geometry stores the geometry imposed on a mesh.
Definition: Geometry.h:39
dolfinx::fem::assemble_matrix
void assemble_matrix(const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const T *)> &mat_add, const Form< T > &a, const std::vector< std::shared_ptr< const DirichletBC< T >>> &bcs)
Assemble bilinear form into a matrix.
Definition: assembler.h:95