DOLFIN-X
DOLFIN-X C++ interface
assemble_vector_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 <dolfinx/common/types.h>
11 #include <functional>
12 #include <memory>
13 #include <petscsys.h>
14 #include <vector>
15 
16 namespace dolfinx
17 {
18 
19 namespace function
20 {
21 class Function;
22 }
23 
24 namespace graph
25 {
26 template <typename T>
27 class AdjacencyList;
28 }
29 
30 namespace mesh
31 {
32 class Mesh;
33 }
34 
35 namespace fem
36 {
37 class DirichletBC;
38 class Form;
39 class DofMap;
40 
42 namespace impl
43 {
44 
49 void assemble_vector(
50  Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b, const Form& L);
51 
53 void assemble_cells(
54  Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b,
55  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_cells,
56  const graph::AdjacencyList<std::int32_t>& dofmap, int num_dofs_per_cell,
57  const std::function<void(PetscScalar*, const PetscScalar*,
58  const PetscScalar*, const double*, const int*,
59  const std::uint8_t*, const std::uint32_t)>& kernel,
60  const Eigen::Array<PetscScalar, Eigen::Dynamic, Eigen::Dynamic,
61  Eigen::RowMajor>& coeffs,
62  const Eigen::Array<PetscScalar, Eigen::Dynamic, 1>& constant_values);
63 
66  Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b,
67  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_facets,
68  const fem::DofMap& dofmap,
69  const std::function<void(PetscScalar*, const PetscScalar*,
70  const PetscScalar*, const double*, const int*,
71  const std::uint8_t*, const std::uint32_t)>& fn,
72  const Eigen::Array<PetscScalar, Eigen::Dynamic, Eigen::Dynamic,
73  Eigen::RowMajor>& coeffs,
74  const Eigen::Array<PetscScalar, Eigen::Dynamic, 1>& constant_values);
75 
78  Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b,
79  const mesh::Mesh& mesh, const std::vector<std::int32_t>& active_facets,
80  const fem::DofMap& dofmap,
81  const std::function<void(PetscScalar*, const PetscScalar*,
82  const PetscScalar*, const double*, const int*,
83  const std::uint8_t*, const std::uint32_t)>& fn,
84  const Eigen::Array<PetscScalar, Eigen::Dynamic, Eigen::Dynamic,
85  Eigen::RowMajor>& coeffs,
86  const std::vector<int>& offsets,
87  const Eigen::Array<PetscScalar, Eigen::Dynamic, 1>& constant_values);
88 
106 void apply_lifting(
107  Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b,
108  const std::vector<std::shared_ptr<const Form>> a,
109  const std::vector<std::vector<std::shared_ptr<const DirichletBC>>>& bcs1,
110  const std::vector<
111  Eigen::Ref<const Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>>>& x0,
112  double scale);
113 
124 void lift_bc(
125  Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b, const Form& a,
126  const Eigen::Ref<const Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>>&
127  bc_values1,
128  const std::vector<bool>& bc_markers1, double scale);
129 
142 void lift_bc(
143  Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b, const Form& a,
144  const Eigen::Ref<const Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>>&
145  bc_values1,
146  const std::vector<bool>& bc_markers1,
147  const Eigen::Ref<const Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>>& x0,
148  double scale);
149 
150 } // namespace impl
151 } // namespace fem
152 } // namespace dolfinx
dolfinx::fem::impl::apply_lifting
void apply_lifting(Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> b, const std::vector< std::shared_ptr< const Form >> a, const std::vector< std::vector< std::shared_ptr< const DirichletBC >>> &bcs1, const std::vector< Eigen::Ref< const Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >>> &x0, double scale)
Modify b such that:
Definition: assemble_vector_impl.cpp:587
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::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::fem::impl::assemble_vector
void assemble_vector(Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> b, const Form &L)
Assemble linear form into an Eigen vector.
Definition: assemble_vector_impl.cpp:295
dolfinx::fem::impl::lift_bc
void lift_bc(Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> b, const Form &a, const Eigen::Ref< const Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> &bc_values1, const std::vector< bool > &bc_markers1, double scale)
Modify RHS vector to account for boundary condition.
Definition: assemble_vector_impl.cpp:637
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