# Linear Complementarity Problems (LCP)¶

## Problem statement¶

The Linear Complementarity problem (LCP) is defined by

Find $$(z,w)$$ such that:

$\begin{split}\begin{equation*} \begin{cases} M \ z + q = w \\ 0 \le w \perp z \ge 0 \end{cases} \end{equation*}\end{split}$

For more details on theory and analysis of LCP, we refer to .

## Implementation in numerics¶

Structure to define the problem: LinearComplementarityProblem.

The generic driver for all LCPs is linearComplementarity_driver() which switches either to lcp_driver_DenseMatrix() or lcp_driver_SparseMatrix() according to the kind of storage used in the problem.

Solvers list LCP_SOLVER

## Error computation¶

The criterion is based on :

\begin{align}\begin{aligned}error = \frac{1}{\|q\| }\sqrt(\sum_{i} (z_i - (z_i - (Mz+q)_i)^{+})^2)\\error = \frac{1}{\|q\| }\sum_{i} [ (z_i*(Mz+q)_i)^{+} + (z_i)^{-} + {(Mz+q)_i}^{-} ]\end{aligned}\end{align}

with $$x^{+} = max(0,x)$$ and $$x^{-} = max(0,-x)$$.

• lcp_compute_error() returns 0 if $$error \leq tolerance$$, else 1.

• A call to this function updates the content of the input vector w with $$Mz + q$$.

## LCP available solvers¶

### Direct solvers¶

#### Lemke (SICONOS_LCP_LEMKE)¶

direct solver for LCP based on pivoting method principle for degenerate problem: the choice of pivot variable is performed via lexicographic ordering.

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 10000

• iparam[SICONOS_LCP_IPARAM_PIVOTING_METHOD_TYPE] = 0

• dparam[SICONOS_DPARAM_TOL] = 1e-6

• dparam = 0.0

• dparam = 0.0

#### Pivot based methods¶

generic solver for pivot-based methods: Bard, Murty and Lemke rules are implemented.

drivers:

and SICONOS_LCP_MURTY, * lcp_pathsearch() for SICONOS_LCP_PATHSEARCH, * lcp_pivot_lumodfor :enumerator:SICONOS_LCP_PIVOT_LUMOD().

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 10000

• iparam[SICONOS_LCP_IPARAM_PIVOTING_METHOD_TYPE] =

• SICONOS_LCP_PIVOT_BARD for SICONOS_LCP_BARD

• SICONOS_LCP_PIVOT_LEAST_INDEX for SICONOS_LCP_MURTY

• SICONOS_LCP_PIVOT_LEMKE for SICONOS_LCP_PIVOT, SICONOS_LCP_PATHSEARCH and SICONOS_LCP_PIVOT_LUMOD

• iparam[SICONOS_IPARAM_PATHSEARCH_STACKSIZE] = 0;

• dparam[SICONOS_DPARAM_TOL] = 100 * epsilon (machine precision)

#### Enumerative solver (SICONOS_LCP_ENUM)¶

Brute-force method which tries every possible solution.

driver: lcp_enum()

parameters:

• iparam[SICONOS_LCP_IPARAM_ENUM_USE_DGELS] = 0 (0 : use dgesv, 1: use dgels)

• iparam[SICONOS_LCP_IPARAM_ENUM_CURRENT_ENUM] (out): key of the solution

• iparam[SICONOS_LCP_IPARAM_ENUM_SEED] = 0, starting key values

• iparam[SICONOS_LCP_IPARAM_ENUM_MULTIPLE_SOLUTIONS] = 0, search for multiple solutions if 1

• iparam[SICONOS_LCP_IPARAM_ENUM_NUMBER_OF_SOLUTIONS] (out): number of solutions

• dparam[SICONOS_DPARAM_TOL] = 1e-6

#### PATH (SICONOS_LCP_PATH)¶

Works only if Siconos has been built with path support (if PathFerris or PathVI has been found, see :ref:siconos_install_guide)

driver: lcp_path()

parameters * iparam[SICONOS_IPARAM_MAX_ITER] = 1000 * dparam[SICONOS_DPARAM_TOL] = 1e-12

### Iterative solvers¶

#### Conjugated Projected Gradient (SICONOS_LCP_CPG)¶

Conjugated Projected Gradient solver for LCP based on quadratic minimization. Reference: “Conjugate gradient type algorithms for frictional multi-contact problems: applications to granular materials”, M. Renouf, P. Alart. doi:10.1016/j.cma.2004.07.009

driver: lcp_cpg()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 10000

• dparam[SICONOS_DPARAM_TOL] = 1e-6

#### Projected Gauss-Seidel (SICONOS_LCP_PGS)¶

driver: lcp_pgs()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 10000

• dparam[SICONOS_DPARAM_TOL] = 1e-6

#### Regularized Projected Gauss-Seidel (SICONOS_LCP_RPGS)¶

Regularized Projected Gauss-Seidel, is a solver for LCP, able to handle matrices with null diagonal terms.

driver: lcp_rpgs()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 10000

• dparam[SICONOS_DPARAM_TOL] = 1e-6

• dparam[SICONOS_LCP_DPARAM_RHO] = 1.0

#### PSOR (SICONOS_LCP_PSOR)¶

Projected Successive over relaxation solver for LCP. See , Chap 5.

driver: lcp_psor()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 1000

• dparam[SICONOS_DPARAM_TOL] = 1e-6

• dparam[SICONOS_LCP_DPARAM_RHO] = 0.1

