DOLFIN-X
DOLFIN-X C++ interface
assembler.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 <Eigen/Sparse>
11 #include <dolfinx/common/types.h>
12 #include <memory>
13 #include <petscmat.h>
14 #include <petscvec.h>
15 #include <vector>
16 
17 namespace dolfinx
18 {
19 namespace function
20 {
21 class FunctionSpace;
22 } // namespace function
23 
24 namespace fem
25 {
26 class DirichletBC;
27 class Form;
28 
29 // -- Scalar ----------------------------------------------------------------
30 
36 PetscScalar assemble_scalar(const Form& M);
37 
38 // -- Vectors ----------------------------------------------------------------
39 
49 void assemble_vector(Vec b, const Form& L);
50 
55 void assemble_vector(
56  Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b, const Form& L);
57 
58 // FIXME: clarify how x0 is used
59 // FIXME: if bcs entries are set
60 
61 // FIXME: need to pass an array of Vec for x0?
62 // FIXME: clarify zeroing of vector
63 
76 void apply_lifting(
77  Vec b, const std::vector<std::shared_ptr<const Form>>& a,
78  const std::vector<std::vector<std::shared_ptr<const DirichletBC>>>& bcs1,
79  const std::vector<Vec>& x0, double scale);
80 
93 void apply_lifting(
94  Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b,
95  const std::vector<std::shared_ptr<const Form>>& a,
96  const std::vector<std::vector<std::shared_ptr<const DirichletBC>>>& bcs1,
97  const std::vector<
98  Eigen::Ref<const Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>>>& x0,
99  double scale);
100 
101 // -- Matrices ---------------------------------------------------------------
102 
103 // Experimental
109 Eigen::SparseMatrix<PetscScalar, Eigen::RowMajor> assemble_matrix_eigen(
110  const Form& a, const std::vector<std::shared_ptr<const DirichletBC>>& bcs);
111 
123 void assemble_matrix(
124  Mat A, const Form& a,
125  const std::vector<std::shared_ptr<const DirichletBC>>& bcs);
126 
138 void assemble_matrix(Mat A, const Form& a, const std::vector<bool>& bc0,
139  const std::vector<bool>& bc1);
140 
157 void add_diagonal(Mat A, const function::FunctionSpace& V,
158  const std::vector<std::shared_ptr<const DirichletBC>>& bcs,
159  PetscScalar diagonal = 1.0);
160 
161 // Developer note: This function calls MatSetValuesLocal and not
162 // MatZeroRowsLocal
163 // (https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatZeroRowsLocal.html)
164 
175 void add_diagonal(
176  Mat A,
177  const Eigen::Ref<const Eigen::Array<std::int32_t, Eigen::Dynamic, 1>>& rows,
178  PetscScalar diagonal = 1.0);
179 
180 // -- Setting bcs ------------------------------------------------------------
181 
182 // FIXME: Move these function elsewhere?
183 
184 // FIXME: clarify x0
185 // FIXME: clarify what happens with ghosts
186 
190 void set_bc(Vec b, const std::vector<std::shared_ptr<const DirichletBC>>& bcs,
191  const Vec x0, double scale = 1.0);
192 
196 void set_bc(
197  Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b,
198  const std::vector<std::shared_ptr<const DirichletBC>>& bcs,
199  const Eigen::Ref<const Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>>& x0,
200  double scale = 1.0);
201 
205 void set_bc(Eigen::Ref<Eigen::Matrix<PetscScalar, Eigen::Dynamic, 1>> b,
206  const std::vector<std::shared_ptr<const DirichletBC>>& bcs,
207  double scale = 1.0);
208 
209 // FIXME: Handle null block
210 // FIXME: Pass function spaces rather than forms
218 std::vector<std::vector<std::shared_ptr<const fem::DirichletBC>>>
219 bcs_rows(const std::vector<const Form*>& L,
220  const std::vector<std::shared_ptr<const fem::DirichletBC>>& bcs);
221 
222 // FIXME: Handle null block
223 // FIXME: Pass function spaces rather than forms
231 std::vector<std::vector<std::vector<std::shared_ptr<const fem::DirichletBC>>>>
232 bcs_cols(const std::vector<std::vector<std::shared_ptr<const Form>>>& a,
233  const std::vector<std::shared_ptr<const DirichletBC>>& bcs);
234 
235 } // namespace fem
236 } // namespace dolfinx
dolfinx::fem::bcs_rows
std::vector< std::vector< std::shared_ptr< const fem::DirichletBC > > > bcs_rows(const std::vector< const Form * > &L, const std::vector< std::shared_ptr< const fem::DirichletBC >> &bcs)
Arrange boundary conditions by block.
Definition: assembler.cpp:324
dolfinx::fem::assemble_vector
void assemble_vector(Vec b, const Form &L)
Assemble linear form into an already allocated PETSc vector. Ghost contributions are not accumulated ...
Definition: assembler.cpp:83
dolfinx::fem::set_bc
void set_bc(Vec b, const std::vector< std::shared_ptr< const DirichletBC >> &bcs, const Vec x0, double scale=1.0)
Set bc values in owned (local) part of the PETScVector, multiplied by 'scale'. The vectors b and x0 m...
Definition: assembler.cpp:281
dolfinx::fem::assemble_scalar
PetscScalar assemble_scalar(const Form &M)
Assemble functional into scalar. Caller is responsible for accumulation across processes.
Definition: assembler.cpp:78
dolfinx::fem::add_diagonal
void add_diagonal(Mat A, const function::FunctionSpace &V, const std::vector< std::shared_ptr< const DirichletBC >> &bcs, PetscScalar diagonal=1.0)
Adds a value to the diagonal of the matrix for rows with a Dirichlet boundary conditions applied....
Definition: assembler.cpp:240
dolfinx::fem::assemble_matrix_eigen
Eigen::SparseMatrix< PetscScalar, Eigen::RowMajor > assemble_matrix_eigen(const Form &a, const std::vector< std::shared_ptr< const DirichletBC >> &bcs)
Assemble bilinear form into an Eigen Sparse matrix.
Definition: assembler.cpp:130
dolfinx::fem::assemble_matrix
void assemble_matrix(Mat A, const Form &a, const std::vector< std::shared_ptr< const DirichletBC >> &bcs)
Assemble bilinear form into a matrix. Matrix must already be initialised. Does not zero or finalise t...
Definition: assembler.cpp:190
dolfinx::fem::bcs_cols
std::vector< std::vector< std::vector< std::shared_ptr< const fem::DirichletBC > > > > bcs_cols(const std::vector< std::vector< std::shared_ptr< const Form >>> &a, const std::vector< std::shared_ptr< const DirichletBC >> &bcs)
Arrange boundary conditions by block.
Definition: assembler.cpp:339
dolfinx::fem::apply_lifting
void apply_lifting(Vec b, const std::vector< std::shared_ptr< const Form >> &a, const std::vector< std::vector< std::shared_ptr< const DirichletBC >>> &bcs1, const std::vector< Vec > &x0, double scale)
Modify b such that:
Definition: assembler.cpp:95