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
105  = a.function_spaces().at(0)->dofmap();
106  std::shared_ptr<const fem::DofMap> dofmap1
107  = a.function_spaces().at(1)->dofmap();
108  assert(dofmap0);
109  assert(dofmap1);
110  const graph::AdjacencyList<std::int32_t>& dofs0 = dofmap0->list();
111  const graph::AdjacencyList<std::int32_t>& dofs1 = dofmap1->list();
112 
113  // Prepare constants
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 bool needs_permutation_data = a.needs_permutation_data();
121  if (needs_permutation_data)
122  mesh->topology_mutable().create_entity_permutations();
123  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info
124  = needs_permutation_data
125  ? mesh->topology().get_cell_permutation_info()
126  : Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>(num_cells);
127 
128  for (int i : a.integral_ids(IntegralType::cell))
129  {
130  const auto& fn = a.kernel(IntegralType::cell, i);
131  const std::vector<std::int32_t>& active_cells
132  = a.domains(IntegralType::cell, i);
133  impl::assemble_cells<T>(mat_set_values, mesh->geometry(), active_cells,
134  dofs0, dofs1, bc0, bc1, fn, coeffs, constants,
135  cell_info);
136  }
137 
138  if (a.num_integrals(IntegralType::exterior_facet) > 0
139  or a.num_integrals(IntegralType::interior_facet) > 0)
140  {
141  // FIXME: cleanup these calls? Some of the happen internally again.
142  mesh->topology_mutable().create_entities(tdim - 1);
143  mesh->topology_mutable().create_connectivity(tdim - 1, tdim);
144 
145  const int facets_per_cell
146  = mesh::cell_num_entities(mesh->topology().cell_type(), tdim - 1);
147  const Eigen::Array<std::uint8_t, Eigen::Dynamic, Eigen::Dynamic>& perms
148  = needs_permutation_data
149  ? mesh->topology().get_facet_permutations()
150  : Eigen::Array<std::uint8_t, Eigen::Dynamic, Eigen::Dynamic>(
151  facets_per_cell, num_cells);
152 
153  for (int i : a.integral_ids(IntegralType::exterior_facet))
154  {
155  const auto& fn = a.kernel(IntegralType::exterior_facet, i);
156  const std::vector<std::int32_t>& active_facets
157  = a.domains(IntegralType::exterior_facet, i);
158  impl::assemble_exterior_facets<T>(mat_set_values, *mesh, active_facets,
159  dofs0, dofs1, bc0, bc1, fn, coeffs,
160  constants, cell_info, perms);
161  }
162 
163  const std::vector<int> c_offsets = a.coefficient_offsets();
164  for (int i : a.integral_ids(IntegralType::interior_facet))
165  {
166  const auto& fn = a.kernel(IntegralType::interior_facet, i);
167  const std::vector<std::int32_t>& active_facets
168  = a.domains(IntegralType::interior_facet, i);
169  impl::assemble_interior_facets<T>(
170  mat_set_values, *mesh, active_facets, *dofmap0, *dofmap1, bc0, bc1,
171  fn, coeffs, c_offsets, constants, cell_info, perms);
172  }
173  }
174 }
175 //-----------------------------------------------------------------------------
176 template <typename T>
177 void assemble_cells(
178  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
179  const std::int32_t*, const T*)>& mat_set,
180  const mesh::Geometry& geometry,
181  const std::vector<std::int32_t>& active_cells,
182  const graph::AdjacencyList<std::int32_t>& dofmap0,
183  const graph::AdjacencyList<std::int32_t>& dofmap1,
184  const std::vector<bool>& bc0, const std::vector<bool>& bc1,
185  const std::function<void(T*, const T*, const T*, const double*, const int*,
186  const std::uint8_t*, const std::uint32_t)>& kernel,
187  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>&
188  coeffs,
189  const Eigen::Array<T, Eigen::Dynamic, 1>& constants,
190  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info)
191 {
192  const int gdim = geometry.dim();
193 
194  // Prepare cell geometry
195  const graph::AdjacencyList<std::int32_t>& x_dofmap = geometry.dofmap();
196 
197  // FIXME: Add proper interface for num coordinate dofs
198  const int num_dofs_g = x_dofmap.num_links(0);
199  const Eigen::Array<double, Eigen::Dynamic, 3, Eigen::RowMajor>& x_g
200  = geometry.x();
201 
202  // Data structures used in assembly
203  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
204  coordinate_dofs(num_dofs_g, gdim);
205  const int num_dofs0 = dofmap0.links(0).size();
206  const int num_dofs1 = dofmap1.links(0).size();
207  Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> Ae(
208  num_dofs0, num_dofs1);
209 
210  // Iterate over active cells
211  for (std::int32_t c : active_cells)
212  {
213  // Get cell coordinates/geometry
214  auto x_dofs = x_dofmap.links(c);
215  for (int i = 0; i < x_dofs.rows(); ++i)
216  for (int j = 0; j < gdim; ++j)
217  coordinate_dofs(i, j) = x_g(x_dofs[i], j);
218 
219  // Tabulate tensor
220  std::fill(Ae.data(), Ae.data() + num_dofs0 * num_dofs1, 0);
221  kernel(Ae.data(), coeffs.row(c).data(), constants.data(),
222  coordinate_dofs.data(), nullptr, nullptr, cell_info[c]);
223 
224  // Zero rows/columns for essential bcs
225  auto dofs0 = dofmap0.links(c);
226  auto dofs1 = dofmap1.links(c);
227  if (!bc0.empty())
228  {
229  for (Eigen::Index i = 0; i < Ae.rows(); ++i)
230  {
231  if (bc0[dofs0[i]])
232  Ae.row(i).setZero();
233  }
234  }
235  if (!bc1.empty())
236  {
237  for (Eigen::Index j = 0; j < Ae.cols(); ++j)
238  {
239  if (bc1[dofs1[j]])
240  Ae.col(j).setZero();
241  }
242  }
243 
244  mat_set(dofs0.size(), dofs0.data(), dofs1.size(), dofs1.data(), Ae.data());
245  }
246 }
247 //-----------------------------------------------------------------------------
248 template <typename T>
249 void assemble_exterior_facets(
250  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
251  const std::int32_t*, const T*)>& mat_set_values,
252  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_facets,
253  const graph::AdjacencyList<std::int32_t>& dofmap0,
254  const graph::AdjacencyList<std::int32_t>& dofmap1,
255  const std::vector<bool>& bc0, const std::vector<bool>& bc1,
256  const std::function<void(T*, const T*, const T*, const double*, const int*,
257  const std::uint8_t*, const std::uint32_t)>& kernel,
258  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>&
259  coeffs,
260  const Eigen::Array<T, Eigen::Dynamic, 1> constants,
261  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info,
262  const Eigen::Array<std::uint8_t, Eigen::Dynamic, Eigen::Dynamic>& perms)
263 {
264  const int gdim = mesh.geometry().dim();
265  const int tdim = mesh.topology().dim();
266 
267  // Prepare cell geometry
268  const graph::AdjacencyList<std::int32_t>& x_dofmap = mesh.geometry().dofmap();
269 
270  // FIXME: Add proper interface for num coordinate dofs
271  const int num_dofs_g = x_dofmap.num_links(0);
272  const Eigen::Array<double, Eigen::Dynamic, 3, Eigen::RowMajor>& x_g
273  = mesh.geometry().x();
274 
275  // Data structures used in assembly
276  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
277  coordinate_dofs(num_dofs_g, gdim);
278  const int num_dofs0 = dofmap0.links(0).size();
279  const int num_dofs1 = dofmap1.links(0).size();
280  Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> Ae(
281  num_dofs0, num_dofs1);
282 
283  // Iterate over all facets
284  auto f_to_c = mesh.topology().connectivity(tdim - 1, tdim);
285  assert(f_to_c);
286  auto c_to_f = mesh.topology().connectivity(tdim, tdim - 1);
287  assert(c_to_f);
288  for (std::int32_t f : active_facets)
289  {
290  auto cells = f_to_c->links(f);
291  assert(cells.rows() == 1);
292 
293  // Get local index of facet with respect to the cell
294  auto facets = c_to_f->links(cells[0]);
295  const auto* it = std::find(facets.data(), facets.data() + facets.rows(), f);
296  assert(it != (facets.data() + facets.rows()));
297  const int local_facet = std::distance(facets.data(), it);
298 
299  // Get cell vertex coordinates
300  auto x_dofs = x_dofmap.links(cells[0]);
301  for (int i = 0; i < num_dofs_g; ++i)
302  for (int j = 0; j < gdim; ++j)
303  coordinate_dofs(i, j) = x_g(x_dofs[i], j);
304 
305  // Tabulate tensor
306  std::fill(Ae.data(), Ae.data() + num_dofs0 * num_dofs1, 0);
307  kernel(Ae.data(), coeffs.row(cells[0]).data(), constants.data(),
308  coordinate_dofs.data(), &local_facet, &perms(local_facet, cells[0]),
309  cell_info[cells[0]]);
310 
311  // Zero rows/columns for essential bcs
312  auto dmap0 = dofmap0.links(cells[0]);
313  auto dmap1 = dofmap1.links(cells[0]);
314  if (!bc0.empty())
315  {
316  for (Eigen::Index i = 0; i < Ae.rows(); ++i)
317  {
318  if (bc0[dmap0[i]])
319  Ae.row(i).setZero();
320  }
321  }
322  if (!bc1.empty())
323  {
324  for (Eigen::Index j = 0; j < Ae.cols(); ++j)
325  {
326  if (bc1[dmap1[j]])
327  Ae.col(j).setZero();
328  }
329  }
330 
331  mat_set_values(dmap0.size(), dmap0.data(), dmap1.size(), dmap1.data(),
332  Ae.data());
333  }
334 }
335 //-----------------------------------------------------------------------------
336 template <typename T>
337 void assemble_interior_facets(
338  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
339  const std::int32_t*, const T*)>& mat_set_values,
340  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_facets,
341  const DofMap& dofmap0, const DofMap& dofmap1, const std::vector<bool>& bc0,
342  const std::vector<bool>& bc1,
343  const std::function<void(T*, const T*, const T*, const double*, const int*,
344  const std::uint8_t*, const std::uint32_t)>& fn,
345  const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>&
346  coeffs,
347  const std::vector<int>& offsets,
348  const Eigen::Array<T, Eigen::Dynamic, 1>& constants,
349  const Eigen::Array<std::uint32_t, Eigen::Dynamic, 1>& cell_info,
350  const Eigen::Array<std::uint8_t, Eigen::Dynamic, Eigen::Dynamic>& perms)
351 {
352  const int gdim = mesh.geometry().dim();
353  const int tdim = mesh.topology().dim();
354 
355  // Prepare cell geometry
356  const graph::AdjacencyList<std::int32_t>& x_dofmap = mesh.geometry().dofmap();
357 
358  // FIXME: Add proper interface for num coordinate dofs
359  const int num_dofs_g = x_dofmap.num_links(0);
360 
361  const Eigen::Array<double, Eigen::Dynamic, 3, Eigen::RowMajor>& x_g
362  = mesh.geometry().x();
363 
364  // Data structures used in assembly
365  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
366  coordinate_dofs(2 * num_dofs_g, gdim);
367  Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> Ae;
368  Eigen::Array<T, Eigen::Dynamic, 1> coeff_array(2 * offsets.back());
369  assert(offsets.back() == coeffs.cols());
370 
371  // Temporaries for joint dofmaps
372  Eigen::Array<std::int32_t, Eigen::Dynamic, 1> dmapjoint0, dmapjoint1;
373 
374  // Iterate over all facets
375  auto c = mesh.topology().connectivity(tdim - 1, tdim);
376  assert(c);
377  auto c_to_f = mesh.topology().connectivity(tdim, tdim - 1);
378  assert(c_to_f);
379  for (std::int32_t facet_index : active_facets)
380  {
381  // Create attached cells
382  auto cells = c->links(facet_index);
383  assert(cells.rows() == 2);
384 
385  // Get local index of facet with respect to the cell
386  auto facets0 = c_to_f->links(cells[0]);
387  const auto* it0 = std::find(facets0.data(), facets0.data() + facets0.rows(),
388  facet_index);
389  assert(it0 != (facets0.data() + facets0.rows()));
390  const int local_facet0 = std::distance(facets0.data(), it0);
391  auto facets1 = c_to_f->links(cells[1]);
392  const auto* it1 = std::find(facets1.data(), facets1.data() + facets1.rows(),
393  facet_index);
394  assert(it1 != (facets1.data() + facets1.rows()));
395  const int local_facet1 = std::distance(facets1.data(), it1);
396 
397  const std::array local_facet{local_facet0, local_facet1};
398 
399  // Get cell geometry
400  auto x_dofs0 = x_dofmap.links(cells[0]);
401  auto x_dofs1 = x_dofmap.links(cells[1]);
402  for (int i = 0; i < num_dofs_g; ++i)
403  {
404  for (int j = 0; j < gdim; ++j)
405  {
406  coordinate_dofs(i, j) = x_g(x_dofs0[i], j);
407  coordinate_dofs(i + num_dofs_g, j) = x_g(x_dofs1[i], j);
408  }
409  }
410 
411  // Get dof maps for cells and pack
412  auto dmap0_cell0 = dofmap0.cell_dofs(cells[0]);
413  auto dmap0_cell1 = dofmap0.cell_dofs(cells[1]);
414  dmapjoint0.resize(dmap0_cell0.size() + dmap0_cell1.size());
415  dmapjoint0.head(dmap0_cell0.size()) = dmap0_cell0;
416  dmapjoint0.tail(dmap0_cell1.size()) = dmap0_cell1;
417 
418  auto dmap1_cell0 = dofmap1.cell_dofs(cells[0]);
419  auto dmap1_cell1 = dofmap1.cell_dofs(cells[1]);
420  dmapjoint1.resize(dmap1_cell0.size() + dmap1_cell1.size());
421  dmapjoint1.head(dmap1_cell0.size()) = dmap1_cell0;
422  dmapjoint1.tail(dmap1_cell1.size()) = dmap1_cell1;
423 
424  // Layout for the restricted coefficients is flattened
425  // w[coefficient][restriction][dof]
426  auto coeff_cell0 = coeffs.row(cells[0]);
427  auto coeff_cell1 = coeffs.row(cells[1]);
428 
429  // Loop over coefficients
430  for (std::size_t i = 0; i < offsets.size() - 1; ++i)
431  {
432  // Loop over entries for coefficient i
433  const int num_entries = offsets[i + 1] - offsets[i];
434  coeff_array.segment(2 * offsets[i], num_entries)
435  = coeff_cell0.segment(offsets[i], num_entries);
436  coeff_array.segment(offsets[i + 1] + offsets[i], num_entries)
437  = coeff_cell1.segment(offsets[i], num_entries);
438  }
439 
440  // Tabulate tensor
441  Ae.setZero(dmapjoint0.size(), dmapjoint1.size());
442  const std::array perm{perms(local_facet[0], cells[0]),
443  perms(local_facet[1], cells[1])};
444  fn(Ae.data(), coeff_array.data(), constants.data(), coordinate_dofs.data(),
445  local_facet.data(), perm.data(), cell_info[cells[0]]);
446 
447  // Zero rows/columns for essential bcs
448  if (!bc0.empty())
449  {
450  for (Eigen::Index i = 0; i < dmapjoint0.size(); ++i)
451  {
452  if (bc0[dmapjoint0[i]])
453  Ae.row(i).setZero();
454  }
455  }
456  if (!bc1.empty())
457  {
458  for (Eigen::Index j = 0; j < dmapjoint1.size(); ++j)
459  {
460  if (bc1[dmapjoint1[j]])
461  Ae.col(j).setZero();
462  }
463  }
464 
465  mat_set_values(dmapjoint0.size(), dmapjoint0.data(), dmapjoint1.size(),
466  dmapjoint1.data(), Ae.data());
467  }
468 }
469 //-----------------------------------------------------------------------------
470 
471 } // namespace dolfinx::fem::impl
Degree-of-freedom map.
Definition: DofMap.h:66
Class for variational forms.
Definition: Form.h:66
std::vector< int > integral_ids(IntegralType type) const
Get the IDs for integrals (kernels) for given integral type. The IDs correspond to the domain IDs whi...
Definition: Form.h:208
std::shared_ptr< const mesh::Mesh > mesh() const
Extract common mesh for the form.
Definition: Form.h:152
std::vector< int > coefficient_offsets() const
Offset for each coefficient expansion array on a cell. Used to pack data for multiple coefficients in...
Definition: Form.h:251
const std::vector< std::int32_t > & domains(IntegralType type, int i) const
Get the list of mesh entity indices for the ith integral (kernel) for the given domain type,...
Definition: Form.h:225
const std::function< void(T *, const T *, const T *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> & kernel(IntegralType type, int i) const
Get the function for 'kernel' for integral i of given type.
Definition: Form.h:169
const std::vector< std::shared_ptr< const function::FunctionSpace > > & function_spaces() const
Return function spaces for all arguments.
Definition: Form.h:157
int num_integrals(IntegralType type) const
Number of integrals of given type.
Definition: Form.h:194
bool needs_permutation_data() const
Get bool indicating whether permutation data needs to be passed into these integrals.
Definition: Form.h:246
This class provides a static adjacency list data structure. It is commonly used to store directed gra...
Definition: AdjacencyList.h:28
Eigen::Array< T, Eigen::Dynamic, 1 >::SegmentReturnType links(int node)
Links (edges) for given node.
Definition: AdjacencyList.h:153
int num_links(int node) const
Number of connections for given node.
Definition: AdjacencyList.h:143
Geometry stores the geometry imposed on a mesh.
Definition: Geometry.h:39
const graph::AdjacencyList< std::int32_t > & dofmap() const
DOF map.
Definition: Geometry.cpp:22
Eigen::Array< double, Eigen::Dynamic, 3, Eigen::RowMajor > & x()
Geometry degrees-of-freedom.
Definition: Geometry.cpp:32
int dim() const
Return Euclidean dimension of coordinate system.
Definition: Geometry.cpp:20
A Mesh consists of a set of connected and numbered mesh topological entities, and geometry data.
Definition: Mesh.h:47
Geometry & geometry()
Get mesh geometry.
Definition: Mesh.cpp:147
Topology & topology()
Get mesh topology.
Definition: Mesh.cpp:141
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
int dim() const
Return topological dimension.
Definition: Topology.cpp:153
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
Eigen::Array< typename U::scalar_type, Eigen::Dynamic, 1 > pack_constants(const U &u)
Pack constants of u of generic type U ready for assembly.
Definition: utils.h:456
Eigen::Array< typename U::scalar_type, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > pack_coefficients(const U &u)
Pack coefficients of u of generic type U ready for assembly.
Definition: utils.h:407
int cell_num_entities(mesh::CellType type, int dim)
Number of entities of dimension dim.
Definition: cell_types.cpp:381