From 3b07369dbe1e57728c570e039d9e204e7d8cbc98 Mon Sep 17 00:00:00 2001 From: Stefano Zampini Date: Sun, 16 Apr 2023 16:35:05 +0300 Subject: [PATCH] PETScWrappers: Add .cc for compatibility this is needed when we don't want to spill "petsc/private" headers into user code --- include/deal.II/lac/petsc_compatibility.h | 56 +++++++++ source/lac/CMakeLists.txt | 1 + source/lac/petsc_compatibility.cc | 118 ++++++++++++++++++ .../lac/petsc_parallel_block_sparse_matrix.cc | 7 +- source/lac/petsc_parallel_block_vector.cc | 7 +- 5 files changed, 178 insertions(+), 11 deletions(-) create mode 100644 source/lac/petsc_compatibility.cc diff --git a/include/deal.II/lac/petsc_compatibility.h b/include/deal.II/lac/petsc_compatibility.h index 400903540b..12e04b00e8 100644 --- a/include/deal.II/lac/petsc_compatibility.h +++ b/include/deal.II/lac/petsc_compatibility.h @@ -30,6 +30,8 @@ # include # include # include +# include +# include # include @@ -132,6 +134,60 @@ namespace PETScWrappers { set_matrix_option(matrix, MAT_KEEP_NONZERO_PATTERN, PETSC_TRUE); } + + /** + * Tell PETSc that the status of the vector has changed. + */ + void + petsc_increment_state_counter(Vec v); + + /** + * Tell PETSc that the status of the matrix has changed. + */ + void + petsc_increment_state_counter(Mat A); + + /** + * Tell PETSc nonlinear solver to use matrix free finite differencing (MFFD). + * + * @p mf_operator indicates to use MFFD for the linear system matrix + * but use a user defined matrix for preconditioning purposed. + * + * @p mf indicates to use MFFD for the both the linear system matrix + * and the preconditioning matrix. + */ + void + set_use_matrix_free(SNES snes, const bool mf_operator, const bool mf); + + /** + * Tell PETSc ODE solver to use matrix free finite differencing (MFFD). + * + * @p mf_operator indicates to use MFFD for the linear system matrix + * but use a user defined matrix for preconditioning purposed. + * + * @p mf indicates to use MFFD for the both the linear system matrix + * and the preconditioning matrix. + */ + void + set_use_matrix_free(TS ts, const bool mf_operator, const bool mf); + + /** + * Set final time for ODE integration. + */ + void + ts_set_max_time(TS ts, const PetscReal maxtime); + + /** + * Set maximum number of steps for ODE integration. + */ + void + ts_set_max_steps(TS ts, const PetscInt maxsteps); + + /** + * Return current step number. + */ + unsigned int + ts_get_step_number(TS ts); } // namespace PETScWrappers DEAL_II_NAMESPACE_CLOSE diff --git a/source/lac/CMakeLists.txt b/source/lac/CMakeLists.txt index 8b6cdd0ded..7be2a079a7 100644 --- a/source/lac/CMakeLists.txt +++ b/source/lac/CMakeLists.txt @@ -99,6 +99,7 @@ if(DEAL_II_WITH_PETSC) set(_unity_include_src ${_unity_include_src} petsc_communication_pattern.cc + petsc_compatibility.cc petsc_full_matrix.cc petsc_matrix_base.cc petsc_matrix_free.cc diff --git a/source/lac/petsc_compatibility.cc b/source/lac/petsc_compatibility.cc new file mode 100644 index 0000000000..4660b5c89d --- /dev/null +++ b/source/lac/petsc_compatibility.cc @@ -0,0 +1,118 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +/* + * Rather than using ifdefs everywhere, try to wrap older versions of PETSc + * functions in one place. This file contains functions that need access to + * internal PETSc headers that we don't want to expose to deal.II users + */ +#include + +#ifdef DEAL_II_WITH_PETSC + +# if DEAL_II_PETSC_VERSION_LT(3, 13, 1) +# include +# endif +# if DEAL_II_PETSC_VERSION_LT(3, 8, 0) +# include +# endif +# include + +// Shorthand notation for PETSc error codes. +# define AssertPETSc(code) \ + do \ + { \ + PetscErrorCode ierr = (code); \ + AssertThrow(ierr == 0, ExcPETScError(ierr)); \ + } \ + while (0) + +DEAL_II_NAMESPACE_OPEN + + +namespace PETScWrappers +{ + void + petsc_increment_state_counter(Vec v) + { + AssertPETSc(PetscObjectStateIncrease(reinterpret_cast(v))); + } + + void + petsc_increment_state_counter(Mat A) + { + AssertPETSc(PetscObjectStateIncrease(reinterpret_cast(A))); + } + + void + set_use_matrix_free(SNES snes, bool mf_operator, bool mf) + { +# if DEAL_II_PETSC_VERSION_LT(3, 13, 1) + snes->mf = mf ? PETSC_TRUE : PETSC_FALSE; + snes->mf_operator = mf_operator ? PETSC_TRUE : PETSC_FALSE; +# else + AssertPETSc(SNESSetUseMatrixFree(snes, + mf_operator ? PETSC_TRUE : PETSC_FALSE, + mf ? PETSC_TRUE : PETSC_FALSE)); +# endif + } + + void + set_use_matrix_free(TS ts, const bool mf_operator, const bool mf) + { + SNES snes; + AssertPETSc(TSGetSNES(ts, &snes)); + set_use_matrix_free(snes, mf_operator, mf); + } + + void + ts_set_max_steps(TS ts, const PetscInt maxsteps) + { +# if DEAL_II_PETSC_VERSION_LT(3, 8, 0) + if (maxsteps >= 0) + ts->max_steps = maxsteps; +# else + AssertPETSc(TSSetMaxSteps(ts, maxsteps)); +# endif + } + + void + ts_set_max_time(TS ts, const PetscReal maxtime) + { +# if DEAL_II_PETSC_VERSION_LT(3, 8, 0) + if (maxtime != PETSC_DEFAULT) + ts->max_time = maxtime; +# else + AssertPETSc(TSSetMaxTime(ts, maxtime)); +# endif + } + + unsigned int + ts_get_step_number(TS ts) + { + PetscInt step; +# if DEAL_II_PETSC_VERSION_LT(3, 8, 0) + step = ts->steps; +# else + AssertPETSc(TSGetStepNumber(ts, &step)); +# endif + return static_cast(step); + } + +} // namespace PETScWrappers + +DEAL_II_NAMESPACE_CLOSE + +#endif // DEAL_II_WITH_PETSC diff --git a/source/lac/petsc_parallel_block_sparse_matrix.cc b/source/lac/petsc_parallel_block_sparse_matrix.cc index 05aed80f87..85e740188d 100644 --- a/source/lac/petsc_parallel_block_sparse_matrix.cc +++ b/source/lac/petsc_parallel_block_sparse_matrix.cc @@ -18,9 +18,6 @@ #ifdef DEAL_II_WITH_PETSC -// For PetscObjectStateIncrease -# include - namespace { // A dummy utility routine to create an empty matrix in case we import @@ -271,9 +268,7 @@ namespace PETScWrappers BlockSparseMatrix::compress(VectorOperation::values operation) { BaseClass::compress(operation); - PetscErrorCode ierr = PetscObjectStateIncrease( - reinterpret_cast(petsc_nest_matrix)); - AssertThrow(ierr == 0, ExcPETScError(ierr)); + petsc_increment_state_counter(petsc_nest_matrix); } diff --git a/source/lac/petsc_parallel_block_vector.cc b/source/lac/petsc_parallel_block_vector.cc index 27ff0a5db2..be7c3549ce 100644 --- a/source/lac/petsc_parallel_block_vector.cc +++ b/source/lac/petsc_parallel_block_vector.cc @@ -17,8 +17,7 @@ #ifdef DEAL_II_WITH_PETSC -// For PetscObjectStateIncrease -# include +# include DEAL_II_NAMESPACE_OPEN @@ -122,9 +121,7 @@ namespace PETScWrappers BlockVector::compress(VectorOperation::values operation) { BlockVectorBase::compress(operation); - PetscErrorCode ierr = PetscObjectStateIncrease( - reinterpret_cast(petsc_nest_vector)); - AssertThrow(ierr == 0, ExcPETScError(ierr)); + petsc_increment_state_counter(petsc_nest_vector); } void -- 2.39.5