#include <stdio.h>
#include <string.h>
#include "include/lss_lib.h"
#include "include/ldp.h"
#include "include/lsi.h"

Include dependency graph for lsi.c:

Functions
int	_invert_ur_matrix (double R, int m, double Ri)
int	lsi (int m_e, int n_e, int m_g, int m_c, double E, double f, double G, double h, double C, double d, double *x)

Function Documentation

int _invert_ur_matrix	(	double *	R,
		int	m,
		double *	Ri
	)

Inverts upper triangular matrix A of m X n sized, and places the A^{-1} to Ai.

This is done by using Gauss elimination for each column of the inverse matrix, to be exact this solves R * Ri_k = I_k, for k = 1...m, and I_k is the k'th column of the Identity matrix.

Returns:: 0 - on success

Definition at line 406 of file lsi.c.

Referenced by lsi().

int lsi	(	int	m_e,
		int	n_e,
		int	m_g,
		int	m_c,
		double *	E,
		double *	f,
		double *	G,
		double *	h,
		double *	C,
		double *	d,
		double *	x
	)

Solves 'Least squares with inequality constrain' problem, aka Minimize ||Ex - f|| subject to Gx >_ h,.

Assumes rank(E) = n_e <= m_e

What the algorithm basically does is that it converts E = Q [R; 0 ]

Then we get || Ex - f || = || QEx - Qf || = || [R ; 0]y - Qf || = || R y - f_1 || + || f_2 ||

We make variable change to z = Rx - f_1 ( => x = R^-1( z + f_1 ) and minimize ||z|| to

Gx >_ h <=> G R^-1 ( z + f_1 ) >_ h <=> GR^-1 z >_ h - GR^-1 f_1

So we just calculate matrix GR^-1 and solve z with LDP, and then we solve x from z: => x = R^-1( z + f_1 );

Also this algorithm can deal with set of equality constrains: assume we have problem LSI: minismize ||Ex - f|| with Gx >_ h and also with Cx = d, let size(C) = (m_c,n_e).

assumes rank(C) = m_c < n , and rank( [C^T; E^T]) = n_e

First we calculate orthogonal matrix K, such that C K = [ C_1 0 ], and size(C_1) = (m_c, m_c) and rank(C_1) = m_c. C_1 is upper triangular.

Represent change of variables, x = K [y_1 ; y_2] = K y

mark E K = [E_1 E_2] G K = [ G_1 G_2 ]

Then we solve C x = d <=> C K y = d <=> [C_1 0 ] [y_1; y_2] = d Which leaves y_2 as non-bind variables. Then this is orignal LSI problem (since now y_1 are bind and used as constats):

Minimize || E_2 y_2 - (f - E_1 y_1) || subject to G_2 y_2 >_ h - G_1 y_1

After that we can calculate the original solution from the y = [y_1 y_2], x = K y

Some words about memory usage:

This algorithm uses Memory Handler for workspace memory allocation.

Returns:: 0 on success 1 on no solution 2 for dimension or memory error

Parameters:

m_e	Matrix E dimension
n_e	Matrix E dimension
m_g	Matrix G dimension
m_c	Matrix C dimension (ignored if C is not defined)
E	Matrix E is a m_e x n_e matrix, defining the E in \|\|Ex - f\|\|. This matrix will be modified on successfull computation
f	Vector f is a m_e x 1 matrix, defining the f in \|\|Ex - f\|\|. This vector will be modified on successfull computation
G	Matrix G is a m_g x (n_e + 1) matrix, defining the G (that is m_g x n_e matrix) in Gx >_ h. Note that the + 1 is because we want to use this G as work space. So G is handled as m_g x n_e matrix, but there must be m_g*(n_e+1) space allocated for it. This matrix will be modified
h	Vector h is a m_g x 1 matrix, defining the h in the Gx >_ h. This vector will be modified.
C	matrix C is a m_c x n matrix representing equality constrains. After successfull computation it will contain lower triangular matrix CK. If NULL is given then no equality constrains are considered. Will be modified if != NULL.
d	Vector of size m_c x 1 defineng the equality constrain constant. Must be != NULL if C != NULL. Will be modified if != NULL
x	Vector x is the solution vector, n_e x 1

Definition at line 108 of file lsi.c.

References _invert_ur_matrix(), _lss_h12(), _lss_matrix_times_vector(), _lss_print_matrix(), ABS, FREE_WORKSPACE, GET_WORKSPACE, ldp(), LSS_TEST, and SMALLEST_NON_ZERO.

Referenced by do_lsi().

Here is the call graph for this function:

Functions

Function Documentation