alpar@9: /* glplpf.h (LP basis factorization, Schur complement version) */ alpar@9: alpar@9: /*********************************************************************** alpar@9: * This code is part of GLPK (GNU Linear Programming Kit). alpar@9: * alpar@9: * Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, alpar@9: * 2009, 2010, 2011 Andrew Makhorin, Department for Applied Informatics, alpar@9: * Moscow Aviation Institute, Moscow, Russia. All rights reserved. alpar@9: * E-mail: . alpar@9: * alpar@9: * GLPK is free software: you can redistribute it and/or modify it alpar@9: * under the terms of the GNU General Public License as published by alpar@9: * the Free Software Foundation, either version 3 of the License, or alpar@9: * (at your option) any later version. alpar@9: * alpar@9: * GLPK is distributed in the hope that it will be useful, but WITHOUT alpar@9: * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY alpar@9: * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public alpar@9: * License for more details. alpar@9: * alpar@9: * You should have received a copy of the GNU General Public License alpar@9: * along with GLPK. If not, see . alpar@9: ***********************************************************************/ alpar@9: alpar@9: #ifndef GLPLPF_H alpar@9: #define GLPLPF_H alpar@9: alpar@9: #include "glpscf.h" alpar@9: #include "glpluf.h" alpar@9: alpar@9: /*********************************************************************** alpar@9: * The structure LPF defines the factorization of the basis mxm matrix alpar@9: * B, where m is the number of rows in corresponding problem instance. alpar@9: * alpar@9: * This factorization is the following septet: alpar@9: * alpar@9: * [B] = (L0, U0, R, S, C, P, Q), (1) alpar@9: * alpar@9: * and is based on the following main equality: alpar@9: * alpar@9: * ( B F^) ( B0 F ) ( L0 0 ) ( U0 R ) alpar@9: * ( ) = P ( ) Q = P ( ) ( ) Q, (2) alpar@9: * ( G^ H^) ( G H ) ( S I ) ( 0 C ) alpar@9: * alpar@9: * where: alpar@9: * alpar@9: * B is the current basis matrix (not stored); alpar@9: * alpar@9: * F^, G^, H^ are some additional matrices (not stored); alpar@9: * alpar@9: * B0 is some initial basis matrix (not stored); alpar@9: * alpar@9: * F, G, H are some additional matrices (not stored); alpar@9: * alpar@9: * P, Q are permutation matrices (stored in both row- and column-like alpar@9: * formats); alpar@9: * alpar@9: * L0, U0 are some matrices that defines a factorization of the initial alpar@9: * basis matrix B0 = L0 * U0 (stored in an invertable form); alpar@9: * alpar@9: * R is a matrix defined from L0 * R = F, so R = inv(L0) * F (stored in alpar@9: * a column-wise sparse format); alpar@9: * alpar@9: * S is a matrix defined from S * U0 = G, so S = G * inv(U0) (stored in alpar@9: * a row-wise sparse format); alpar@9: * alpar@9: * C is the Schur complement for matrix (B0 F G H). It is defined from alpar@9: * S * R + C = H, so C = H - S * R = H - G * inv(U0) * inv(L0) * F = alpar@9: * = H - G * inv(B0) * F. Matrix C is stored in an invertable form. alpar@9: * alpar@9: * REFERENCES alpar@9: * alpar@9: * 1. M.A.Saunders, "LUSOL: A basis package for constrained optimiza- alpar@9: * tion," SCCM, Stanford University, 2006. alpar@9: * alpar@9: * 2. M.A.Saunders, "Notes 5: Basis Updates," CME 318, Stanford Univer- alpar@9: * sity, Spring 2006. alpar@9: * alpar@9: * 3. M.A.Saunders, "Notes 6: LUSOL---a Basis Factorization Package," alpar@9: * ibid. */ alpar@9: alpar@9: typedef struct LPF LPF; alpar@9: alpar@9: struct LPF alpar@9: { /* LP basis factorization */ alpar@9: int valid; alpar@9: /* the factorization is valid only if this flag is set */ alpar@9: /*--------------------------------------------------------------*/ alpar@9: /* initial basis matrix B0 */ alpar@9: int m0_max; alpar@9: /* maximal value of m0 (increased automatically, if necessary) */ alpar@9: int m0; alpar@9: /* the order of B0 */ alpar@9: LUF *luf; alpar@9: /* LU-factorization of B0 */ alpar@9: /*--------------------------------------------------------------*/ alpar@9: /* current basis matrix B */ alpar@9: int m; alpar@9: /* the order of B */ alpar@9: double *B; /* double B[1+m*m]; */ alpar@9: /* B in dense format stored by rows and used only for debugging; alpar@9: normally this array is not allocated */ alpar@9: /*--------------------------------------------------------------*/ alpar@9: /* augmented matrix (B0 F G H) of the order m0+n */ alpar@9: int n_max; alpar@9: /* maximal number of additional rows and columns */ alpar@9: int n; alpar@9: /* current number of additional rows and columns */ alpar@9: /*--------------------------------------------------------------*/ alpar@9: /* m0xn matrix R in column-wise format */ alpar@9: int *R_ptr; /* int R_ptr[1+n_max]; */ alpar@9: /* R_ptr[j], 1 <= j <= n, is a pointer to j-th column */ alpar@9: int *R_len; /* int R_len[1+n_max]; */ alpar@9: /* R_len[j], 1 <= j <= n, is the length of j-th column */ alpar@9: /*--------------------------------------------------------------*/ alpar@9: /* nxm0 matrix S in row-wise format */ alpar@9: int *S_ptr; /* int S_ptr[1+n_max]; */ alpar@9: /* S_ptr[i], 1 <= i <= n, is a pointer to i-th row */ alpar@9: int *S_len; /* int S_len[1+n_max]; */ alpar@9: /* S_len[i], 1 <= i <= n, is the length of i-th row */ alpar@9: /*--------------------------------------------------------------*/ alpar@9: /* Schur complement C of the order n */ alpar@9: SCF *scf; /* SCF scf[1:n_max]; */ alpar@9: /* factorization of the Schur complement */ alpar@9: /*--------------------------------------------------------------*/ alpar@9: /* matrix P of the order m0+n */ alpar@9: int *P_row; /* int P_row[1+m0_max+n_max]; */ alpar@9: /* P_row[i] = j means that P[i,j] = 1 */ alpar@9: int *P_col; /* int P_col[1+m0_max+n_max]; */ alpar@9: /* P_col[j] = i means that P[i,j] = 1 */ alpar@9: /*--------------------------------------------------------------*/ alpar@9: /* matrix Q of the order m0+n */ alpar@9: int *Q_row; /* int Q_row[1+m0_max+n_max]; */ alpar@9: /* Q_row[i] = j means that Q[i,j] = 1 */ alpar@9: int *Q_col; /* int Q_col[1+m0_max+n_max]; */ alpar@9: /* Q_col[j] = i means that Q[i,j] = 1 */ alpar@9: /*--------------------------------------------------------------*/ alpar@9: /* Sparse Vector Area (SVA) is a set of locations intended to alpar@9: store sparse vectors which represent columns of matrix R and alpar@9: rows of matrix S; each location is a doublet (ind, val), where alpar@9: ind is an index, val is a numerical value of a sparse vector alpar@9: element; in the whole each sparse vector is a set of adjacent alpar@9: locations defined by a pointer to its first element and its alpar@9: length, i.e. the number of its elements */ alpar@9: int v_size; alpar@9: /* the SVA size, in locations; locations are numbered by integers alpar@9: 1, 2, ..., v_size, and location 0 is not used */ alpar@9: int v_ptr; alpar@9: /* pointer to the first available location */ alpar@9: int *v_ind; /* int v_ind[1+v_size]; */ alpar@9: /* v_ind[k], 1 <= k <= v_size, is the index field of location k */ alpar@9: double *v_val; /* double v_val[1+v_size]; */ alpar@9: /* v_val[k], 1 <= k <= v_size, is the value field of location k */ alpar@9: /*--------------------------------------------------------------*/ alpar@9: double *work1; /* double work1[1+m0+n_max]; */ alpar@9: /* working array */ alpar@9: double *work2; /* double work2[1+m0+n_max]; */ alpar@9: /* working array */ alpar@9: }; alpar@9: alpar@9: /* return codes: */ alpar@9: #define LPF_ESING 1 /* singular matrix */ alpar@9: #define LPF_ECOND 2 /* ill-conditioned matrix */ alpar@9: #define LPF_ELIMIT 3 /* update limit reached */ alpar@9: alpar@9: #define lpf_create_it _glp_lpf_create_it alpar@9: LPF *lpf_create_it(void); alpar@9: /* create LP basis factorization */ alpar@9: alpar@9: #define lpf_factorize _glp_lpf_factorize alpar@9: int lpf_factorize(LPF *lpf, int m, const int bh[], int (*col) alpar@9: (void *info, int j, int ind[], double val[]), void *info); alpar@9: /* compute LP basis factorization */ alpar@9: alpar@9: #define lpf_ftran _glp_lpf_ftran alpar@9: void lpf_ftran(LPF *lpf, double x[]); alpar@9: /* perform forward transformation (solve system B*x = b) */ alpar@9: alpar@9: #define lpf_btran _glp_lpf_btran alpar@9: void lpf_btran(LPF *lpf, double x[]); alpar@9: /* perform backward transformation (solve system B'*x = b) */ alpar@9: alpar@9: #define lpf_update_it _glp_lpf_update_it alpar@9: int lpf_update_it(LPF *lpf, int j, int bh, int len, const int ind[], alpar@9: const double val[]); alpar@9: /* update LP basis factorization */ alpar@9: alpar@9: #define lpf_delete_it _glp_lpf_delete_it alpar@9: void lpf_delete_it(LPF *lpf); alpar@9: /* delete LP basis factorization */ alpar@9: alpar@9: #endif alpar@9: alpar@9: /* eof */