TAOLinearBoundSolver.h

Note

The documentation on this page was automatically extracted from the DOLFIN C++ code and may need to be edited or expanded.

class TAOLinearBoundSolver

Parent class(es)

This class provides a bound constrained solver for a linear variational inequality defined by a matrix A and a vector b. It solves the problem:

Find \(x_l\leq x\leq x_u\) such that \((Ax-b)\cdot (y-x)\geq 0,\; \forall x_l\leq y\leq x_u\)

It is a wrapper for the TAO bound constrained solver.

@code{.py}

# Assemble the linear system A, b = assemble_system(a, L, bc) # Define the constraints constraint_u = Constant(1.) constraint_l = Constant(0.) u_min = interpolate(constraint_l, V) u_max = interpolate(constraint_u, V) # Define the function to store the solution usol=Function(V) # Create the TAOLinearBoundSolver solver=TAOLinearBoundSolver(“tao_gpcg”,”gmres”) # Set some parameters solver.parameters[“monitor_convergence”]=True solver.parameters[“report”]=True # Solve the problem solver.solve(A, usol.vector(), b , u_min.vector(), u_max.vector()) info(solver.parameters,True)

@endcode

explicit TAOLinearBoundSolver(MPI_Comm comm)

Create TAO bound constrained solver

TAOLinearBoundSolver(const std::string method = "default", const std::string ksp_type = "default", const std::string pc_type = "default")

Create TAO bound constrained solver

std::size_t solve(const GenericMatrix &A, GenericVector &x, const GenericVector &b, const GenericVector &xl, const GenericVector &xu)

Solve the linear variational inequality defined by A and b with xl =< x <= xu

std::size_t solve(const PETScMatrix &A, PETScVector &x, const PETScVector &b, const PETScVector &xl, const PETScVector &xu)

Solve the linear variational inequality defined by A and b with xl =< x <= xu

void set_solver(const std::string&)

Set the TAO solver type

void set_ksp(const std::string ksp_type = "default")

Set PETSC Krylov Solver (ksp) used by TAO

Tao tao() const

Return TAO solver pointer

static std::map<std::string, std::string> methods()

Return a list of available Tao solver methods

static std::map<std::string, std::string> krylov_solvers()

Return a list of available krylov solvers

static std::map<std::string, std::string> preconditioners()

Return a list of available preconditioners

static Parameters default_parameters()

Default parameter values

std::shared_ptr<const PETScMatrix> get_matrix() const

Return Matrix shared pointer

std::shared_ptr<const PETScVector> get_vector() const

Return load vector shared pointer