Parent classes:
- dolfin::GenericLinearSolver
- dolfin::PETScObject
This class provides an easy to use multigrid solver for Dolfin. Example usage for 3 refinements:
// ... create coarse bilinearform a, linearform L, solution u and boundary conditions bcs
LinearVariationalProblem problem(a, L, u, bcs);
// create multigrid solver
fmg::MultigridSolver mg_solver(problem, 3);
// set some parameters (optional)
mg_solver.parameters["prolongation_assembly_mode"] = "uniform";
mg_solver.parameters["coarse_solver_type"] = "lu";
mg_solver.parameters("lu_solver_parameters")["same_nonzero_pattern"] = true;
// solve problem
mg_solver.solve();
// plot solution
plot(mg_solver.solution());
Note
For performance reasons it is mostly recommended to use MultigridPreconditionedKrylovSolver instead of using this solver directly.
The following parameters may be configured for the solver:
Parameter Default Description pre_smoother “fsor” Specifies the pre-smoother for each level. See SmoothingOperator for a list of available smoothers. pre_smoother_relax 1.0 Specifies the pre-smoother relaxation factor for each level. pre_smoother_reps 1 Specifies the number of pre-smoother repetitions for each level. post_smoother “bsor” Specifies the post-smoother for each level. See SmoothingOperator for a list of available smoothers. post_smoother_relax 1.0 Specifies the post-smoother relaxation factor for each level. post_smoother_reps 1 Specifies the number of post-smoother repetitions for each level. cycle_pattern “V” Specifies the cycle pattern. See VectorCyclePattern for valid pattern values. coarse_solver_type “lu” Specifies the type of the coarse solver. Valid values are “lu” and “krylov”. lu_solver_method “petsc” Specifies the name of the coarse solver when coarse_solver_type is set to “lu”. lu_solver_parameters Parameters Specifies the parameters of the coarse solver when coarse_solver_type is set to “lu”. See Dolfin manual for valid parameters. krylov_solver_method “cg” for symmetric problems, “gmres” else Specifies the name of the Krylov solver when coarse_solver_type is set to “krylov”. krylov_solver_preconditioner “none” Specifies the name of the Krylov solver preconditioner when coarse_solver_type is set to “krylov”. lu_krylov_parameters Parameters Specifies the parameters of the Krylov solver when coarse_solver_type is set to “krylov”. See Dolfin manual for valid parameters. prolongation_assembly_method “uniform” for uniform problems, “cellwise” else See ProlongationAssembler. prolongation_assembly_eps 1e-13 See ProlongationAssembler prolongation_assembly_space_grouping true See ProlongationAssembler prolongation_assembly_hashing false See ProlongationAssembler dirichlet_projection_method “identify” If set to “remove”, Dirichlet boundary conditions will be removed from the linear system. If set to “identify”, the linear system is modified such that Dirichlet boundary conditions are identifies with the right hand side. dirichlet_dof_identification_method “topological” Method of Dirichlet dof identification. Valid values are: “topological”, “pointwise”, “geometric”. Refer to the Dolfin documentation for more information. tolerance_checking true If set to true the absolute_tolerance and relative_tolerance values are used as stopping criterion. absolute_tolerance 1e-15 The absolute tolerance used as stopping criterion. If tolerance_checking is enabled the solver will stop if the current defect in the defect_norm norm is less than absolute_tolerance. relative_tolerance 1e-6 The relative tolerance used as stopping criterion. If tolerance_checking is enabled the solver will stop if the current defect in the defect_norm norm is less than relative_tolerance times the initial defect. maximum_iterations 10000 The maximum number of iterations used as failure criterion. If error_on_nonconvergence is enabled the solver will throw an exception if the number of iterations exceeds this value. If error_on_nonconvergence is not enabled the solver will return without an error. monitor_convergence false Print residuals for each iteration. error_on_nonconvergence true See maximum_iterations parameter. galerkin_lhs_assembly true If set to true the lhs of the problems on the coarser levels are not assembled but computed from the lhs of the finest level using .
galerkin_rhs_assembly true If set to true the rhs of the problems on the coarser levels are not assembled but computed from the rhs of the finest level using .
coarse_levels_rhs_assembly false If set to true the rhs of the problems on the coarser levels are assembled. Usually the rhs on the coarser levels are not needed. level_offset 0 Specifies which level is considered as the finest level for which the problem is solved. keep_A0 false If set to true the A0 matrix will be kept for each level. keep_P0 false If set to true the P0 matrix will be kept for each level.
Methods:
- MultigridSolver(MultigridProblem& problem)¶
Construct from given multigird problem.
- MultigridSolver(dolfin::LinearVariationalProblem& problem, uint nrefine=0, bool symmetric=true, bool uniform=true, bool interpolate_coef=false, bool interpolate_solution=false)¶
Construct from given linear problem (adaptively refined). If nrefine is greater zero, the problem will be refined uniformly nrefine times. For non-symmetric problems set symmetric explicitly to false. If the problem was adapted non-uniformly in advance you have to set uniform to false.
- MultigridSolver(dolfin::NonlinearVariationalProblem& problem, uint nrefine=0, bool symmetric=false, bool uniform=true, bool interpolate_coef=false, bool interpolate_solution=false)¶
Construct from given nonlinear problem (adaptively refined). If nrefine is greater zero, the problem will be refined uniformly nrefine times. For symmetric problems set symmetric explicitly to true. If the problem was adapted non-uniformly in advance you have to set uniform to false.
- void init()¶
Initialize the solver.
- bool initialized() const¶
Return true if solver is initialized.
- void update(UpdateFlag flags=UpdateFlags::all)¶
Update some aspects of the solver or the problem. This needs to be called to make changes of the parameters object take effect. The flags parameter may be a combination of bitwise ORed UpdateFlags. This is also used when solving a nonlinear problem to update the system each Newton step.
- void update_level(uint i, UpdateFlag flags=UpdateFlags::all)¶
Update single level.
- void update_solution()¶
Update the solution given by get_solution() to all coarser levels (not implemented yet).
- void update_solution(const dolfin::Function& u)¶
Update the solution of the problem to function u to all coarser levels (not implemented yet).
- void update_solution(const dolfin::GenericVector& xr)¶
Update the solution of the problem to reduced vector xr to all coarser levels (not implemented yet).
- uint solve(bool x_assumed_zero=true)¶
Solve problem. The solution may be obtained by get_solution().
- uint solve(dolfin::Function& u, bool u_assumed_zero=true)¶
Solve problem. The solution will be stored in u.
- uint solve(dolfin::GenericVector& xr, bool x_assumed_zero=true)¶
Solve problem. The solution will be stored reduced vector xr.
- uint solve(dolfin::GenericVector& xr, const dolfin::GenericVector& br, bool x_assumed_zero)¶
Solve problem using reduced right hand side br. The solution will be stored reduced vector xr.
Not supported.
Not supported.
- uint solve(const dolfin::GenericLinearOperator& A, dolfin::GenericVector& x, const dolfin::GenericVector& b)¶
Not supported.
- const boost::shared_ptr<const dolfin::GenericMatrix> problem_operator() const¶
Return operator A of the problem.
- const boost::shared_ptr<const dolfin::GenericVector> problem_rhs() const¶
Return vector b of the problem.
- const boost::shared_ptr<dolfin::Function> solution() const¶
Return solution u of the problem.
- dolfin::GenericLinearAlgebraFactory& linear_algebra_factory() const¶
User may return a reference to the linear algebra factory being used. This is currently PETSc and should not be changed.
- bool reduce_required() const¶
Returns true if solve() requires a reduced vector.
- boost::shared_ptr<dolfin::GenericVector> create_reduced_vector() const¶
Create a vector suitable for solve() (vector without boundary dofs)
- void reduce_vector(const dolfin::GenericVector& x, dolfin::GenericVector& xr) const¶
Copy vector x to xr with boundary dofs removed.
- void expand_vector(const dolfin::GenericVector& xr, dolfin::GenericVector& x) const¶
Copy vector xr to x with boundary values inserted.
- const MultigridProblem& problem() const¶
Return multigrid problem object.
- MultigridLevel& finest_level() const¶
Return finest level of the problem.
- static const dolfin::Parameters& default_parameters(bool linear=true, bool symmetric=true, bool uniform=true)¶
Return default parameters.
- void before_solve()¶
This is function is called at the beginning of solve.
- void after_solve()¶
This is function is called before return of solve.
- bool test_residual_break_condition(const MultigridLevel& level, const GenericCyclePattern& cycle_pattern)¶
User needs to return true if converged, false to continue computation or throw an exception if not converged. This is function is called only on the finest level after computing the defect.
- bool test_loop_break_condition(const GenericCyclePattern& cycle_pattern)¶
User needs to return true if converged, false to continue computation or throw an exception if not converged. This is function is called at the end of every loop.
- void init_smoothers(MultigridLevel& level)¶
User may set level.pre_smoother and level.post_smoother to modify smoothers.
- void init_coarse_solver(MultigridLevel& level)¶
User may set level.coarse_solver to modify coarse solver. This requires also to call level.coarse_solver->set_operator(level.A)
- MultigridLevel* create_level()¶
User may return a new instance of a custom subclass of MultigridLevel.
User may set the shared pointer to a custom cycle pattern.
- void init_user(MultigridLevel& level)¶
User defined initialization routine will be called last (does currently nothing).
- void before_cycle_action(MultigridLevel& level, GenericCyclePattern& cycle)¶
Called each time before a cycle action (does currently nothing).
- void next_cycle_action(MultigridLevel& level, GenericCyclePattern& cycle)¶
Called to increment the cycle action (calls cycle.next()).
- void pr_transition(MultigridLevel& level)¶
Called on transition between prolongation and restriction (does currently nothing).