Reference documentation for deal.II version GIT relicensing-487-ge9eb5ab491 2024-04-25 07:20:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
petsc_ts.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_petsc_ts_h
16#define dealii_petsc_ts_h
17
18#include <deal.II/base/config.h>
19
20#ifdef DEAL_II_WITH_PETSC
21# include <deal.II/base/mpi.h>
24
28
29# include <petscts.h>
30
31# include <exception>
32
33# if defined(DEAL_II_HAVE_CXX20)
34# include <concepts>
35# endif
36
38
39namespace PETScWrappers
40{
45 {
46 public:
52 using real_type = PetscReal;
53
89 // Running parameters
90 const std::string &options_prefix = "",
91 const std::string &ts_type = "",
92 const real_type initial_time = 0.0,
93 const real_type final_time = 0.0,
94 const real_type initial_step_size = 0.0,
95 const int max_steps = -1,
96 const bool match_step = false,
97 const bool restart_if_remesh = false,
98 // Error parameters
99 const std::string &ts_adapt_type = "none",
100 const real_type minimum_step_size = -1.0,
101 const real_type maximum_step_size = -1.0,
102 const real_type absolute_tolerance = -1.0,
103 const real_type relative_tolerance = -1.0,
104 const bool ignore_algebraic_lte = true)
119 {}
120
124 void
126
130 std::string options_prefix;
131
135 std::string ts_type;
136
141
146
153
160
165
170
174 std::string ts_adapt_type;
175
182
189
196
203
208 };
209
322 template <typename VectorType = PETScWrappers::VectorBase,
323 typename PMatrixType = PETScWrappers::MatrixBase,
324 typename AMatrixType = PMatrixType>
325# if defined(DEAL_II_HAVE_CXX20) && \
326 !defined(DEAL_II_DOXYGEN_DO_NOT_PARSE_REQUIRES_CLAUSES)
328 std::constructible_from<VectorType, Vec>) &&
330 std::constructible_from<PMatrixType, Mat>) &&
332 std::constructible_from<AMatrixType, Mat>))
333# endif
335 {
336 public:
342 using real_type = PetscReal;
343
348 const MPI_Comm mpi_comm = PETSC_COMM_WORLD);
349
354
360 operator TS() const;
361
365 TS
367
373
377 void
379
384 void
386
399 void
400 set_matrix(PMatrixType &P);
401
407 void
408 set_matrices(AMatrixType &A, PMatrixType &P);
409
415
421
425 unsigned int
427
435 unsigned int
436 solve(VectorType &y);
437
447 unsigned int
448 solve(VectorType &y, PMatrixType &P);
449
460 unsigned int
461 solve(VectorType &y, AMatrixType &A, PMatrixType &P);
462
471 std::function<void(const real_type t,
472 const VectorType &y,
473 const VectorType &y_dot,
474 VectorType &res)>
476
491 std::function<void(const real_type t,
492 const VectorType &y,
493 const VectorType &y_dot,
494 const real_type alpha,
495 AMatrixType &A,
496 PMatrixType &P)>
498
507 std::function<void(const real_type t, const VectorType &y, VectorType &res)>
509
519 std::function<void(const real_type t,
520 const VectorType &y,
521 AMatrixType &A,
522 PMatrixType &P)>
524
536 std::function<void(const real_type t,
537 const VectorType &y,
538 const unsigned int step_number)>
540
560 std::function<void(const real_type t,
561 const VectorType &y,
562 const VectorType &ydot,
563 const real_type alpha)>
565
577 std::function<void(const VectorType &src, VectorType &dst)>
579
589
603 std::function<void(const real_type t, VectorType &y)> distribute;
604
621 std::function<void(const real_type t,
622 const unsigned int step,
623 const VectorType &y,
624 bool &resize)>
626
641 std::function<void(const std::vector<VectorType> &all_in,
642 std::vector<VectorType> &all_out)>
644
645 private:
649 TS ts;
650
656
661
667
672
677
683 mutable std::exception_ptr pending_exception;
684
689
697 void
699
707 void
708 setup_algebraic_constraints(const VectorType &y);
709 };
710
711} // namespace PETScWrappers
712
714
715#endif // DEAL_II_WITH_PETSC
716
717#endif
void add_parameters(ParameterHandler &prm)
void reinit(const TimeStepperData &data)
std::exception_ptr pending_exception
Definition petsc_ts.h:683
unsigned int solve(VectorType &y, PMatrixType &P)
std::function< void(const real_type t, const unsigned int step, const VectorType &y, bool &resize)> decide_for_coarsening_and_refinement
Definition petsc_ts.h:625
void setup_algebraic_constraints(const VectorType &y)
std::function< void(const real_type t, const VectorType &y, const VectorType &y_dot, VectorType &res)> implicit_function
Definition petsc_ts.h:475
MPI_Comm get_mpi_communicator() const
std::function< void(const real_type t, const VectorType &y, AMatrixType &A, PMatrixType &P)> explicit_jacobian
Definition petsc_ts.h:523
unsigned int get_step_number()
TimeStepper(const TimeStepperData &data=TimeStepperData(), const MPI_Comm mpi_comm=PETSC_COMM_WORLD)
SmartPointer< AMatrixType, TimeStepper > A
Definition petsc_ts.h:654
unsigned int solve(VectorType &y, AMatrixType &A, PMatrixType &P)
std::function< void(const real_type t, const VectorType &y, const unsigned int step_number)> monitor
Definition petsc_ts.h:539
std::function< void(const real_type t, const VectorType &y, const VectorType &y_dot, const real_type alpha, AMatrixType &A, PMatrixType &P)> implicit_jacobian
Definition petsc_ts.h:497
std::function< IndexSet()> algebraic_components
Definition petsc_ts.h:588
std::function< void(const real_type t, VectorType &y)> distribute
Definition petsc_ts.h:603
std::function< void(const real_type t, const VectorType &y, const VectorType &ydot, const real_type alpha)> setup_jacobian
Definition petsc_ts.h:564
SmartPointer< PMatrixType, TimeStepper > P
Definition petsc_ts.h:655
std::function< void(const real_type t, const VectorType &y, VectorType &res)> explicit_function
Definition petsc_ts.h:508
unsigned int solve(VectorType &y)
std::function< void(const VectorType &src, VectorType &dst)> solve_with_jacobian
Definition petsc_ts.h:578
void set_matrices(AMatrixType &A, PMatrixType &P)
PreconditionShell solve_with_jacobian_pc
Definition petsc_ts.h:660
std::function< void(const std::vector< VectorType > &all_in, std::vector< VectorType > &all_out)> interpolate
Definition petsc_ts.h:643
void set_matrix(PMatrixType &P)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
TimeStepperData(const std::string &options_prefix="", const std::string &ts_type="", const real_type initial_time=0.0, const real_type final_time=0.0, const real_type initial_step_size=0.0, const int max_steps=-1, const bool match_step=false, const bool restart_if_remesh=false, const std::string &ts_adapt_type="none", const real_type minimum_step_size=-1.0, const real_type maximum_step_size=-1.0, const real_type absolute_tolerance=-1.0, const real_type relative_tolerance=-1.0, const bool ignore_algebraic_lte=true)
Definition petsc_ts.h:88