ASPECT
stokes_matrix_free_local_smoothing.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2025 by the authors of the ASPECT code.
3 
4  This file is part of ASPECT.
5 
6  ASPECT is free software; you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation; either version 2, or (at your option)
9  any later version.
10 
11  ASPECT is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with ASPECT; see the file LICENSE. If not see
18  <http://www.gnu.org/licenses/>.
19 */
20 
21 #ifndef _aspect_simulator_solver_stokes_matrix_free_local_smoothing_h
22 #define _aspect_simulator_solver_stokes_matrix_free_local_smoothing_h
23 
24 #include <aspect/global.h>
27 
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>
31 
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>
40 
41 namespace aspect
42 {
53  template <int dim, int velocity_degree>
54  class StokesMatrixFreeHandlerLocalSmoothingImplementation: public StokesMatrixFreeHandler<dim>
55  {
56  public:
63  const Parameters<dim> &parameters);
64 
69 
73  void initialize() override;
74 
78  std::string name() const override;
79 
102  StokesSolver::SolverOutputs
103  solve(const LinearAlgebra::BlockSparseMatrix &system_matrix,
104  const LinearAlgebra::BlockVector &system_rhs,
105  const bool solve_newton_system,
106  const double last_pressure_normalization_adjustment,
107  LinearAlgebra::BlockVector &solution_vector) override;
108 
113  void setup_dofs() override;
114 
122  void assemble() override;
123 
128  void build_preconditioner() override;
129 
133  static
134  void declare_parameters (ParameterHandler &prm);
135 
139  void parse_parameters (ParameterHandler &prm) override;
140 
144  std::size_t get_dof_handler_memory_consumption() const override;
145 
149  std::size_t get_mg_transfer_memory_consumption() const override;
150 
154  std::size_t get_constraint_memory_consumption() const override;
155 
160  std::size_t get_cell_data_memory_consumption() const override;
161 
162  private:
169 
174  void correct_stokes_rhs();
175 
176 
178 
180 
186 
192 
193  DoFHandler<dim> dof_handler_v;
194  DoFHandler<dim> dof_handler_p;
195  DoFHandler<dim> dof_handler_projection;
196 
197  FESystem<dim> fe_v;
198  FESystem<dim> fe_p;
199  FESystem<dim> fe_projection;
200 
205 
209  MGLevelObject<MatrixFreeStokesOperators::OperatorCellData<dim, GMGNumberType>> level_cell_data;
210 
217 
222 
223  AffineConstraints<double> constraints_v;
224  AffineConstraints<double> constraints_p;
225 
226  MGLevelObject<GMGABlockMatrixType> mg_matrices_A_block;
227  MGLevelObject<GMGSchurComplementMatrixType> mg_matrices_Schur_complement;
228 
229  MGConstrainedDoFs mg_constrained_dofs_A_block;
232 
235 
236  std::vector<std::shared_ptr<MatrixFree<dim,double>>> matrix_free_objects;
237  };
238 }
239 
240 #endif
const Simulator< dim > * simulator
::TrilinosWrappers::MPI::BlockVector BlockVector
Definition: global.h:269
MGLevelObject< GMGSchurComplementMatrixType > mg_matrices_Schur_complement
MatrixFreeStokesOperators::OperatorCellData< dim, GMGNumberType > active_cell_data
::TrilinosWrappers::BlockSparseMatrix BlockSparseMatrix
Definition: global.h:280
void parse_parameters(ParameterHandler &prm) override
double GMGNumberType
std::vector< std::shared_ptr< MatrixFree< dim, double > > > matrix_free_objects
::MGTransferMatrixFree< dim, NumberType > MGTransferMF
Definition: compat.h:64
static void declare_parameters(ParameterHandler &prm)
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 > &parameters)
MGLevelObject< MatrixFreeStokesOperators::OperatorCellData< dim, GMGNumberType > > level_cell_data