]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Also annotate the implementation of functions with 'requires' clauses. 15132/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Sun, 23 Apr 2023 22:05:12 +0000 (16:05 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Sun, 23 Apr 2023 22:21:23 +0000 (16:21 -0600)
include/deal.II/lac/petsc_snes.templates.h
include/deal.II/lac/petsc_ts.templates.h

index dff02f56234acd77a1b1263e6cf1c6d12e3cde44..9d272094b90fef2690c232f1d0cd5314a0c64ded 100644 (file)
@@ -51,6 +51,8 @@ DEAL_II_NAMESPACE_OPEN
 namespace PETScWrappers
 {
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
   NonlinearSolver<VectorType, PMatrixType, AMatrixType>::NonlinearSolver(
     const NonlinearSolverData &data,
     const MPI_Comm &           mpi_comm)
@@ -63,6 +65,8 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
   NonlinearSolver<VectorType, PMatrixType, AMatrixType>::operator SNES() const
   {
     return snes;
@@ -71,8 +75,9 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  SNES
-  NonlinearSolver<VectorType, PMatrixType, AMatrixType>::petsc_snes()
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  SNES NonlinearSolver<VectorType, PMatrixType, AMatrixType>::petsc_snes()
   {
     return snes;
   }
@@ -80,9 +85,10 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  MPI_Comm
-  NonlinearSolver<VectorType, PMatrixType, AMatrixType>::get_mpi_communicator()
-    const
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  MPI_Comm NonlinearSolver<VectorType, PMatrixType, AMatrixType>::
+    get_mpi_communicator() const
   {
     return PetscObjectComm(reinterpret_cast<PetscObject>(snes));
   }
@@ -90,6 +96,8 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
   NonlinearSolver<VectorType, PMatrixType, AMatrixType>::~NonlinearSolver()
   {
     AssertPETSc(SNESDestroy(&snes));
@@ -98,8 +106,9 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  void
-  NonlinearSolver<VectorType, PMatrixType, AMatrixType>::reinit()
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  void NonlinearSolver<VectorType, PMatrixType, AMatrixType>::reinit()
   {
     AssertPETSc(SNESReset(snes));
   }
@@ -107,8 +116,9 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  void
-  NonlinearSolver<VectorType, PMatrixType, AMatrixType>::reinit(
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  void NonlinearSolver<VectorType, PMatrixType, AMatrixType>::reinit(
     const NonlinearSolverData &data)
   {
     reinit();
@@ -154,8 +164,9 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  void
-  NonlinearSolver<VectorType, PMatrixType, AMatrixType>::set_matrix(
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  void NonlinearSolver<VectorType, PMatrixType, AMatrixType>::set_matrix(
     PMatrixType &P)
   {
     this->A = nullptr;
@@ -165,8 +176,9 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  void
-  NonlinearSolver<VectorType, PMatrixType, AMatrixType>::set_matrices(
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  void NonlinearSolver<VectorType, PMatrixType, AMatrixType>::set_matrices(
     AMatrixType &A,
     PMatrixType &P)
   {
@@ -177,8 +189,10 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  unsigned int
-  NonlinearSolver<VectorType, PMatrixType, AMatrixType>::solve(VectorType &x)
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  unsigned int NonlinearSolver<VectorType, PMatrixType, AMatrixType>::solve(
+    VectorType &x)
   {
     const auto snes_function =
       [](SNES, Vec x, Vec f, void *ctx) -> PetscErrorCode {
@@ -393,9 +407,11 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  unsigned int
-  NonlinearSolver<VectorType, PMatrixType, AMatrixType>::solve(VectorType & x,
-                                                               PMatrixType &P)
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  unsigned int NonlinearSolver<VectorType, PMatrixType, AMatrixType>::solve(
+    VectorType & x,
+    PMatrixType &P)
   {
     set_matrix(P);
     return solve(x);
@@ -404,10 +420,12 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  unsigned int
-  NonlinearSolver<VectorType, PMatrixType, AMatrixType>::solve(VectorType & x,
-                                                               AMatrixType &A,
-                                                               PMatrixType &P)
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  unsigned int NonlinearSolver<VectorType, PMatrixType, AMatrixType>::solve(
+    VectorType & x,
+    AMatrixType &A,
+    PMatrixType &P)
   {
     set_matrices(A, P);
     return solve(x);
index 91d3268c64aef26b14604ec6f89f61d4b0f87919..75ae0bb2fd8ec4ccca0f66bf9f8d81e1490c36ee 100644 (file)
@@ -51,6 +51,8 @@ DEAL_II_NAMESPACE_OPEN
 namespace PETScWrappers
 {
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
   TimeStepper<VectorType, PMatrixType, AMatrixType>::TimeStepper(
     const TimeStepperData &data,
     const MPI_Comm &       mpi_comm)
@@ -63,6 +65,8 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
   TimeStepper<VectorType, PMatrixType, AMatrixType>::~TimeStepper()
   {
     AssertPETSc(TSDestroy(&ts));
@@ -71,6 +75,8 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
   TimeStepper<VectorType, PMatrixType, AMatrixType>::operator TS() const
   {
     return ts;
@@ -79,8 +85,9 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  TS
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::petsc_ts()
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  TS TimeStepper<VectorType, PMatrixType, AMatrixType>::petsc_ts()
   {
     return ts;
   }
@@ -88,9 +95,11 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
   inline MPI_Comm
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::get_mpi_communicator()
-    const
+    TimeStepper<VectorType, PMatrixType, AMatrixType>::get_mpi_communicator()
+      const
   {
     return PetscObjectComm(reinterpret_cast<PetscObject>(ts));
   }
@@ -98,8 +107,10 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
   typename TimeStepper<VectorType, PMatrixType, AMatrixType>::real_type
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::get_time()
+    TimeStepper<VectorType, PMatrixType, AMatrixType>::get_time()
   {
     PetscReal      t;
     PetscErrorCode ierr = TSGetTime(ts, &t);
@@ -110,8 +121,10 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
   typename TimeStepper<VectorType, PMatrixType, AMatrixType>::real_type
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::get_time_step()
+    TimeStepper<VectorType, PMatrixType, AMatrixType>::get_time_step()
   {
     PetscReal      dt;
     PetscErrorCode ierr = TSGetTimeStep(ts, &dt);
@@ -122,8 +135,9 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  void
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::reinit()
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  void TimeStepper<VectorType, PMatrixType, AMatrixType>::reinit()
   {
     AssertPETSc(TSReset(ts));
   }
@@ -131,8 +145,9 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  void
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::reinit(
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  void TimeStepper<VectorType, PMatrixType, AMatrixType>::reinit(
     const TimeStepperData &data)
   {
     reinit();
@@ -199,8 +214,10 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  void
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::set_matrix(PMatrixType &P)
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  void TimeStepper<VectorType, PMatrixType, AMatrixType>::set_matrix(
+    PMatrixType &P)
   {
     this->A = nullptr;
     this->P = &P;
@@ -209,8 +226,9 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  void
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::set_matrices(
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  void TimeStepper<VectorType, PMatrixType, AMatrixType>::set_matrices(
     AMatrixType &A,
     PMatrixType &P)
   {
@@ -221,8 +239,10 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  unsigned int
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::solve(VectorType &y)
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  unsigned int TimeStepper<VectorType, PMatrixType, AMatrixType>::solve(
+    VectorType &y)
   {
     const auto ts_ifunction =
       [](TS, PetscReal t, Vec x, Vec xdot, Vec f, void *ctx) -> PetscErrorCode {
@@ -529,9 +549,11 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  unsigned int
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::solve(VectorType & y,
-                                                           PMatrixType &P)
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  unsigned int TimeStepper<VectorType, PMatrixType, AMatrixType>::solve(
+    VectorType & y,
+    PMatrixType &P)
   {
     set_matrix(P);
     return solve(y);
@@ -540,10 +562,12 @@ namespace PETScWrappers
 
 
   template <typename VectorType, typename PMatrixType, typename AMatrixType>
-  unsigned int
-  TimeStepper<VectorType, PMatrixType, AMatrixType>::solve(VectorType & y,
-                                                           AMatrixType &A,
-                                                           PMatrixType &P)
+  DEAL_II_CXX20_REQUIRES((concepts::is_dealii_petsc_vector_type<VectorType> ||
+                          std::constructible_from<VectorType, Vec>))
+  unsigned int TimeStepper<VectorType, PMatrixType, AMatrixType>::solve(
+    VectorType & y,
+    AMatrixType &A,
+    PMatrixType &P)
   {
     set_matrices(A, P);
     return solve(y);

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.