Return to Question

added 23 characters in body

Source Link

edited Oct 6, 2011 at 20:38

I'm trying to implement the "Dantzig; van de Panne and Whinston" principal pivoting algorithm for solving symmetric positive semi-definite LCPs from "The Linear Complementarity Problem" book (Algorithm 4.2.11 on page 248 of my edition).

The algorithm looks like this:

Given $(q,M)$ with $M$ being square, symmetric and positive semi-definite and of size $n$, try to find $z$ s.t. $Mz+q = w$ and $z\perp w, w \ge 0 \space z \ge 0$, or indicate if no such $z$ exists.

...

Initialize:

$\alpha = \emptyset$, $\beta = \{ 1, ..., n\}$

Step 1: Breaking ties arbitrarily, find: $r \in arg \space min \{ q_i : i \in \beta \}$

If $q_r \ge 0$, stop. A solution of $(q,M)$ is given by $z_\alpha = q_\alpha, z_\beta = 0$

If $M_{rr} = 0$ and $M_{ir} \ge 0$ for all $i \in \alpha$ stop. There is no solution

Step 2: Use the minimum ratio test to find the index $s \in \alpha \cup \{r\}$ of a blocking variable to $r$. Break ties arbitrarily, but prefer $r$ if there's a choice.

Step 3: Using a simple principal pivot on row/column s to update $M$ and $q$

If $s = r$, transfer $r$ from $\beta$ to $\alpha$ and goto Step 1

If $s \ne r$, transfer $s$ from $\alpha$ to $\beta$ and goto Step 2

...

This algorithm works most of the time, but in some cases where I perform $n$ pivots, I'm getting to step 1 and $\beta$ is empty, and I'm not sure how to proceed. Some detective work shows that often one more pivot will arrive at the correct answer, but I have a number of different cases and they all require slight differences once they run out of elements in $\beta$ in order to arrive at the correct answer.

One such case is: $M = \begin{bmatrix} 60 & 26 & 37 \\\\ 26 & 36 & 29 \\\\ 37 & 29 & 30 \end{bmatrix}$ and $q = \begin{bmatrix} -1 \\\\ -1 \\\\ -1 \end{bmatrix}$

For which I know the answer is: $z = \begin{bmatrix} 0 \\\\ .004184 \\\\ .029289 \end{bmatrix}$

The matrix M is symmetric PSD as near as I can tell (all the eigen values are positive), so it should be something that Dantzig's algorithm can process. Stepping through the code the above problem gets pivoted 3 times, and I find that a 4th pivot will arrive at a correct answer.

I'm getting tired of hack and check so I'm wondering if there's something I'm just missing about the algorithm's termination condition, or if this case implies a bug in my implementation.

I can provide definitions for "minimum ratio test", "arg min", and "simple principal pivot" if they aren't familiar terms (though they seem to be pretty universal to all pivoting methods for LCPs).

The algorithm looks like this:

Given $(q,M)$ with $M$ being square, symmetric and positive semi-definite and of size $n$, try to find $z$ s.t. $Mz+q = w$ and $z\perp w, w \ge 0 \space z \ge 0$, or indicate if no such $z$ exists.

...

Initialize:

$\alpha = \emptyset$, $\beta = \{ 1, ..., n\}$

Step 1: Breaking ties arbitrarily, find: $r \in arg \space min \{ q_i : i \in \beta \}$

If $q_r \ge 0$, stop. A solution of $(q,M)$ is given by $z_\alpha = q_\alpha, z_\beta = 0$

If $M_{rr} = 0$ and $M_{ir} \ge 0$ for all $i \in \alpha$ stop. There is no solution

Step 2: Use the minimum ratio test to find the index $s \in \alpha \cup \{r\}$ of a blocking variable to $r$. Break ties arbitrarily, but prefer $r$ if there's a choice.

Step 3: Using a simple principal pivot on row/column s to update $M$ and $q$

If $s = r$, transfer $r$ from $\beta$ to $\alpha$ and goto Step 1

If $s \ne r$, transfer $s$ from $\alpha$ to $\beta$ and goto Step 2

...

This algorithm works most of the time, but in some cases I'm getting to step 1 and $\beta$ is empty, and I'm not sure how to proceed. Some detective work shows that often one more pivot will arrive at the correct answer, but I have a number of different cases and they all require slight differences once they run out of elements in $\beta$ in order to arrive at the correct answer.

One such case is: $M = \begin{bmatrix} 60 & 26 & 37 \\\\ 26 & 36 & 29 \\\\ 37 & 29 & 30 \end{bmatrix}$ and $q = \begin{bmatrix} -1 \\\\ -1 \\\\ -1 \end{bmatrix}$

For which I know the answer is: $z = \begin{bmatrix} 0 \\\\ .004184 \\\\ .029289 \end{bmatrix}$

