[1] | 1 | /* glpscf.h (Schur complement factorization) */ |
---|
| 2 | |
---|
| 3 | /*********************************************************************** |
---|
| 4 | * This code is part of GLPK (GNU Linear Programming Kit). |
---|
| 5 | * |
---|
| 6 | * Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, |
---|
| 7 | * 2009, 2010 Andrew Makhorin, Department for Applied Informatics, |
---|
| 8 | * Moscow Aviation Institute, Moscow, Russia. All rights reserved. |
---|
| 9 | * E-mail: <mao@gnu.org>. |
---|
| 10 | * |
---|
| 11 | * GLPK is free software: you can redistribute it and/or modify it |
---|
| 12 | * under the terms of the GNU General Public License as published by |
---|
| 13 | * the Free Software Foundation, either version 3 of the License, or |
---|
| 14 | * (at your option) any later version. |
---|
| 15 | * |
---|
| 16 | * GLPK is distributed in the hope that it will be useful, but WITHOUT |
---|
| 17 | * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY |
---|
| 18 | * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public |
---|
| 19 | * License for more details. |
---|
| 20 | * |
---|
| 21 | * You should have received a copy of the GNU General Public License |
---|
| 22 | * along with GLPK. If not, see <http://www.gnu.org/licenses/>. |
---|
| 23 | ***********************************************************************/ |
---|
| 24 | |
---|
| 25 | #ifndef GLPSCF_H |
---|
| 26 | #define GLPSCF_H |
---|
| 27 | |
---|
| 28 | /*********************************************************************** |
---|
| 29 | * The structure SCF defines the following factorization of a square |
---|
| 30 | * nxn matrix C (which is the Schur complement): |
---|
| 31 | * |
---|
| 32 | * F * C = U * P, |
---|
| 33 | * |
---|
| 34 | * where F is a square transforming matrix, U is an upper triangular |
---|
| 35 | * matrix, P is a permutation matrix. |
---|
| 36 | * |
---|
| 37 | * It is assumed that matrix C is small and dense, so matrices F and U |
---|
| 38 | * are stored in the dense format by rows as follows: |
---|
| 39 | * |
---|
| 40 | * 1 n n_max 1 n n_max |
---|
| 41 | * 1 * * * * * * x x x x 1 * * * * * * x x x x |
---|
| 42 | * * * * * * * x x x x . * * * * * x x x x |
---|
| 43 | * * * * * * * x x x x . . * * * * x x x x |
---|
| 44 | * * * * * * * x x x x . . . * * * x x x x |
---|
| 45 | * * * * * * * x x x x . . . . * * x x x x |
---|
| 46 | * n * * * * * * x x x x n . . . . . * x x x x |
---|
| 47 | * x x x x x x x x x x . . . . . . x x x x |
---|
| 48 | * x x x x x x x x x x . . . . . . . x x x |
---|
| 49 | * x x x x x x x x x x . . . . . . . . x x |
---|
| 50 | * n_max x x x x x x x x x x n_max . . . . . . . . . x |
---|
| 51 | * |
---|
| 52 | * matrix F matrix U |
---|
| 53 | * |
---|
| 54 | * where '*' are matrix elements, 'x' are reserved locations. |
---|
| 55 | * |
---|
| 56 | * Permutation matrix P is stored in row-like format. |
---|
| 57 | * |
---|
| 58 | * Matrix C normally is not stored. |
---|
| 59 | * |
---|
| 60 | * REFERENCES |
---|
| 61 | * |
---|
| 62 | * 1. M.A.Saunders, "LUSOL: A basis package for constrained optimiza- |
---|
| 63 | * tion," SCCM, Stanford University, 2006. |
---|
| 64 | * |
---|
| 65 | * 2. M.A.Saunders, "Notes 5: Basis Updates," CME 318, Stanford Univer- |
---|
| 66 | * sity, Spring 2006. |
---|
| 67 | * |
---|
| 68 | * 3. M.A.Saunders, "Notes 6: LUSOL---a Basis Factorization Package," |
---|
| 69 | * ibid. */ |
---|
| 70 | |
---|
| 71 | typedef struct SCF SCF; |
---|
| 72 | |
---|
| 73 | struct SCF |
---|
| 74 | { /* Schur complement factorization */ |
---|
| 75 | int n_max; |
---|
| 76 | /* maximal order of matrices C, F, U, P; n_max >= 1 */ |
---|
| 77 | int n; |
---|
| 78 | /* current order of matrices C, F, U, P; n >= 0 */ |
---|
| 79 | double *f; /* double f[1+n_max*n_max]; */ |
---|
| 80 | /* matrix F stored by rows */ |
---|
| 81 | double *u; /* double u[1+n_max*(n_max+1)/2]; */ |
---|
| 82 | /* upper triangle of matrix U stored by rows */ |
---|
| 83 | int *p; /* int p[1+n_max]; */ |
---|
| 84 | /* matrix P; p[i] = j means that P[i,j] = 1 */ |
---|
| 85 | int t_opt; |
---|
| 86 | /* type of transformation used to restore triangular structure of |
---|
| 87 | matrix U: */ |
---|
| 88 | #define SCF_TBG 1 /* Bartels-Golub elimination */ |
---|
| 89 | #define SCF_TGR 2 /* Givens plane rotation */ |
---|
| 90 | int rank; |
---|
| 91 | /* estimated rank of matrices C and U */ |
---|
| 92 | double *c; /* double c[1+n_max*n_max]; */ |
---|
| 93 | /* matrix C stored in the same format as matrix F and used only |
---|
| 94 | for debugging; normally this array is not allocated */ |
---|
| 95 | double *w; /* double w[1+n_max]; */ |
---|
| 96 | /* working array */ |
---|
| 97 | }; |
---|
| 98 | |
---|
| 99 | /* return codes: */ |
---|
| 100 | #define SCF_ESING 1 /* singular matrix */ |
---|
| 101 | #define SCF_ELIMIT 2 /* update limit reached */ |
---|
| 102 | |
---|
| 103 | #define scf_create_it _glp_scf_create_it |
---|
| 104 | SCF *scf_create_it(int n_max); |
---|
| 105 | /* create Schur complement factorization */ |
---|
| 106 | |
---|
| 107 | #define scf_update_exp _glp_scf_update_exp |
---|
| 108 | int scf_update_exp(SCF *scf, const double x[], const double y[], |
---|
| 109 | double z); |
---|
| 110 | /* update factorization on expanding C */ |
---|
| 111 | |
---|
| 112 | #define scf_solve_it _glp_scf_solve_it |
---|
| 113 | void scf_solve_it(SCF *scf, int tr, double x[]); |
---|
| 114 | /* solve either system C * x = b or C' * x = b */ |
---|
| 115 | |
---|
| 116 | #define scf_reset_it _glp_scf_reset_it |
---|
| 117 | void scf_reset_it(SCF *scf); |
---|
| 118 | /* reset factorization for empty matrix C */ |
---|
| 119 | |
---|
| 120 | #define scf_delete_it _glp_scf_delete_it |
---|
| 121 | void scf_delete_it(SCF *scf); |
---|
| 122 | /* delete Schur complement factorization */ |
---|
| 123 | |
---|
| 124 | #endif |
---|
| 125 | |
---|
| 126 | /* eof */ |
---|