Reference documentation for deal.II version GIT relicensing-1062-gc06da148b8 2024-07-15 19: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 - 2024 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
325 template <typename VectorType = PETScWrappers::VectorBase,
326 typename PMatrixType = PETScWrappers::MatrixBase,
327 typename AMatrixType = PMatrixType>
328# if defined(DEAL_II_HAVE_CXX20) && \
329 !defined(DEAL_II_DOXYGEN_DO_NOT_PARSE_REQUIRES_CLAUSES)
331 std::constructible_from<VectorType, Vec>) &&
333 std::constructible_from<PMatrixType, Mat>) &&
335 std::constructible_from<AMatrixType, Mat>))
336# endif
338 {
339 public:
345 using real_type = PetscReal;
346
351 const MPI_Comm mpi_comm = PETSC_COMM_WORLD);
352
357
363 operator TS() const;
364
368 TS
370
376
380 void
382
387 void
389
402 void
403 set_matrix(PMatrixType &P);
404
410 void
411 set_matrices(AMatrixType &A, PMatrixType &P);
412
418
424
428 unsigned int
430
438 unsigned int
439 solve(VectorType &y);
440
450 unsigned int
451 solve(VectorType &y, PMatrixType &P);
452
463 unsigned int
464 solve(VectorType &y, AMatrixType &A, PMatrixType &P);
465
474 std::function<void(const real_type t,
475 const VectorType &y,
476 const VectorType &y_dot,
477 VectorType &res)>
479
494 std::function<void(const real_type t,
495 const VectorType &y,
496 const VectorType &y_dot,
497 const real_type alpha,
498 AMatrixType &A,
499 PMatrixType &P)>
501
510 std::function<void(const real_type t, const VectorType &y, VectorType &res)>
512
522 std::function<void(const real_type t,
523 const VectorType &y,
524 AMatrixType &A,
525 PMatrixType &P)>
527
539 std::function<void(const real_type t,
540 const VectorType &y,
541 const unsigned int step_number)>
543
563 std::function<void(const real_type t,
564 const VectorType &y,
565 const VectorType &ydot,
566 const real_type alpha)>
568
580 std::function<void(const VectorType &src, VectorType &dst)>
582
614
619 DEAL_II_DEPRECATED_EARLY
620 std::function<void(const real_type t, VectorType &y)> distribute;
621
640 std::function<void(const real_type t, VectorType &y)>
642
650 DEAL_II_DEPRECATED_EARLY
651 std::function<void(const real_type t,
652 const unsigned int step,
653 const VectorType &y,
654 bool &resize)>
656
676 std::function<
677 bool(const real_type t, const unsigned int step, const VectorType &y)>
679
694 std::function<void(const std::vector<VectorType> &all_in,
695 std::vector<VectorType> &all_out)>
697
698 private:
702 TS ts;
703
709
714
720
725
730
736 mutable std::exception_ptr pending_exception;
737
742
750 void
752
760 void
761 setup_algebraic_constraints(const VectorType &y);
762 };
763
764} // namespace PETScWrappers
765
767
768#endif // DEAL_II_WITH_PETSC
769
770#endif
void add_parameters(ParameterHandler &prm)
void reinit(const TimeStepperData &data)
std::exception_ptr pending_exception
Definition petsc_ts.h:736
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:655
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:478
MPI_Comm get_mpi_communicator() const
std::function< void(const real_type t, VectorType &y)> update_constrained_components
Definition petsc_ts.h:641
std::function< void(const real_type t, const VectorType &y, AMatrixType &A, PMatrixType &P)> explicit_jacobian
Definition petsc_ts.h:526
std::function< bool(const real_type t, const unsigned int step, const VectorType &y)> decide_and_prepare_for_remeshing
Definition petsc_ts.h:678
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:707
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:542
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:500
std::function< IndexSet()> algebraic_components
Definition petsc_ts.h:613
std::function< void(const real_type t, VectorType &y)> distribute
Definition petsc_ts.h:620
std::function< void(const real_type t, const VectorType &y, const VectorType &ydot, const real_type alpha)> setup_jacobian
Definition petsc_ts.h:567
SmartPointer< PMatrixType, TimeStepper > P
Definition petsc_ts.h:708
std::function< void(const real_type t, const VectorType &y, VectorType &res)> explicit_function
Definition petsc_ts.h:511
unsigned int solve(VectorType &y)
std::function< void(const VectorType &src, VectorType &dst)> solve_with_jacobian
Definition petsc_ts.h:581
void set_matrices(AMatrixType &A, PMatrixType &P)
PreconditionShell solve_with_jacobian_pc
Definition petsc_ts.h:713
std::function< void(const std::vector< VectorType > &all_in, std::vector< VectorType > &all_out)> interpolate
Definition petsc_ts.h:696
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