A Stokes preconditioner (saddle-point problem) for the problem \(\begin{pmatrix} A & B \\ C & 0 \end{pmatrix} \begin{pmatrix} u \\ p \end{pmatrix} = \begin{pmatrix} f \\ g \end{pmatrix}, \). More...
#include <dumux/linear/stokes_solver.hh>
where A is the discrete velocity operator and B and C are discrete gradient and divergence operators. This preconditioner is especially suited for solving the incompressible Stokes equations.
| M | Type of the matrix. |
| X | Type of the update. |
| Y | Type of the defect. |
| l | Preconditioner block level |
Public Types | |
| using | matrix_type = M |
| The matrix type the preconditioner is for. | |
| using | domain_type = X |
| The domain type of the preconditioner. | |
| using | range_type = Y |
| The range type of the preconditioner. | |
| using | field_type = typename X::field_type |
| The field type of the preconditioner. | |
| using | scalar_field_type = Dune::Simd::Scalar<field_type> |
| Scalar type underlying the field_type. | |
| using | PressureLinearOperator = Dune::MatrixAdapter<P,V,V> |
| the type of the pressure operator | |
Public Member Functions | |
| StokesPreconditioner (const std::shared_ptr< const Dune::AssembledLinearOperator< M, X, Y > > &fullLinearOperator, const std::shared_ptr< const PressureLinearOperator > &pressureLinearOperator, const Dune::ParameterTree ¶ms) | |
| Constructor. | |
| void | pre (X &update, Y ¤tDefect) override |
| Prepare the preconditioner. | |
| void | apply (X &update, const Y ¤tDefect) override |
| Apply the preconditioner. | |
| void | post (X &update) override |
| Clean up. | |
| Dune::SolverCategory::Category | category () const override |
| Category of the preconditioner (see SolverCategory::Category) | |
| using Dumux::Detail::StokesPreconditioner< M, X, Y, l >::domain_type = X |
| using Dumux::Detail::StokesPreconditioner< M, X, Y, l >::field_type = typename X::field_type |
| using Dumux::Detail::StokesPreconditioner< M, X, Y, l >::matrix_type = M |
| using Dumux::Detail::StokesPreconditioner< M, X, Y, l >::PressureLinearOperator = Dune::MatrixAdapter<P,V,V> |
| using Dumux::Detail::StokesPreconditioner< M, X, Y, l >::range_type = Y |
| using Dumux::Detail::StokesPreconditioner< M, X, Y, l >::scalar_field_type = Dune::Simd::Scalar<field_type> |
|
inline |
| fullLinearOperator | the Stokes linear operator |
| pressureLinearOperator | the linear operator to be used for the pressure space preconditioner |
| params | a parameter tree for the preconditioner configuration |
|
inlineoverride |
| update | The update to be computed. |
| currentDefect | The current defect. |
The currentDefect has be be in a consistent representation, Definition 2.3 Blatt and Bastian (2009) https://doi.org/10.1504/IJCSE.2008.021112 The update is initially zero. At exit the update has to be in a consistent representation. This usually requires communication.
|
inlineoverride |
|
inlineoverride |
|
inlineoverride |