21 #ifndef _aspect_simulator_solver_stokes_matrix_free_local_smoothing_h 22 #define _aspect_simulator_solver_stokes_matrix_free_local_smoothing_h 28 #include <deal.II/matrix_free/matrix_free.h> 29 #include <deal.II/matrix_free/operators.h> 30 #include <deal.II/matrix_free/fe_evaluation.h> 32 #include <deal.II/multigrid/mg_constrained_dofs.h> 33 #include <deal.II/multigrid/multigrid.h> 34 #include <deal.II/multigrid/mg_transfer_matrix_free.h> 35 #include <deal.II/multigrid/mg_transfer_global_coarsening.templates.h> 36 #include <deal.II/multigrid/mg_tools.h> 37 #include <deal.II/multigrid/mg_coarse.h> 38 #include <deal.II/multigrid/mg_smoother.h> 39 #include <deal.II/multigrid/mg_matrix.h> 53 template <
int dim,
int velocity_degree>
54 class StokesMatrixFreeHandlerLocalSmoothingImplementation:
public StokesMatrixFreeHandler<dim>
63 const Parameters<dim> ¶meters);
78 std::string
name()
const override;
102 StokesSolver::SolverOutputs
105 const bool solve_newton_system,
106 const double last_pressure_normalization_adjustment,
209 MGLevelObject<MatrixFreeStokesOperators::OperatorCellData<dim, GMGNumberType>>
level_cell_data;
void correct_stokes_rhs()
const Simulator< dim > * simulator
MGTransferMF< dim, GMGNumberType > mg_transfer_Schur_complement
::TrilinosWrappers::MPI::BlockVector BlockVector
ABlockMatrixType A_block_matrix
DoFHandler< dim > dof_handler_v
void evaluate_material_model()
MGLevelObject< GMGSchurComplementMatrixType > mg_matrices_Schur_complement
MatrixFreeStokesOperators::OperatorCellData< dim, GMGNumberType > active_cell_data
MGConstrainedDoFs mg_constrained_dofs_Schur_complement
BTBlockOperatorType BT_block
AffineConstraints< double > constraints_v
void build_preconditioner() override
MGConstrainedDoFs mg_constrained_dofs_A_block
::TrilinosWrappers::BlockSparseMatrix BlockSparseMatrix
void parse_parameters(ParameterHandler &prm) override
DoFHandler< dim > dof_handler_p
std::size_t get_constraint_memory_consumption() const override
std::size_t get_dof_handler_memory_consumption() const override
MGLevelObject< GMGABlockMatrixType > mg_matrices_A_block
MGTransferMF< dim, GMGNumberType > mg_transfer_A_block
std::vector< std::shared_ptr< MatrixFree< dim, double > > > matrix_free_objects
std::size_t get_cell_data_memory_consumption() const override
void setup_dofs() override
::MGTransferMatrixFree< dim, NumberType > MGTransferMF
std::string name() const override
void initialize() override
FESystem< dim > fe_projection
AffineConstraints< double > constraints_p
~StokesMatrixFreeHandlerLocalSmoothingImplementation() override=default
MGConstrainedDoFs mg_constrained_dofs_projection
static void declare_parameters(ParameterHandler &prm)
std::size_t get_mg_transfer_memory_consumption() const override
DoFHandler< dim > dof_handler_projection
StokesMatrixType stokes_matrix
StokesSolver::SolverOutputs solve(const LinearAlgebra::BlockSparseMatrix &system_matrix, const LinearAlgebra::BlockVector &system_rhs, const bool solve_newton_system, const double last_pressure_normalization_adjustment, LinearAlgebra::BlockVector &solution_vector) override
StokesMatrixFreeHandlerLocalSmoothingImplementation(Simulator< dim > &simulator, const Parameters< dim > ¶meters)
SchurComplementMatrixType Schur_complement_block_matrix
MGLevelObject< MatrixFreeStokesOperators::OperatorCellData< dim, GMGNumberType > > level_cell_data