LCPsolve.jl
A solver for linear complementarity problems.
This package provides a Julia implementation of the Matlab routine written by Yuval Tassa. The solver is particularly useful when the problem to be solved is ill-conditioned. This is often the case, for example, when the linear system arises from a discretization of a Hamilton-Jacobi-Bellman variational inequality in the process of value function iteration. Illustrative applications in economics for solving optimal stopping problems can be found here.
Please see NLsolve.jl for solvers suitable for alternative scenarios.
Usage
An object of type LCP
is used to specify the problem.
Passing the object to solve!
yields the results,
which is returned in an object of type SolverResults
.
The solution is stored in the field sol
.
using LCPsolve, SparseArrays
n = 10000
M = spdiagm(0=>[9;fill(17,n-2);9], -1=>fill(-8,n-1), 1=>fill(-8,n-1))
q = -log.(collect(LinRange(0.05,5,n)))
result = solve!(LCP(M,q))
Setting up the Problem
The complementarity problem can be specified
by constructing an object of type LCP
.
By default, calling LCP(M, q)
creates a problem for finding a vector x
such that
for a square matrix M
and a vector q
the following three equations hold simultaneously:
x >= 0;
Mx + q >= 0;
x'(Mx + q) = 0.
It is advisable to construct M
as a sparse matrix of type SparseMatrixCSC
.
By providing LCP
the keyword arguments l
and u
,
which represent constraints as vectors,
one can solve the more general problem
in which each element of x
indexed by i
satisfies the following relations:
l[i] < x[i] < u[i] => (Mx)[i] + q[i] = 0;
(Mx)[i] + q[i] < 0 => x[i] = u[i];
(Mx)[i] + q[i] > 0 => x[i] = l[i].
Solving the Problem
To solve the problem, pass the LCP
object to solve!
.
The method solve!
accepts an optional argument x0
as the initial starting point for the solver.
Keyword arguments can be passed for adjusting the behavior of the solver.
For details, please use the help mode in REPL.
An object of type SolverResults
will be returned by solve!
.
The solution is stored in the field sol
.
The field converged
indicates whether convergence has been reached.
Applications in Optimal Stopping Problems
The examples
folder contains more illustrations for the usage of the solver.
References
Fischer, A. (1995). A Newton-type method for positive-semidefinite linear complementarity problems. Journal of Optimization Theory and Applications, 86(3), 585-608.
Bazaraa, M. S., Sherali, H. D., & Shetty, C. M. (2013). Nonlinear programming: theory and algorithms. John Wiley & Sons.
Tassa, Y. (2008). LCP / MCP solver (Newton-based). MATLAB Central File Exchange. Retrieved from https://www.mathworks.com/matlabcentral/fileexchange/20952-lcp-mcp-solver-newton-based.