#### Latin method (SICONOS_LCP_LATIN and SICONOS_LCP_LATIN_W)¶

Latin stands for LArge Time INcrements. ‘w’ version is the Latin solver with relaxation.

drivers: lcp_latin() and lcp_latin_w()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 1000

• dparam[SICONOS_DPARAM_TOL] = 1e-6

• dparam[SICONOS_LCP_DPARAM_LATIN_PARAMETER] = 0.3

• dparam[SICONOS_LCP_DPARAM_RHO] = 1.0 (only useful for solver with relaxation)

#### Sparse-block Gauss-Seidel (SICONOS_LCP_NSGS)¶

Gauss-Seidel solver based on a Sparse-Block storage for the matrix M of the LCP.

Matrix M of the LCP must be a SparseBlockStructuredMatrix.

This solver first build a local problem for each row of blocks and then call any of the other solvers through lcp_driver().

driver: lcp_nsgs_SBM()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 1000

• iparam[SICONOS_LCP_IPARAM_NSGS_ITERATIONS_SUM] (out): sum of all local number of iterations (if it has sense for the local solver)

• dparam[SICONOS_DPARAM_TOL] = 1e-6

• dparam[SICONOS_LCP_DPARAM_NSGS_LOCAL_ERROR_SUM] (in): sum of all local error values

internal solver : SICONOS_LCP_PSOR

### Equation-based solvers¶

#### Nonsmooth Newton, min formulation (SICONOS_LCP_NEWTONMIN)¶

Nonsmooth Newton method based on the min formulation of the LCP.

\begin{align}\begin{aligned}0 \le z \perp w \ge 0 \Longrightarrow \min(w,\rho z)=0 \Longrightarrow w = \max(0,w - \rho z)\\\begin{split}H(z) = H(\left[ \begin{array}{c} z \\ w \end{array}\right])= \left[ \begin{array}{c} w-Mz-q \\ min(w,\rho z) \end{array}\right] =0\end{split}\end{aligned}\end{align}

References: Alart & Curnier 1990, Pang 1990

driver: lcp_newton_min()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 1000

• dparam[SICONOS_DPARAM_TOL] = 1e-6

#### Nonsmooth Newton, Fisher-Burmeister (SICONOS_LCP_NEWTON_FB_FBLSA)¶

Nonsmooth Newton method based on the Fischer-Bursmeister NCP function. It uses a variant of line search algorithm (VFBLSA in Facchinei-Pang 2003).

\begin{align}\begin{aligned}0 \le z \perp w \ge 0 \Longrightarrow \phi(z,w)=\sqrt{z^2+w^2}-(z+w)=0\\\begin{split}\Phi(z) = \left[ \begin{array}{c} \phi(z_1,w_1) \\ \phi(z_1,w_1) \\ \vdots \\ \phi(z_n,w_n) \end{array}\right] =0\end{split}\end{aligned}\end{align}

References: Alart & Curnier 1990, Pang 1990

driver: lcp_newton_FB()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 1000

• iparam[SICONOS_IPARAM_LSA_NONMONOTONE_LS] = 0 (in): if > 0. use a non-monotone linear search

• iparam[SICONOS_IPARAM_LSA_NONMONOTONE_LS_M] = 0 (in): if a non-monotone linear search is used, specify the number of merit values to remember

• iparam[SICONOS_IPARAM_STOPPING_CRITERION] = SICONOS_STOPPING_CRITERION_USER_ROUTINE;

• dparam[SICONOS_DPARAM_TOL] = 1e-10

• dparam[SICONOS_DPARAM_LSA_ALPHA_MIN] = 1e-16;

#### Nonsmooth Newton, Fisher-Burmeister (SICONOS_LCP_NEWTON_MIN_FBLSA)¶

A nonsmooth Newton method based based on the minFBLSA algorithm : the descent direction is given by a min reformulation but the linesearch is done with Fischer-Burmeister (and if needed the gradient direction).

driver: lcp_newton_minFB()

parameters: same as SICONOS_LCP_NEWTON_FB_FBLSA.

#### GAMS solver (SICONOS_LCP_GAMS)¶

Optimization solvers from GAMS.

Works only if Siconos has been built with GAMS support (see :ref:siconos_install_guide)*

driver: lcp_gams()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 10000

• dparam[SICONOS_DPARAM_TOL] = 1e-12

### QP-reformulation¶

#### quadratic program formulation (SICONOS_LCP_QP)¶

Quadratic program formulation for solving a LCP

The QP we solve is

Minimize: $$z^T (M z + q)$$ subject to $$Mz + q \geq 0$$

which is the classical reformulation that can be found in Cottle, Pang and Stone (2009).

If the symmetry condition is not fulfilled, use the NSQP Solver.

driver: lcp_qp()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 0

• dparam[SICONOS_DPARAM_TOL] = 1e-6

#### quadratic program formulation (SICONOS_LCP_NSQP)¶

Non symmetric (and not nonsmooth as one could have thought in a plateform dedicated to nonsmooth problems) quadratic program formulation for solving an LCP with a non symmetric matrix.

driver: lcp_nsqp()

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 0

• dparam[SICONOS_DPARAM_TOL] = 1e-6

### AVI reformulation¶

#### AVI with Cao/Ferris solver (SICONOS_AVI_CAOFERRIS)¶

Reformulates the LCP as an Affine Variational Inequalities (AVI), then uses the solver by Cao and Ferris.

parameters:

• iparam[SICONOS_IPARAM_MAX_ITER] = 10000

• dparam[SICONOS_DPARAM_TOL] = 1e-12