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 <Eigen/Dense>
10 #include <functional>
11 #include <petscmat.h>
12 #include <petscsys.h>
13 #include <vector>
14 
15 namespace dolfinx
16 {
17 
18 namespace function
19 {
20 class Function;
21 }
22 
23 namespace graph
24 {
25 template <typename T>
27 }
28 
29 namespace mesh
30 {
31 class Mesh;
32 }
33 
34 namespace fem
35 {
36 class Form;
37 class DofMap;
38 
39 namespace impl
40 {
41 
47 
48 template <typename ScalarType>
49 void assemble_matrix(
50  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
51  const std::int32_t*, const ScalarType*)>&
52  mat_set_values_local,
53  const Form& a, const std::vector<bool>& bc0, const std::vector<bool>& bc1);
54 
56 template <typename ScalarType>
57 void assemble_cells(
58  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
59  const std::int32_t*, const ScalarType*)>&
60  mat_set_values_local,
61  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_cells,
62  const graph::AdjacencyList<std::int32_t>& dofmap0, int num_dofs_per_cell0,
63  const graph::AdjacencyList<std::int32_t>& dofmap1, int num_dofs_per_cell1,
64  const std::vector<bool>& bc0, const std::vector<bool>& bc1,
65  const std::function<void(ScalarType*, const ScalarType*, const ScalarType*,
66  const double*, const int*, const std::uint8_t*,
67  const std::uint32_t)>& kernel,
68  const Eigen::Array<ScalarType, Eigen::Dynamic, Eigen::Dynamic,
69  Eigen::RowMajor>& coeffs,
70  const Eigen::Array<ScalarType, Eigen::Dynamic, 1>& constant_values);
71 
73 template <typename ScalarType>
75  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
76  const std::int32_t*, const ScalarType*)>&
77  mat_set_values_local,
78  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_facets,
79  const DofMap& dofmap0, const DofMap& dofmap1, const std::vector<bool>& bc0,
80  const std::vector<bool>& bc1,
81  const std::function<void(ScalarType*, const ScalarType*, const ScalarType*,
82  const double*, const int*, const std::uint8_t*,
83  const std::uint32_t)>& fn,
84  const Eigen::Array<ScalarType, Eigen::Dynamic, Eigen::Dynamic,
85  Eigen::RowMajor>& coeffs,
86  const Eigen::Array<ScalarType, Eigen::Dynamic, 1> constant_values);
87 
89 template <typename ScalarType>
91  const std::function<int(std::int32_t, const std::int32_t*, std::int32_t,
92  const std::int32_t*, const ScalarType*)>&
93  mat_set_values_local,
94  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_facets,
95  const DofMap& dofmap0, const DofMap& dofmap1, const std::vector<bool>& bc0,
96  const std::vector<bool>& bc1,
97  const std::function<void(ScalarType*, const ScalarType*, const ScalarType*,
98  const double*, const int*, const std::uint8_t*,
99  const std::uint32_t)>& kernel,
100  const Eigen::Array<ScalarType, Eigen::Dynamic, Eigen::Dynamic,
101  Eigen::RowMajor>& coeffs,
102  const std::vector<int>& offsets,
103  const Eigen::Array<ScalarType, Eigen::Dynamic, 1>& constant_values);
104 
105 } // namespace impl
106 } // namespace fem
107 } // namespace dolfinx
dolfinx::fem::Form
Base class for variational forms.
Definition: Form.h:65
dolfinx::graph::AdjacencyList
This class provides a static adjacency list data structure. It is commonly used to store directed gra...
Definition: assemble_matrix_impl.h:26
dolfinx::fem::DofMap
Degree-of-freedom map.
Definition: DofMap.h:40
dolfinx::fem::impl::assemble_matrix
void assemble_matrix(const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const ScalarType *)> &mat_set_values_local, const Form &a, const std::vector< bool > &bc0, const std::vector< bool > &bc1)
The matrix A must already be initialised. The matrix may be a proxy, i.e. a view into a larger matrix...
Definition: assemble_matrix_impl.cpp:24
dolfinx::fem::impl::assemble_interior_facets
void assemble_interior_facets(const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const ScalarType *)> &mat_set_values_local, const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_facets, const DofMap &dofmap0, const DofMap &dofmap1, const std::vector< bool > &bc0, const std::vector< bool > &bc1, const std::function< void(ScalarType *, const ScalarType *, const ScalarType *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &kernel, const Eigen::Array< ScalarType, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const std::vector< int > &offsets, const Eigen::Array< ScalarType, Eigen::Dynamic, 1 > &constant_values)
Execute kernel over interior facets and accumulate result in Mat.
Definition: assemble_matrix_impl.cpp:274
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::impl::assemble_cells
void assemble_cells(const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const ScalarType *)> &mat_set_values_local, const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_cells, const graph::AdjacencyList< std::int32_t > &dofmap0, int num_dofs_per_cell0, const graph::AdjacencyList< std::int32_t > &dofmap1, int num_dofs_per_cell1, const std::vector< bool > &bc0, const std::vector< bool > &bc1, const std::function< void(ScalarType *, const ScalarType *, const ScalarType *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &kernel, const Eigen::Array< ScalarType, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const Eigen::Array< ScalarType, Eigen::Dynamic, 1 > &constant_values)
Execute kernel over cells and accumulate result in Mat.
Definition: assemble_matrix_impl.cpp:101
dolfinx::function::Function
This class represents a function in a finite element function space , given by.
Definition: Function.h:41
dolfinx::fem::impl::assemble_exterior_facets
void assemble_exterior_facets(const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const ScalarType *)> &mat_set_values_local, const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_facets, const DofMap &dofmap0, const DofMap &dofmap1, const std::vector< bool > &bc0, const std::vector< bool > &bc1, const std::function< void(ScalarType *, const ScalarType *, const ScalarType *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &fn, const Eigen::Array< ScalarType, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const Eigen::Array< ScalarType, Eigen::Dynamic, 1 > constant_values)
Execute kernel over exterior facets and accumulate result in Mat.
Definition: assemble_matrix_impl.cpp:178