14#include <boost/lexical_cast.hpp>
19#include <petscoptions.h>
38void error(
int error_code, std::string filename, std::string petsc_function);
48create_vectors(MPI_Comm comm,
49 const std::vector<std::span<const PetscScalar>>& x);
56Vec create_vector(
const common::IndexMap& map,
int bs);
66Vec create_vector(MPI_Comm comm, std::array<std::int64_t, 2> range,
67 std::span<const std::int64_t> ghosts,
int bs);
78Vec create_vector_wrap(
const common::IndexMap& map,
int bs,
79 std::span<const PetscScalar> x);
85Vec create_vector_wrap(
const la::Vector<V>& x)
87 assert(x.index_map());
88 return create_vector_wrap(*x.index_map(), x.bs(), x.array());
103std::vector<IS> create_index_sets(
105 std::pair<std::reference_wrapper<const common::IndexMap>,
int>>& maps);
108std::vector<std::vector<PetscScalar>> get_local_vectors(
111 std::pair<std::reference_wrapper<const common::IndexMap>,
int>>& maps);
114void scatter_local_vectors(
115 Vec x,
const std::vector<std::span<const PetscScalar>>& x_b,
117 std::pair<std::reference_wrapper<const common::IndexMap>,
int>>& maps);
121Mat create_matrix(MPI_Comm comm,
const SparsityPattern& sp,
122 std::optional<std::string> type = std::nullopt);
129MatNullSpace create_nullspace(MPI_Comm comm, std::span<const Vec> basis);
139void set(std::string option);
143void set(std::string option,
const T value)
145 if (option[0] !=
'-')
146 option =
'-' + option;
149 ierr = PetscOptionsSetValue(
nullptr, option.c_str(),
150 boost::lexical_cast<std::string>(value).c_str());
152 petsc::error(ierr, __FILE__,
"PetscOptionsSetValue");
156void clear(std::string option);
174 Vector(
const common::IndexMap& map,
int bs);
177 Vector(
const Vector& x) =
delete;
193 Vector(Vec x,
bool inc_ref_count);
199 Vector& operator=(
const Vector& x) =
delete;
202 Vector& operator=(Vector&& x);
209 std::int64_t
size()
const;
212 std::int32_t local_size()
const;
218 MPI_Comm comm()
const;
221 void set_options_prefix(std::string options_prefix);
225 std::string get_options_prefix()
const;
228 void set_from_options();
244 Operator(Mat A,
bool inc_ref_count);
247 Operator(
const Operator& A) =
delete;
250 Operator(Operator&& A);
256 Operator& operator=(
const Operator& A) =
delete;
259 Operator& operator=(Operator&& A);
263 std::array<std::int64_t, 2>
size()
const;
270 Vec create_vector(std::size_t dim)
const;
285class Matrix :
public Operator
292 static auto set_fn(Mat A, InsertMode mode)
294 return [A, mode, cache = std::vector<PetscInt>()](
295 std::span<const std::int32_t> rows,
296 std::span<const std::int32_t> cols,
297 std::span<const PetscScalar> vals)
mutable ->
int
300#ifdef PETSC_USE_64BIT_INDICES
301 cache.resize(rows.size() + cols.size());
302 std::ranges::copy(rows, cache.begin());
303 std::ranges::copy(cols, std::next(cache.begin(), rows.size()));
304 const PetscInt* _rows = cache.data();
305 const PetscInt* _cols = cache.data() + rows.size();
306 ierr = MatSetValuesLocal(A, rows.size(), _rows, cols.size(), _cols,
309 ierr = MatSetValuesLocal(A, rows.size(), rows.data(), cols.size(),
310 cols.data(), vals.data(), mode);
315 petsc::error(ierr, __FILE__,
"MatSetValuesLocal");
326 static auto set_block_fn(Mat A, InsertMode mode)
328 return [A, mode, cache = std::vector<PetscInt>()](
329 std::span<const std::int32_t> rows,
330 std::span<const std::int32_t> cols,
331 std::span<const PetscScalar> vals)
mutable ->
int
334#ifdef PETSC_USE_64BIT_INDICES
335 cache.resize(rows.size() + cols.size());
336 std::ranges::copy(rows, cache.begin());
337 std::ranges::copy(cols, std::next(cache.begin(), rows.size()));
338 const PetscInt* _rows = cache.data();
339 const PetscInt* _cols = cache.data() + rows.size();
340 ierr = MatSetValuesBlockedLocal(A, rows.size(), _rows, cols.size(), _cols,
343 ierr = MatSetValuesBlockedLocal(A, rows.size(), rows.data(), cols.size(),
344 cols.data(), vals.data(), mode);
349 petsc::error(ierr, __FILE__,
"MatSetValuesBlockedLocal");
363 static auto set_block_expand_fn(Mat A,
int bs0,
int bs1, InsertMode mode)
365 return [A, bs0, bs1, mode, cache0 = std::vector<PetscInt>(),
366 cache1 = std::vector<PetscInt>()](
367 std::span<const std::int32_t> rows,
368 std::span<const std::int32_t> cols,
369 std::span<const PetscScalar> vals)
mutable ->
int
372 cache0.resize(bs0 * rows.size());
373 cache1.resize(bs1 * cols.size());
374 for (std::size_t i = 0; i < rows.size(); ++i)
375 for (
int k = 0; k < bs0; ++k)
376 cache0[bs0 * i + k] = bs0 * rows[i] + k;
378 for (std::size_t i = 0; i < cols.size(); ++i)
379 for (
int k = 0; k < bs1; ++k)
380 cache1[bs1 * i + k] = bs1 * cols[i] + k;
382 ierr = MatSetValuesLocal(A, cache0.size(), cache0.data(), cache1.size(),
383 cache1.data(), vals.data(), mode);
386 petsc::error(ierr, __FILE__,
"MatSetValuesLocal");
393 Matrix(MPI_Comm comm,
const SparsityPattern& sp,
394 std::optional<std::string> type = std::nullopt);
400 Matrix(Mat A,
bool inc_ref_count);
403 Matrix(
const Matrix& A) =
delete;
406 Matrix(Matrix&& A) =
default;
412 Matrix& operator=(
const Matrix& A) =
delete;
415 Matrix& operator=(Matrix&& A) =
default;
420 enum class AssemblyType : std::int32_t
431 void apply(AssemblyType type);
434 double norm(Norm norm_type)
const;
440 void set_options_prefix(std::string options_prefix);
444 std::string get_options_prefix()
const;
447 void set_from_options();
457 explicit KrylovSolver(MPI_Comm comm);
462 KrylovSolver(KSP ksp,
bool inc_ref_count);
465 KrylovSolver(
const KrylovSolver& solver) =
delete;
468 KrylovSolver(KrylovSolver&& solver);
474 KrylovSolver& operator=(
const KrylovSolver&) =
delete;
477 KrylovSolver& operator=(KrylovSolver&& solver);
480 void set_operator(
const Mat A);
483 void set_operators(
const Mat A,
const Mat P);
487 int solve(Vec x,
const Vec b,
bool transpose =
false)
const;
491 void set_options_prefix(std::string options_prefix);
495 std::string get_options_prefix()
const;
498 void set_from_options()
const;
Definition SparsityPattern.h:26
int size(MPI_Comm comm)
Definition MPI.cpp:72
constexpr std::array< std::int64_t, 2 > local_range(int rank, std::int64_t N, int size)
Return local range for the calling process, partitioning the global [0, N - 1] range across all ranks...
Definition MPI.h:90
Miscellaneous classes, functions and types.
Definition dolfinx_common.h:8
Linear algebra interface.
Definition sparsitybuild.h:15
auto norm(const V &x, Norm type=Norm::l2)
Definition Vector.h:268