I'm getting tired of hack and check so I'm wondering if there's something I'm just missing about the algorithm's termination condition, or if this case implies a bug in my implementation.

I can provide definitions for "minimum ratio test", "arg min", and "simple principal pivot" if they aren't familiar terms (though they seem to be pretty universal to all pivoting methods for LCPs).

The algorithm looks like this:

Given $(q,M)$ with $M$ being square, symmetric and positive semi-definite and of size $n$, try to find $z$ s.t. $Mz+q = w$ and $z\perp w, w \ge 0 \space z \ge 0$, or indicate if no such $z$ exists.

...

Initialize:

$\alpha = \emptyset$, $\beta = \{ 1, ..., n\}$

Step 1: Breaking ties arbitrarily, find: $r \in arg \space min \{ q_i : i \in \beta \}$

If $q_r \ge 0$, stop. A solution of $(q,M)$ is given by $z_\alpha = q_\alpha, z_\beta = 0$

If $M_{rr} = 0$ and $M_{ir} \ge 0$ for all $i \in \alpha$ stop. There is no solution

Step 2: Use the minimum ratio test to find the index $s \in \alpha \cup \{r\}$ of a blocking variable to $r$. Break ties arbitrarily, but prefer $r$ if there's a choice.

Step 3: Using a simple principal pivot on row/column s to update $M$ and $q$

If $s = r$, transfer $r$ from $\beta$ to $\alpha$ and goto Step 1

If $s \ne r$, transfer $s$ from $\alpha$ to $\beta$ and goto Step 2

...

This algorithm works most of the time, but in cases where I perform $n$ pivots, I'm getting to step 1 and $\beta$ is empty, and I'm not sure how to proceed. Some detective work shows that often one more pivot will arrive at the correct answer, but I have a number of different cases and they all require slight differences once they run out of elements in $\beta$ in order to arrive at the correct answer.

One such case is: $M = \begin{bmatrix} 60 & 26 & 37 \\\\ 26 & 36 & 29 \\\\ 37 & 29 & 30 \end{bmatrix}$ and $q = \begin{bmatrix} -1 \\\\ -1 \\\\ -1 \end{bmatrix}$

For which I know the answer is: $z = \begin{bmatrix} 0 \\\\ .004184 \\\\ .029289 \end{bmatrix}$

I'm getting tired of hack and check so I'm wondering if there's something I'm just missing about the algorithm's termination condition, or if this case implies a bug in my implementation.

I can provide definitions for "minimum ratio test", "arg min", and "simple principal pivot" if they aren't familiar terms (though they seem to be pretty universal to all pivoting methods for LCPs).

Source Link

asked Oct 6, 2011 at 19:43

Jay Lemmon

Linear complementarity problem: principal pivoting algorithm

The algorithm looks like this:

Given $(q,M)$ with $M$ being square, symmetric and positive semi-definite and of size $n$, try to find $z$ s.t. $Mz+q = w$ and $z\perp w, w \ge 0 \space z \ge 0$, or indicate if no such $z$ exists.

...

Initialize:

$\alpha = \emptyset$, $\beta = \{ 1, ..., n\}$

Step 1: Breaking ties arbitrarily, find: $r \in arg \space min \{ q_i : i \in \beta \}$

If $q_r \ge 0$, stop. A solution of $(q,M)$ is given by $z_\alpha = q_\alpha, z_\beta = 0$

If $M_{rr} = 0$ and $M_{ir} \ge 0$ for all $i \in \alpha$ stop. There is no solution

Step 2: Use the minimum ratio test to find the index $s \in \alpha \cup \{r\}$ of a blocking variable to $r$. Break ties arbitrarily, but prefer $r$ if there's a choice.

Step 3: Using a simple principal pivot on row/column s to update $M$ and $q$

If $s = r$, transfer $r$ from $\beta$ to $\alpha$ and goto Step 1

If $s \ne r$, transfer $s$ from $\alpha$ to $\beta$ and goto Step 2

...

One such case is: $M = \begin{bmatrix} 60 & 26 & 37 \\\\ 26 & 36 & 29 \\\\ 37 & 29 & 30 \end{bmatrix}$ and $q = \begin{bmatrix} -1 \\\\ -1 \\\\ -1 \end{bmatrix}$

For which I know the answer is: $z = \begin{bmatrix} 0 \\\\ .004184 \\\\ .029289 \end{bmatrix}$

I'm getting tired of hack and check so I'm wondering if there's something I'm just missing about the algorithm's termination condition, or if this case implies a bug in my implementation.

I can provide definitions for "minimum ratio test", "arg min", and "simple principal pivot" if they aren't familiar terms (though they seem to be pretty universal to all pivoting methods for LCPs).

linear-algebra linear-programming