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