DyLP trunk
|
00001 /* glpinv.h */ 00002 00003 /*---------------------------------------------------------------------- 00004 -- Copyright (C) 2000, 2001, 2002, 2003 Andrew Makhorin, Department 00005 -- for Applied Informatics, Moscow Aviation Institute, Moscow, Russia. 00006 -- All rights reserved. E-mail: <mao@mai2.rcnet.ru>. 00007 -- 00008 -- This file is a part of GLPK (GNU Linear Programming Kit). 00009 -- 00010 -- Licensed under the Eclipse Public License (EPL) by permission of the 00011 -- author for inclusion in the DyLP LP distribution. 00012 ----------------------------------------------------------------------*/ 00013 00014 /* 00015 @(#)glpinv.h 1.3 06/22/04 00016 svn/cvs: $Id$ 00017 */ 00018 00019 00020 #ifndef _GLPINV_H 00021 #define _GLPINV_H 00022 00023 #include "glpluf.h" 00024 00025 #define inv_create dy_glp_inv_create 00026 #define inv_decomp dy_glp_inv_decomp 00027 #define inv_h_solve dy_glp_inv_h_solve 00028 #define inv_ftran dy_glp_inv_ftran 00029 #define inv_btran dy_glp_inv_btran 00030 #define inv_update dy_glp_inv_update 00031 #define inv_delete dy_glp_inv_delete 00032 00033 /*---------------------------------------------------------------------- 00034 -- The structure INV defines an invertable form of the basis matrix B, 00035 -- which is based on LU-factorization and is the following sextet: 00036 -- 00037 -- [B] = (F, H, V, P0, P, Q), (1) 00038 -- 00039 -- where F, H, and V are such matrices that 00040 -- 00041 -- B = F * H * V, (2) 00042 -- 00043 -- and P0, P, and Q are such permutation matrices that the matrix 00044 -- 00045 -- L = P0 * F * inv(P0) (3) 00046 -- 00047 -- is lower triangular with unity diagonal, and the matrix 00048 -- 00049 -- U = P * V * Q (4) 00050 -- 00051 -- is upper triangular. All the matrices have the same order m, which 00052 -- is the order of the basis matrix B. 00053 -- 00054 -- The matrices F, V, P, and Q are stored in the structure LUF (see the 00055 -- section GLPLUF), which is a member of the structure INV. 00056 -- 00057 -- The matrix H is stored in the form of eta file using row-like format 00058 -- as follows: 00059 -- 00060 -- H = H[1] * H[2] * ... * H[nfs], (5) 00061 -- 00062 -- where H[k], k = 1, 2, ..., nfs, is a row-like factor, which differs 00063 -- from the unity matrix only by one row, nfs is current number of row- 00064 -- like factors. After the factorization has been built for some given 00065 -- basis matrix B the matrix H has no factors and thus it is the unity 00066 -- matrix. Then each time when the factorization is recomputed for an 00067 -- adjacent basis matrix, the next factor H[k], k = 1, 2, ... is built 00068 -- and added to the end of the eta file H. 00069 -- 00070 -- Being sparse vectors non-trivial rows of the factors H[k] are stored 00071 -- in the right part of the sparse vector area (SVA) in the same manner 00072 -- as rows and columns of the matrix F. 00073 -- 00074 -- For more details see the program documentation. */ 00075 00076 typedef struct INV INV; 00077 00078 struct INV 00079 { /* invertable (factorized) form of the basis matrix */ 00080 int m; 00081 /* order of the matrices B, F, H, V, P0, P, Q */ 00082 int valid; 00083 /* if this flag is not set, the invertable form is invalid and 00084 can't be updated nor used in ftran and btran operations */ 00085 LUF *luf; 00086 /* LU-factorization (holds the matrices F, V, P, Q) */ 00087 /*--------------------------------------------------------------*/ 00088 /* matrix H in the form of eta file */ 00089 int hh_max; 00090 /* maximal number of row-like factors (that limits maximal number 00091 of updates of the factorization) */ 00092 int hh_nfs; 00093 /* current number of row-like factors (0 <= hh_nfs <= hh_max) */ 00094 int *hh_ndx; /* int hh_ndx[1+hh_max]; */ 00095 /* hh_ndx[0] is not used; 00096 hh_ndx[k], k = 1, ..., nfs, is number of a non-trivial row of 00097 the factor H[k] */ 00098 int *hh_ptr; /* int hh_ptr[1+hh_max]; */ 00099 /* hh_ptr[0] is not used; 00100 hh_ptr[k], k = 1, ..., nfs, is a pointer to the first element 00101 of the non-trivial row of the factor H[k] in the sparse vector 00102 area */ 00103 int *hh_len; /* int hh_len[1+hh_max]; */ 00104 /* hh_len[0] is not used; 00105 hh_len[k], k = 1, ..., nfs, is total number of elements in the 00106 non-trivial row of the factor H[k] */ 00107 /*--------------------------------------------------------------*/ 00108 /* matrix P0 */ 00109 int *p0_row; /* int p0_row[1+n]; */ 00110 /* p0_row[0] is not used; p0_row[i] = j means that p0[i,j] = 1 */ 00111 int *p0_col; /* int p0_col[1+n]; */ 00112 /* p0_col[0] is not used; p0_col[j] = i means that p0[i,j] = 1 */ 00113 /* if i-th row or column of the matrix F corresponds to i'-th row 00114 or column of the matrix L = P0*F*inv(P0), then p0_row[i'] = i 00115 and p0_col[i] = i' */ 00116 /*--------------------------------------------------------------*/ 00117 /* partially transformed column is inv(F*H)*B[j], where B[j] is a 00118 new column, which will replace the existing j-th column of the 00119 basis matrix */ 00120 int cc_len; 00121 /* number of (non-zero) elements in the partially transformed 00122 column; if cc_len < 0, the column has been not prepared yet */ 00123 int *cc_ndx; /* int cc_ndx[1+m]; */ 00124 /* cc_ndx[0] is not used; 00125 cc_ndx[k], k = 1, ..., cc_len, is a row index of the partially 00126 transformed column element */ 00127 double *cc_val; /* double cc_val[1+m]; */ 00128 /* cc_val[0] is not used; 00129 cc_val[k], k = 1, ..., cc_len, is a numerical (non-zero) value 00130 of the column element */ 00131 /*--------------------------------------------------------------*/ 00132 /* control parameters */ 00133 double upd_tol; 00134 #if 0 00135 /* update tolerance; if on updating the factorization absolute 00136 value of some diagonal element of the matrix U = P*V*Q is less 00137 than upd_tol, the factorization is considered as inaccurate */ 00138 #else 00139 /* update tolerance; if after the factorization has been updated 00140 absolute value of some diagonal element u[k,k] of the matrix 00141 U = P*V*Q is less than upd_tol * max(|u[k,*]|, |u[*,k]|), the 00142 factorization is considered as inaccurate */ 00143 #endif 00144 /*--------------------------------------------------------------*/ 00145 /* some statistics */ 00146 int nnz_h; 00147 /* current number of non-zeros in all factors of the matrix H */ 00148 double min_vrratio ; 00149 /* minimum value of u[k,k]/(upd_tol)*max(|u[k,*]|, |u[*,k]|) since 00150 last pivot. */ 00151 }; 00152 00153 INV *inv_create(int m, int max_upd); 00154 /* create factorization of the basis matrix */ 00155 00156 int inv_decomp(INV *inv, 00157 void *info, int (*col)(void *info, int j, int rn[], double bj[])); 00158 /* compute factorization of the basis matrix */ 00159 00160 void inv_h_solve(INV *inv, int tr, double x[]); 00161 /* solve system H*x = b or H'*x = b */ 00162 00163 void inv_ftran(INV *inv, double x[], int save); 00164 /* perform forward transformation (FTRAN) */ 00165 00166 void inv_btran(INV *inv, double x[]); 00167 /* perform backward transformation (BTRAN) */ 00168 00169 int inv_update(INV *inv, int j); 00170 /* update factorization for adjacent basis matrix */ 00171 00172 void inv_delete(INV *inv); 00173 /* delete factorization of the basis matrix */ 00174 00175 #endif 00176 00177 /* eof */