MultigridPreconditionedKrylovSolver.h

class MultigridPreconditionedKrylovSolver

Parent classes:

This class provides a PETScKrylovSolver which uses a MultigridSolver as preconditioner. 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["direct_solver_type"] = "lu";

// create multigrid pcg solver
fmg::MultigridPreconditionedKrylovSolver mgpcg_solver(mg_solver, "cg");

// solve problem
mgpcg_solver.solve();

// plot solution
plot(mgpcg_solver.solution());

Methods:

MultigridPreconditionedKrylovSolver(MultigridSolver& mg_solver, const std::string& method="default")

Constructor. Note the passed mg_sover is modified to work as a preconditioner i.e.:

mg_solver.parameters["error_on_nonconvergence"] = false;
mg_solver.parameters["tolerance_checking"] = false;
mg_solver.parameters["maximum_iterations"] = 1;

are set.

The method sting can be any valid value as for PETScKrylovSolver. The “default” method is “cg” for symmetric problems and “gmres” for non-symmetric problems.

MultigridPreconditionedKrylovSolver(MultigridProblem& problem, const dolfin::Parameters& mg_params=dolfin::empty_parameters, const std::string& method="default")

Construct from given multigird problem info.

MultigridPreconditionedKrylovSolver(dolfin::LinearVariationalProblem& problem, uint nrefine=0, bool symmetric=true, bool uniform=true, bool interpolate_coef=false, bool interpolate_solution=false, const dolfin::Parameters& mg_params=dolfin::empty_parameters, const std::string& method="default")

Construct from given linear problem (adaptively refined). If nrefine is greater zero, the problem will be refined uniformly nrefine times. If the problem was adapted non-uniformly in advance you have to set uniform to false. Parameters for the multigrid solver can be passed with mg_params. For non-symmetric problems set symmetric explicitly to false. The “default” method is “cg” for symmetric problems and “gmres” for non-symmetric problems.

MultigridPreconditionedKrylovSolver(dolfin::NonlinearVariationalProblem& problem, uint nrefine=0, bool symmetric=false, bool uniform=true, bool interpolate_coef=false, bool interpolate_solution=false, const dolfin::Parameters& mg_params=dolfin::empty_parameters, const std::string& method="default")

Construct from given nonlinear problem (adaptively refined). If nrefine is greater zero, the problem will be refined uniformly nrefine times. If the problem was adapted non-uniformly in advance you have to set uniform to false. Parameters for the multigrid solver can be passed with mg_params. For symmetric problems set symmetric explicitly to true. The “default” method is “cg” for symmetric problems and “gmres” for non-symmetric problems.

void init()

Initialize the solver.

bool initialized() const

Return true if solver is initialized.

uint solve(bool x_assumed_zero=true)

Solve the problem (as defined in the mg_solver) The solution may be obtained with mg_solver().solution()

uint solve(dolfin::GenericVector& x, bool x_assumed_zero=true)

Solve the problem (as defined in the mg_solver) The solution is stored to x

uint solve(dolfin::GenericVector& x, const dolfin::GenericVector& b, bool x_assumed_zero)

Solve the problem (as defined in the mg_solver) The solution is stored to x

const boost::shared_ptr<dolfin::Function> solution() const

Return solution u of the problem, when solved without passing a vector to solve().

Previous topic

MultigridLevel.h

Next topic

MultigridPreconditioner.h

This Page