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.