21 #ifndef _aspect_simulator_solver_stokes_matrix_free_h 22 #define _aspect_simulator_solver_stokes_matrix_free_h 30 #include <deal.II/multigrid/mg_transfer_matrix_free.h> 31 #include <deal.II/multigrid/mg_transfer_global_coarsening.templates.h> 33 #include <deal.II/lac/solver_bicgstab.h> 113 template <
class StokesMatrixType,
class ABlockMatrixType,
class BTBlockOperatorType,
114 class SchurComplementMatrixType,
class ABlockPreconditionerType,
115 class SchurComplementPreconditionerType>
139 const ABlockMatrixType &A_block,
140 const BTBlockOperatorType &BT_block,
141 const SchurComplementMatrixType &Schur_complement_block,
142 const ABlockPreconditionerType &A_block_preconditioner,
143 const SchurComplementPreconditionerType &Schur_complement_preconditioner,
144 const bool do_solve_A,
145 const bool do_solve_Schur_complement,
146 const bool A_block_is_symmetric,
147 const double A_block_tolerance,
148 const double Schur_complement_tolerance);
153 void vmult (::LinearAlgebra::distributed::BlockVector<double> &dst,
154 const ::LinearAlgebra::distributed::BlockVector<double> &src)
const;
156 unsigned int n_iterations_A_block()
const;
157 unsigned int n_iterations_Schur_complement()
const;
182 mutable ::LinearAlgebra::distributed::BlockVector<double>
utmp;
185 template <
class StokesMatrixType,
class ABlockMatrixType,
class BTBlockOperatorType,
class SchurComplementMatrixType,
186 class ABlockPreconditionerType,
class SchurComplementPreconditionerType>
188 ABlockPreconditionerType, SchurComplementPreconditionerType>::
189 BlockSchurGMGPreconditioner (
const StokesMatrixType &Stokes_matrix,
190 const ABlockMatrixType &A_block,
191 const BTBlockOperatorType &BT_block,
192 const SchurComplementMatrixType &Schur_complement_block,
193 const ABlockPreconditionerType &A_block_preconditioner,
194 const SchurComplementPreconditionerType &Schur_complement_preconditioner,
195 const bool do_solve_A,
196 const bool do_solve_Schur_complement,
197 const bool A_block_symmetric,
198 const double A_block_tolerance,
199 const double Schur_complement_tolerance)
201 stokes_matrix (Stokes_matrix),
204 Schur_complement_block (Schur_complement_block),
205 A_block_preconditioner (A_block_preconditioner),
206 Schur_complement_preconditioner (Schur_complement_preconditioner),
207 do_solve_A (do_solve_A),
208 do_solve_Schur_complement (do_solve_Schur_complement),
209 A_block_is_symmetric (A_block_symmetric),
211 n_iterations_Schur_complement_ (0),
212 A_block_tolerance (A_block_tolerance),
213 Schur_complement_tolerance (Schur_complement_tolerance)
216 template <
class StokesMatrixType,
class ABlockMatrixType,
class BTBlockOperatorType,
217 class SchurComplementMatrixType,
class ABlockPreconditionerType,
218 class SchurComplementPreconditionerType>
221 ABlockPreconditionerType, SchurComplementPreconditionerType>
:: 227 template <
class StokesMatrixType,
class ABlockMatrixType,
class BTBlockOperatorType,
class SchurComplementMatrixType,
228 class ABlockPreconditionerType,
class SchurComplementPreconditionerType>
231 ABlockPreconditionerType, SchurComplementPreconditionerType>
:: 237 template <
class StokesMatrixType,
class ABlockMatrixType,
class BTBlockOperatorType,
class SchurComplementMatrixType,
238 class ABlockPreconditionerType,
class SchurComplementPreconditionerType>
241 ABlockPreconditionerType, SchurComplementPreconditionerType>
:: 242 vmult (::LinearAlgebra::distributed::BlockVector<double> &dst,
243 const ::LinearAlgebra::distributed::BlockVector<double> &src)
const 262 SolverCG<::LinearAlgebra::distributed::Vector<double>> solver(solver_control);
269 if (src.block(1).l2_norm() > 1e-50)
274 dst.block(1), src.block(1),
281 catch (
const std::exception &exc)
284 "BlockSchurGMGPreconditioner::vmult",
285 std::vector<SolverControl> {solver_control},
287 src.block(0).get_mpi_communicator());
297 dst.block(1) *= -1.0;
301 utmp.block(0) *= -1.0;
302 utmp.block(0) += src.block(0);
311 PrimitiveVectorMemory<::LinearAlgebra::distributed::Vector<double>> mem;
317 SolverCG<::LinearAlgebra::distributed::Vector<double>> solver(solver_control,mem);
327 SolverBicgstab<::LinearAlgebra::distributed::Vector<double>>
328 solver(solver_control,
330 SolverBicgstab<::LinearAlgebra::distributed::Vector<double>>::AdditionalData(
false));
340 catch (
const std::exception &exc)
343 "BlockSchurGMGPreconditioner::vmult",
344 std::vector<SolverControl> {solver_control},
346 src.block(0).get_mpi_communicator());
virtual std::size_t get_dof_handler_memory_consumption() const =0
const SchurComplementMatrixType & Schur_complement_block
const Simulator< dim > * simulator
void throw_linear_solver_failure_exception(const std::string &solver_name, const std::string &function_name, const std::vector< SolverControl > &solver_controls, const std::exception &exc, const MPI_Comm mpi_communicator, const std::string &output_filename="")
const bool A_block_is_symmetric
const SchurComplementPreconditionerType & Schur_complement_preconditioner
virtual std::size_t get_mg_transfer_memory_consumption() const =0
void vmult(::LinearAlgebra::distributed::BlockVector< double > &dst, const ::LinearAlgebra::distributed::BlockVector< double > &src) const
virtual std::size_t get_cell_data_memory_consumption() const =0
std::unique_ptr< StokesMatrixFreeHandler< dim > > create_matrix_free_solver(Simulator< dim > &simulator, const Parameters< dim > ¶meters)
unsigned int n_iterations_A_block() const
mutable ::LinearAlgebra::distributed::BlockVector< double > utmp
unsigned int n_iterations_Schur_complement_
virtual void assemble()=0
const ABlockMatrixType & A_block
virtual void build_preconditioner()=0
const ABlockPreconditionerType & A_block_preconditioner
unsigned int n_iterations_Schur_complement() const
const double Schur_complement_tolerance
const StokesMatrixType & stokes_matrix
const bool do_solve_Schur_complement
unsigned int n_iterations_A_
const double A_block_tolerance
static void declare_parameters(ParameterHandler &prm)
virtual void setup_dofs()=0
virtual std::size_t get_constraint_memory_consumption() const =0
const BTBlockOperatorType & BT_block