1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/src/glpscf.h Mon Dec 06 13:09:21 2010 +0100
1.3 @@ -0,0 +1,126 @@
1.4 +/* glpscf.h (Schur complement factorization) */
1.5 +
1.6 +/***********************************************************************
1.7 +* This code is part of GLPK (GNU Linear Programming Kit).
1.8 +*
1.9 +* Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008,
1.10 +* 2009, 2010 Andrew Makhorin, Department for Applied Informatics,
1.11 +* Moscow Aviation Institute, Moscow, Russia. All rights reserved.
1.12 +* E-mail: <mao@gnu.org>.
1.13 +*
1.14 +* GLPK is free software: you can redistribute it and/or modify it
1.15 +* under the terms of the GNU General Public License as published by
1.16 +* the Free Software Foundation, either version 3 of the License, or
1.17 +* (at your option) any later version.
1.18 +*
1.19 +* GLPK is distributed in the hope that it will be useful, but WITHOUT
1.20 +* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1.21 +* or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
1.22 +* License for more details.
1.23 +*
1.24 +* You should have received a copy of the GNU General Public License
1.25 +* along with GLPK. If not, see <http://www.gnu.org/licenses/>.
1.26 +***********************************************************************/
1.27 +
1.28 +#ifndef GLPSCF_H
1.29 +#define GLPSCF_H
1.30 +
1.31 +/***********************************************************************
1.32 +* The structure SCF defines the following factorization of a square
1.33 +* nxn matrix C (which is the Schur complement):
1.34 +*
1.35 +* F * C = U * P,
1.36 +*
1.37 +* where F is a square transforming matrix, U is an upper triangular
1.38 +* matrix, P is a permutation matrix.
1.39 +*
1.40 +* It is assumed that matrix C is small and dense, so matrices F and U
1.41 +* are stored in the dense format by rows as follows:
1.42 +*
1.43 +* 1 n n_max 1 n n_max
1.44 +* 1 * * * * * * x x x x 1 * * * * * * x x x x
1.45 +* * * * * * * x x x x . * * * * * x x x x
1.46 +* * * * * * * x x x x . . * * * * x x x x
1.47 +* * * * * * * x x x x . . . * * * x x x x
1.48 +* * * * * * * x x x x . . . . * * x x x x
1.49 +* n * * * * * * x x x x n . . . . . * x x x x
1.50 +* x x x x x x x x x x . . . . . . x x x x
1.51 +* x x x x x x x x x x . . . . . . . x x x
1.52 +* x x x x x x x x x x . . . . . . . . x x
1.53 +* n_max x x x x x x x x x x n_max . . . . . . . . . x
1.54 +*
1.55 +* matrix F matrix U
1.56 +*
1.57 +* where '*' are matrix elements, 'x' are reserved locations.
1.58 +*
1.59 +* Permutation matrix P is stored in row-like format.
1.60 +*
1.61 +* Matrix C normally is not stored.
1.62 +*
1.63 +* REFERENCES
1.64 +*
1.65 +* 1. M.A.Saunders, "LUSOL: A basis package for constrained optimiza-
1.66 +* tion," SCCM, Stanford University, 2006.
1.67 +*
1.68 +* 2. M.A.Saunders, "Notes 5: Basis Updates," CME 318, Stanford Univer-
1.69 +* sity, Spring 2006.
1.70 +*
1.71 +* 3. M.A.Saunders, "Notes 6: LUSOL---a Basis Factorization Package,"
1.72 +* ibid. */
1.73 +
1.74 +typedef struct SCF SCF;
1.75 +
1.76 +struct SCF
1.77 +{ /* Schur complement factorization */
1.78 + int n_max;
1.79 + /* maximal order of matrices C, F, U, P; n_max >= 1 */
1.80 + int n;
1.81 + /* current order of matrices C, F, U, P; n >= 0 */
1.82 + double *f; /* double f[1+n_max*n_max]; */
1.83 + /* matrix F stored by rows */
1.84 + double *u; /* double u[1+n_max*(n_max+1)/2]; */
1.85 + /* upper triangle of matrix U stored by rows */
1.86 + int *p; /* int p[1+n_max]; */
1.87 + /* matrix P; p[i] = j means that P[i,j] = 1 */
1.88 + int t_opt;
1.89 + /* type of transformation used to restore triangular structure of
1.90 + matrix U: */
1.91 +#define SCF_TBG 1 /* Bartels-Golub elimination */
1.92 +#define SCF_TGR 2 /* Givens plane rotation */
1.93 + int rank;
1.94 + /* estimated rank of matrices C and U */
1.95 + double *c; /* double c[1+n_max*n_max]; */
1.96 + /* matrix C stored in the same format as matrix F and used only
1.97 + for debugging; normally this array is not allocated */
1.98 + double *w; /* double w[1+n_max]; */
1.99 + /* working array */
1.100 +};
1.101 +
1.102 +/* return codes: */
1.103 +#define SCF_ESING 1 /* singular matrix */
1.104 +#define SCF_ELIMIT 2 /* update limit reached */
1.105 +
1.106 +#define scf_create_it _glp_scf_create_it
1.107 +SCF *scf_create_it(int n_max);
1.108 +/* create Schur complement factorization */
1.109 +
1.110 +#define scf_update_exp _glp_scf_update_exp
1.111 +int scf_update_exp(SCF *scf, const double x[], const double y[],
1.112 + double z);
1.113 +/* update factorization on expanding C */
1.114 +
1.115 +#define scf_solve_it _glp_scf_solve_it
1.116 +void scf_solve_it(SCF *scf, int tr, double x[]);
1.117 +/* solve either system C * x = b or C' * x = b */
1.118 +
1.119 +#define scf_reset_it _glp_scf_reset_it
1.120 +void scf_reset_it(SCF *scf);
1.121 +/* reset factorization for empty matrix C */
1.122 +
1.123 +#define scf_delete_it _glp_scf_delete_it
1.124 +void scf_delete_it(SCF *scf);
1.125 +/* delete Schur complement factorization */
1.126 +
1.127 +#endif
1.128 +
1.129 +/* eof */