lemon-project-template-glpk

comparison deps/glpk/doc/glpk04.tex @ 9:33de93886c88

Import GLPK 4.47
author Alpar Juttner <alpar@cs.elte.hu>
date Sun, 06 Nov 2011 20:59:10 +0100
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:726f2b9b30a1
1 %* glpk04.tex *%
2
3 \chapter{Advanced API Routines}
4
5 \section{Background}
6 \label{basbgd}
7
8 Using vector and matrix notations LP problem (1.1)---(1.3) (see Section
9 \ref{seclp}, page \pageref{seclp}) can be stated as follows:
10
11 \medskip
12
13 \noindent
14 \hspace{.5in} minimize (or maximize)
15 $$z=c^Tx_S+c_0\eqno(3.1)$$
16 \hspace{.5in} subject to linear constraints
17 $$x_R=Ax_S\eqno(3.2)$$
18 \hspace{.5in} and bounds of variables
19 $$
20 \begin{array}{l@{\ }c@{\ }l@{\ }c@{\ }l}
21 l_R&\leq&x_R&\leq&u_R\\
22 l_S&\leq&x_S&\leq&u_S\\
23 \end{array}\eqno(3.3)
24 $$
25 where:
26
27 \noindent
28 $x_R=(x_1,\dots,x_m)$ is the vector of auxiliary variables;
29
30 \noindent
31 $x_S=(x_{m+1},\dots,x_{m+n})$ is the vector of structural
32 variables;
33
34 \noindent
35 $z$ is the objective function;
36
37 \noindent
38 $c=(c_1,\dots,c_n)$ is the vector of objective coefficients;
39
40 \noindent
41 $c_0$ is the constant term (``shift'') of the objective function;
42
43 \noindent
44 $A=(a_{11},\dots,a_{mn})$ is the constraint matrix;
45
46 \noindent
47 $l_R=(l_1,\dots,l_m)$ is the vector of lower bounds of auxiliary
48 variables;
49
50 \noindent
51 $u_R=(u_1,\dots,u_m)$ is the vector of upper bounds of auxiliary
52 variables;
53
54 \noindent
55 $l_S=(l_{m+1},\dots,l_{m+n})$ is the vector of lower bounds of
56 structural variables;
57
58 \noindent
59 $u_S=(u_{m+1},\dots,u_{m+n})$ is the vector of upper bounds of
60 structural variables.
61
62 \medskip
63
64 From the simplex method's standpoint there is no difference between
65 auxiliary and structural variables. This allows combining all these
66 variables into one vector that leads to the following problem statement:
67
68 \medskip
69
70 \noindent
71 \hspace{.5in} minimize (or maximize)
72 $$z=(0\ |\ c)^Tx+c_0\eqno(3.4)$$
73 \hspace{.5in} subject to linear constraints
74 $$(I\ |-\!A)x=0\eqno(3.5)$$
75 \hspace{.5in} and bounds of variables
76 $$l\leq x\leq u\eqno(3.6)$$
77 where:
78
79 \noindent
80 $x=(x_R\ |\ x_S)$ is the $(m+n)$-vector of (all) variables;
81
82 \noindent
83 $(0\ |\ c)$ is the $(m+n)$-vector of objective
84 coefficients;\footnote{Subvector 0 corresponds to objective coefficients
85 at auxiliary variables.}
86
87 \noindent
88 $(I\ |-\!A)$ is the {\it augmented} constraint
89 $m\times(m+n)$-matrix;\footnote{Note that due to auxiliary variables
90 matrix $(I\ |-\!A)$ contains the unity submatrix and therefore has full
91 rank. This means, in particular, that the system (3.5) has no linearly
92 dependent constraints.}
93
94 \noindent
95 $l=(l_R\ |\ l_S)$ is the $(m+n)$-vector of lower bounds of (all)
96 variables;
97
98 \noindent
99 $u=(u_R\ |\ u_S)$ is the $(m+n)$-vector of upper bounds of (all)
100 variables.
101
102 \medskip
103
104 By definition an {\it LP basic solution} geometrically is a point in
105 the space of all variables, which is the intersection of planes
106 corresponding to active constraints\footnote{A constraint is called
107 {\it active} if in a given point it is satisfied as equality, otherwise
108 it is called {\it inactive}.}. The space of all variables has the
109 dimension $m+n$, therefore, to define some basic solution we have to
110 define $m+n$ active constraints. Note that $m$ constraints (3.5) being
111 linearly independent equalities are always active, so remaining $n$
112 active constraints can be chosen only from bound constraints (3.6).
113
114 A variable is called {\it non-basic}, if its (lower or upper) bound is
115 active, otherwise it is called {\it basic}. Since, as was said above,
116 exactly $n$ bound constraints must be active, in any basic solution
117 there are always $n$ non-basic variables and $m$ basic variables.
118 (Note that a free variable also can be non-basic. Although such
119 variable has no bounds, we can think it as the difference between two
120 non-negative variables, which both are non-basic in this case.)
121
122 Now consider how to determine numeric values of all variables for a
123 given basic solution.
124
125 Let $\Pi$ be an appropriate permutation matrix of the order $(m+n)$.
126 Then we can write:
127 $$\left(\begin{array}{@{}c@{}}x_B\\x_N\\\end{array}\right)=
128 \Pi\left(\begin{array}{@{}c@{}}x_R\\x_S\\\end{array}\right)=\Pi x,
129 \eqno(3.7)$$
130 where $x_B$ is the vector of basic variables, $x_N$ is the vector of
131 non-basic variables, $x=(x_R\ |\ x_S)$ is the vector of all variables
132 in the original order. In this case the system of linear constraints
133 (3.5) can be rewritten as follows:
134 $$(I\ |-\!A)\Pi^T\Pi x=0\ \ \ \Rightarrow\ \ \ (B\ |\ N)
135 \left(\begin{array}{@{}c@{}}x_B\\x_N\\\end{array}\right)=0,\eqno(3.8)$$
136 where
137 $$(B\ |\ N)=(I\ |-\!A)\Pi^T.\eqno(3.9)$$
138 Matrix $B$ is a square non-singular $m\times m$-matrix, which is
139 composed from columns of the augmented constraint matrix corresponding
140 to basic variables. It is called the {\it basis matrix} or simply the
141 {\it basis}. Matrix $N$ is a rectangular $m\times n$-matrix, which is
142 composed from columns of the augmented constraint matrix corresponding
143 to non-basic variables.
144
145 From (3.8) it follows that:
146 $$Bx_B+Nx_N=0,\eqno(3.10)$$
147 therefore,
148 $$x_B=-B^{-1}Nx_N.\eqno(3.11)$$
149 Thus, the formula (3.11) shows how to determine numeric values of basic
150 variables $x_B$ assuming that non-basic variables $x_N$ are fixed on
151 their active bounds.
152
153 The $m\times n$-matrix
154 $$\Xi=-B^{-1}N,\eqno(3.12)$$
155 which appears in (3.11), is called the {\it simplex
156 tableau}.\footnote{This definition corresponds to the GLPK
157 implementation.} It shows how basic variables depend on non-basic
158 variables:
159 $$x_B=\Xi x_N.\eqno(3.13)$$
160
161 The system (3.13) is equivalent to the system (3.5) in the sense that
162 they both define the same set of points in the space of (primal)
163 variables, which satisfy to these systems. If, moreover, values of all
164 basic variables satisfy to their bound constraints (3.3), the
165 corresponding basic solution is called {\it (primal) feasible},
166 otherwise {\it (primal) infeasible}. It is understood that any (primal)
167 feasible basic solution satisfy to all constraints (3.2) and (3.3).
168
169 The LP theory says that if LP has optimal solution, it has (at least
170 one) basic feasible solution, which corresponds to the optimum. And the
171 most natural way to determine whether a given basic solution is optimal
172 or not is to use the Karush---Kuhn---Tucker optimality conditions.
173
174 \def\arraystretch{1.5}
175
176 For the problem statement (3.4)---(3.6) the optimality conditions are
177 the following:\footnote{These conditions can be appiled to any solution,
178 not only to a basic solution.}
179 $$(I\ |-\!A)x=0\eqno(3.14)$$
180 $$(I\ |-\!A)^T\pi+\lambda_l+\lambda_u=\nabla z=(0\ |\ c)^T\eqno(3.15)$$
181 $$l\leq x\leq u\eqno(3.16)$$
182 $$\lambda_l\geq 0,\ \ \lambda_u\leq 0\ \ \mbox{(minimization)}
183 \eqno(3.17)$$
184 $$\lambda_l\leq 0,\ \ \lambda_u\geq 0\ \ \mbox{(maximization)}
185 \eqno(3.18)$$
186 $$(\lambda_l)_k(x_k-l_k)=0,\ \ (\lambda_u)_k(x_k-u_k)=0,\ \ k=1,2,\dots,
187 m+n\eqno(3.19)$$
188 where:
189 $\pi=(\pi_1,\pi_2,\dots,\pi_m)$ is a $m$-vector of Lagrange
190 multipliers for equality constraints (3.5);
191 $\lambda_l=[(\lambda_l)_1,(\lambda_l)_2,\dots,(\lambda_l)_n]$ is a
192 $n$-vector of Lagrange multipliers for lower bound constraints (3.6);
193 $\lambda_u=[(\lambda_u)_1,(\lambda_u)_2,\dots,(\lambda_u)_n]$ is a
194 $n$-vector of Lagrange multipliers for upper bound constraints (3.6).
195
196 Condition (3.14) is the {\it primal} (original) system of equality
197 constraints (3.5).
198
199 Condition (3.15) is the {\it dual} system of equality constraints.
200 It requires the gradient of the objective function to be a linear
201 combination of normals to the planes defined by constraints of the
202 original problem.
203
204 Condition (3.16) is the primal (original) system of bound constraints
205 (3.6).
206
207 Condition (3.17) (or (3.18) in case of maximization) is the dual system
208 of bound constraints.
209
210 Condition (3.19) is the {\it complementary slackness condition}. It
211 requires, for each original (auxiliary or structural) variable $x_k$,
212 that either its (lower or upper) bound must be active, or zero bound of
213 the corresponding Lagrange multiplier ($(\lambda_l)_k$ or
214 $(\lambda_u)_k$) must be active.
215
216 In GLPK two multipliers $(\lambda_l)_k$ and $(\lambda_u)_k$ for each
217 primal (original) variable $x_k$, $k=1,2,\dots,m+n$, are combined into
218 one multiplier:
219 $$\lambda_k=(\lambda_l)_k+(\lambda_u)_k,\eqno(3.20)$$
220 which is called a {\it dual variable} for $x_k$. This {\it cannot} lead
221 to the ambiguity, because both lower and upper bounds of $x_k$ cannot be
222 active at the same time,\footnote{If $x_k$ is a fixed variable, we can
223 think it as double-bounded variable $l_k\leq x_k\leq u_k$, where
224 $l_k=u_k.$} so at least one of $(\lambda_l)_k$ and $(\lambda_u)_k$ must
225 be equal to zero, and because these multipliers have different signs,
226 the combined multiplier, which is their sum, uniquely defines each of
227 them.
228
229 \def\arraystretch{1}
230
231 Using dual variables $\lambda_k$ the dual system of bound constraints
232 (3.17) and (3.18) can be written in the form of so called {\it ``rule of
233 signs''} as follows:
234
235 \begin{center}
236 \begin{tabular}{|@{\,}c@{$\,$}|@{$\,$}c@{$\,$}|@{$\,$}c@{$\,$}|
237 @{$\,$}c|c@{$\,$}|@{$\,$}c@{$\,$}|@{$\,$}c@{$\,$}|}
238 \hline
239 Original bound&\multicolumn{3}{c|}{Minimization}&\multicolumn{3}{c|}
240 {Maximization}\\
241 \cline{2-7}
242 constraint&$(\lambda_l)_k$&$(\lambda_u)_k$&$(\lambda_l)_k+
243 (\lambda_u)_k$&$(\lambda_l)_k$&$(\lambda_u)_k$&$(\lambda_l)_k+
244 (\lambda_u)_k$\\
245 \hline
246 $-\infty<x_k<+\infty$&$=0$&$=0$&$\lambda_k=0$&$=0$&$=0$&$\lambda_k=0$\\
247 $x_k\geq l_k$&$\geq 0$&$=0$&$\lambda_k\geq 0$&$\leq 0$&$=0$&$\lambda_k
248 \leq0$\\
249 $x_k\leq u_k$&$=0$&$\leq 0$&$\lambda_k\leq 0$&$=0$&$\geq 0$&$\lambda_k
250 \geq0$\\
251 $l_k\leq x_k\leq u_k$&$\geq 0$& $\leq 0$& $-\infty\!<\!\lambda_k\!<
252 \!+\infty$
253 &$\leq 0$& $\geq 0$& $-\infty\!<\!\lambda_k\!<\!+\infty$\\
254 $x_k=l_k=u_k$&$\geq 0$& $\leq 0$& $-\infty\!<\!\lambda_k\!<\!+\infty$&
255 $\leq 0$&
256 $\geq 0$& $-\infty\!<\!\lambda_k\!<\!+\infty$\\
257 \hline
258 \end{tabular}
259 \end{center}
260
261 May note that each primal variable $x_k$ has its dual counterpart
262 $\lambda_k$ and vice versa. This allows applying the same partition for
263 the vector of dual variables as (3.7):
264 $$\left(\begin{array}{@{}c@{}}\lambda_B\\\lambda_N\\\end{array}\right)=
265 \Pi\lambda,\eqno(3.21)$$
266 where $\lambda_B$ is a vector of dual variables for basic variables
267 $x_B$, $\lambda_N$ is a vector of dual variables for non-basic variables
268 $x_N$.
269
270 By definition, bounds of basic variables are inactive constraints, so in
271 any basic solution $\lambda_B=0$. Corresponding values of dual variables
272 $\lambda_N$ for non-basic variables $x_N$ can be determined in the
273 following way. From the dual system (3.15) we have:
274 $$(I\ |-\!A)^T\pi+\lambda=(0\ |\ c)^T,\eqno(3.22)$$
275 so multiplying both sides of (3.22) by matrix $\Pi$ gives:
276 $$\Pi(I\ |-\!A)^T\pi+\Pi\lambda=\Pi(0\ |\ c)^T.\eqno(3.23)$$
277 From (3.9) it follows that
278 $$\Pi(I\ |-\!A)^T=[(I\ |-\!A)\Pi^T]^T=(B\ |\ N)^T.\eqno(3.24)$$
279 Further, we can apply the partition (3.7) also to the vector of
280 objective coefficients (see (3.4)):
281 $$\left(\begin{array}{@{}c@{}}c_B\\c_N\\\end{array}\right)=
282 \Pi\left(\begin{array}{@{}c@{}}0\\c\\\end{array}\right),\eqno(3.25)$$
283 where $c_B$ is a vector of objective coefficients at basic variables,
284 $c_N$ is a vector of objective coefficients at non-basic variables.
285 Now, substituting (3.24), (3.21), and (3.25) into (3.23), leads to:
286 $$(B\ |\ N)^T\pi+(\lambda_B\ |\ \lambda_N)^T=(c_B\ |\ c_N)^T,
287 \eqno(3.26)$$
288 and transposing both sides of (3.26) gives the system:
289 $$\left(\begin{array}{@{}c@{}}B^T\\N^T\\\end{array}\right)\pi+
290 \left(\begin{array}{@{}c@{}}\lambda_B\\\lambda_N\\\end{array}\right)=
291 \left(\begin{array}{@{}c@{}}c_B\\c_T\\\end{array}\right),\eqno(3.27)$$
292 which can be written as follows:
293 $$\left\{
294 \begin{array}{@{\ }r@{\ }c@{\ }r@{\ }c@{\ }l@{\ }}
295 B^T\pi&+&\lambda_B&=&c_B\\
296 N^T\pi&+&\lambda_N&=&c_N\\
297 \end{array}
298 \right.\eqno(3.28)
299 $$
300 Lagrange multipliers $\pi=(\pi_i)$ correspond to equality constraints
301 (3.5) and therefore can have any sign. This allows resolving the first
302 subsystem of (3.28) as follows:\footnote{$B^{-T}$ means $(B^T)^{-1}=
303 (B^{-1})^T$.}
304 $$\pi=B^{-T}(c_B-\lambda_B)=-B^{-T}\lambda_B+B^{-T}c_B,\eqno(3.29)$$
305 and substitution of $\pi$ from (3.29) into the second subsystem of
306 (3.28) gives:
307 $$\lambda_N=-N^T\pi+c_N=N^TB^{-T}\lambda_B+(c_N-N^TB^{-T}c_B).
308 \eqno(3.30)$$
309 The latter system can be written in the following final form:
310 $$\lambda_N=-\Xi^T\lambda_B+d,\eqno(3.31)$$
311 where $\Xi$ is the simplex tableau (see (3.12)), and
312 $$d=c_N-N^TB^{-T}c_B=c_N+\Xi^Tc_B\eqno(3.32)$$
313 is the vector of so called {\it reduced costs} of non-basic variables.
314
315 \pagebreak
316
317 Above it was said that in any basic solution $\lambda_B=0$, so
318 $\lambda_N=d$ as it follows from (3.31).
319
320 The system (3.31) is equivalent to the system (3.15) in the sense that
321 they both define the same set of points in the space of dual variables
322 $\lambda$, which satisfy to these systems. If, moreover, values of all
323 dual variables $\lambda_N$ (i.e. reduced costs $d$) satisfy to their
324 bound constraints (i.e. to the ``rule of signs''; see the table above),
325 the corresponding basic solution is called {\it dual feasible},
326 otherwise {\it dual infeasible}. It is understood that any dual feasible
327 solution satisfy to all constraints (3.15) and (3.17) (or (3.18) in case
328 of maximization).
329
330 It can be easily shown that the complementary slackness condition
331 (3.19) is always satisfied for {\it any} basic solution. Therefore,
332 a basic solution\footnote{It is assumed that a complete basic solution
333 has the form $(x,\lambda)$, i.e. it includes primal as well as dual
334 variables.} is {\it optimal} if and only if it is primal and dual
335 feasible, because in this case it satifies to all the optimality
336 conditions (3.14)---(3.19).
337
338 \def\arraystretch{1.5}
339
340 The meaning of reduced costs $d=(d_j)$ of non-basic variables can be
341 explained in the following way. From (3.4), (3.7), and (3.25) it follows
342 that:
343 $$z=c_B^Tx_B+c_N^Tx_N+c_0.\eqno(3.33)$$
344 Substituting $x_B$ from (3.11) into (3.33) we can eliminate basic
345 variables and express the objective only through non-basic variables:
346 $$
347 \begin{array}{r@{\ }c@{\ }l}
348 z&=&c_B^T(-B^{-1}Nx_N)+c_N^Tx_N+c_0=\\
349 &=&(c_N^T-c_B^TB^{-1}N)x_N+c_0=\\
350 &=&(c_N-N^TB^{-T}c_B)^Tx_N+c_0=\\
351 &=&d^Tx_N+c_0.
352 \end{array}\eqno(3.34)
353 $$
354 From (3.34) it is seen that reduced cost $d_j$ shows how the objective
355 function $z$ depends on non-basic variable $(x_N)_j$ in the neighborhood
356 of the current basic solution, i.e. while the current basis remains
357 unchanged.
358
359 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
360
361 \newpage
362
363 \section{LP basis routines}
364 \label{lpbasis}
365
366 \subsection{glp\_bf\_exists---check if the basis factorization exists}
367
368 \subsubsection*{Synopsis}
369
370 \begin{verbatim}
371 int glp_bf_exists(glp_prob *lp);
372 \end{verbatim}
373
374 \subsubsection*{Returns}
375
376 If the basis factorization for the current basis associated with the
377 specified problem object exists and therefore is available for
378 computations, the routine \verb|glp_bf_exists| returns non-zero.
379 Otherwise the routine returns zero.
380
381 \subsubsection*{Comments}
382
383 Let the problem object have $m$ rows and $n$ columns. In GLPK the
384 {\it basis matrix} $B$ is a square non-singular matrix of the order $m$,
385 whose columns correspond to basic (auxiliary and/or structural)
386 variables. It is defined by the following main
387 equality:\footnote{For more details see Subsection \ref{basbgd},
388 page \pageref{basbgd}.}
389 $$(B\ |\ N)=(I\ |-\!A)\Pi^T,$$
390 where $I$ is the unity matrix of the order $m$, whose columns correspond
391 to auxiliary variables; $A$ is the original constraint
392 $m\times n$-matrix, whose columns correspond to structural variables;
393 $(I\ |-\!A)$ is the augmented constraint\linebreak
394 $m\times(m+n)$-matrix, whose columns correspond to all (auxiliary and
395 structural) variables following in the original order; $\Pi$ is a
396 permutation matrix of the order $m+n$; and $N$ is a rectangular
397 $m\times n$-matrix, whose columns correspond to non-basic (auxiliary
398 and/or structural) variables.
399
400 For various reasons it may be necessary to solve linear systems with
401 matrix $B$. To provide this possibility the GLPK implementation
402 maintains an invertable form of $B$ (that is, some representation of
403 $B^{-1}$) called the {\it basis factorization}, which is an internal
404 component of the problem object. Typically, the basis factorization is
405 computed by the simplex solver, which keeps it in the problem object
406 to be available for other computations.
407
408 Should note that any changes in the problem object, which affects the
409 basis matrix (e.g. changing the status of a row or column, changing
410 a basic column of the constraint matrix, removing an active constraint,
411 etc.), invalidates the basis factorization. So before calling any API
412 routine, which uses the basis factorization, the application program
413 must make sure (using the routine \verb|glp_bf_exists|) that the
414 factorization exists and therefore available for computations.
415
416 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
417
418 \subsection{glp\_factorize---compute the basis factorization}
419
420 \subsubsection*{Synopsis}
421
422 \begin{verbatim}
423 int glp_factorize(glp_prob *lp);
424 \end{verbatim}
425
426 \subsubsection*{Description}
427
428 The routine \verb|glp_factorize| computes the basis factorization for
429 the current basis associated with the specified problem
430 object.\footnote{The current basis is defined by the current statuses
431 of rows (auxiliary variables) and columns (structural variables).}
432
433 The basis factorization is computed from ``scratch'' even if it exists,
434 so the application program may use the routine \verb|glp_bf_exists|,
435 and, if the basis factorization already exists, not to call the routine
436 \verb|glp_factorize| to prevent an extra work.
437
438 The routine \verb|glp_factorize| {\it does not} compute components of
439 the basic solution (i.e. primal and dual values).
440
441 \subsubsection*{Returns}
442
443 \begin{tabular}{@{}p{25mm}p{97.3mm}@{}}
444 0 & The basis factorization has been successfully computed.\\
445 \verb|GLP_EBADB| & The basis matrix is invalid, because the number of
446 basic (auxiliary and structural) variables is not the same as the number
447 of rows in the problem object.\\
448 \verb|GLP_ESING| & The basis matrix is singular within the working
449 precision.\\
450 \verb|GLP_ECOND| & The basis matrix is ill-conditioned, i.e. its
451 condition number is too large.\\
452 \end{tabular}
453
454 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
455
456 \newpage
457
458 \subsection{glp\_bf\_updated---check if the basis factorization has\\
459 been updated}
460
461 \subsubsection*{Synopsis}
462
463 \begin{verbatim}
464 int glp_bf_updated(glp_prob *lp);
465 \end{verbatim}
466
467 \subsubsection*{Returns}
468
469 If the basis factorization has been just computed from ``scratch'', the
470 routine \verb|glp_bf_updated| returns zero. Otherwise, if the
471 factorization has been updated at least once, the routine returns
472 non-zero.
473
474 \subsubsection*{Comments}
475
476 {\it Updating} the basis factorization means recomputing it to reflect
477 changes in the basis matrix. For example, on every iteration of the
478 simplex method some column of the current basis matrix is replaced by a
479 new column that gives a new basis matrix corresponding to the adjacent
480 basis. In this case computing the basis factorization for the adjacent
481 basis from ``scratch'' (as the routine \verb|glp_factorize| does) would
482 be too time-consuming.
483
484 On the other hand, since the basis factorization update is a numeric
485 computational procedure, applying it many times may lead to accumulating
486 round-off errors. Therefore the basis is periodically refactorized
487 (reinverted) from ``scratch'' (with the routine \verb|glp_factorize|)
488 that allows improving its numerical properties.
489
490 The routine \verb|glp_bf_updated| allows determining if the basis
491 factorization has been updated at least once since it was computed from
492 ``scratch''.
493
494 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
495
496 \newpage
497
498 \subsection{glp\_get\_bfcp---retrieve basis factorization control
499 parameters}
500
501 \subsubsection*{Synopsis}
502
503 \begin{verbatim}
504 void glp_get_bfcp(glp_prob *lp, glp_bfcp *parm);
505 \end{verbatim}
506
507 \subsubsection*{Description}
508
509 The routine \verb|glp_get_bfcp| retrieves control parameters, which are
510 used on computing and updating the basis factorization associated with
511 the specified problem object.
512
513 Current values of the control parameters are stored in a \verb|glp_bfcp|
514 structure, which the parameter \verb|parm| points to. For a detailed
515 description of the structure \verb|glp_bfcp| see comments to the routine
516 \verb|glp_set_bfcp| in the next subsection.
517
518 \subsubsection*{Comments}
519
520 The purpose of the routine \verb|glp_get_bfcp| is two-fold. First, it
521 allows the application program obtaining current values of control
522 parameters used by internal GLPK routines, which compute and update the
523 basis factorization.
524
525 The second purpose of this routine is to provide proper values for all
526 fields of the structure \verb|glp_bfcp| in the case when the application
527 program needs to change some control parameters.
528
529 \subsection{glp\_set\_bfcp---change basis factorization control
530 parameters}
531
532 \subsubsection*{Synopsis}
533
534 \begin{verbatim}
535 void glp_set_bfcp(glp_prob *lp, const glp_bfcp *parm);
536 \end{verbatim}
537
538 \subsubsection*{Description}
539
540 The routine \verb|glp_set_bfcp| changes control parameters, which are
541 used by internal GLPK routines on computing and updating the basis
542 factorization associated with the specified problem object.
543
544 New values of the control parameters should be passed in a structure
545 \verb|glp_bfcp|, which the parameter \verb|parm| points to. For a
546 detailed description of the structure \verb|glp_bfcp| see paragraph
547 ``Control parameters'' below.
548
549 The parameter \verb|parm| can be specified as \verb|NULL|, in which case
550 all control parameters are reset to their default values.
551
552 \subsubsection*{Comments}
553
554 Before changing some control parameters with the routine
555 \verb|glp_set_bfcp| the application program should retrieve current
556 values of all control parameters with the routine \verb|glp_get_bfcp|.
557 This is needed for backward compatibility, because in the future there
558 may appear new members in the structure \verb|glp_bfcp|.
559
560 Note that new values of control parameters come into effect on a next
561 computation of the basis factorization, not immediately.
562
563 \subsubsection*{Example}
564
565 \begin{verbatim}
566 glp_prob *lp;
567 glp_bfcp parm;
568 . . .
569 /* retrieve current values of control parameters */
570 glp_get_bfcp(lp, &parm);
571 /* change the threshold pivoting tolerance */
572 parm.piv_tol = 0.05;
573 /* set new values of control parameters */
574 glp_set_bfcp(lp, &parm);
575 . . .
576 \end{verbatim}
577
578 \subsubsection*{Control parameters}
579
580 This paragraph describes all basis factorization control parameters
581 currently used in the package. Symbolic names of control parameters are
582 names of corresponding members in the structure \verb|glp_bfcp|.
583
584 \def\arraystretch{1}
585
586 \medskip
587
588 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
589 \multicolumn{2}{@{}l}{{\tt int type} (default: {\tt GLP\_BF\_FT})} \\
590 &Basis factorization type:\\
591 &\verb|GLP_BF_FT|---$LU$ + Forrest--Tomlin update;\\
592 &\verb|GLP_BF_BG|---$LU$ + Schur complement + Bartels--Golub update;\\
593 &\verb|GLP_BF_GR|---$LU$ + Schur complement + Givens rotation update.
594 \\
595 &In case of \verb|GLP_BF_FT| the update is applied to matrix $U$, while
596 in cases of \verb|GLP_BF_BG| and \verb|GLP_BF_GR| the update is applied
597 to the Schur complement.
598 \end{tabular}
599
600 \medskip
601
602 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
603 \multicolumn{2}{@{}l}{{\tt int lu\_size} (default: {\tt 0})} \\
604 &The initial size of the Sparse Vector Area, in non-zeros, used on
605 computing $LU$-factorization of the basis matrix for the first time.
606 If this parameter is set to 0, the initial SVA size is determined
607 automatically.\\
608 \end{tabular}
609
610 \medskip
611
612 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
613 \multicolumn{2}{@{}l}{{\tt double piv\_tol} (default: {\tt 0.10})} \\
614 &Threshold pivoting (Markowitz) tolerance, 0 $<$ \verb|piv_tol| $<$ 1,
615 used on computing $LU$-factorization of the basis matrix. Element
616 $u_{ij}$ of the active submatrix of factor $U$ fits to be pivot if it
617 satisfies to the stability criterion
618 $|u_{ij}| >= {\tt piv\_tol}\cdot\max|u_{i*}|$, i.e. if it is not very
619 small in the magnitude among other elements in the same row. Decreasing
620 this parameter may lead to better sparsity at the expense of numerical
621 accuracy, and vice versa.\\
622 \end{tabular}
623
624 \medskip
625
626 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
627 \multicolumn{2}{@{}l}{{\tt int piv\_lim} (default: {\tt 4})} \\
628 &This parameter is used on computing $LU$-factorization of the basis
629 matrix and specifies how many pivot candidates needs to be considered
630 on choosing a pivot element, \verb|piv_lim| $\geq$ 1. If \verb|piv_lim|
631 candidates have been considered, the pivoting routine prematurely
632 terminates the search with the best candidate found.\\
633 \end{tabular}
634
635 \medskip
636
637 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
638 \multicolumn{2}{@{}l}{{\tt int suhl} (default: {\tt GLP\_ON})} \\
639 &This parameter is used on computing $LU$-factorization of the basis
640 matrix. Being set to {\tt GLP\_ON} it enables applying the following
641 heuristic proposed by Uwe Suhl: if a column of the active submatrix has
642 no eligible pivot candidates, it is no more considered until it becomes
643 a column singleton. In many cases this allows reducing the time needed
644 for pivot searching. To disable this heuristic the parameter \verb|suhl|
645 should be set to {\tt GLP\_OFF}.\\
646 \end{tabular}
647
648 \medskip
649
650 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
651 \multicolumn{2}{@{}l}{{\tt double eps\_tol} (default: {\tt 1e-15})} \\
652 &Epsilon tolerance, \verb|eps_tol| $\geq$ 0, used on computing
653 $LU$-factorization of the basis matrix. If an element of the active
654 submatrix of factor $U$ is less than \verb|eps_tol| in the magnitude,
655 it is replaced by exact zero.\\
656 \end{tabular}
657
658 \medskip
659
660 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
661 \multicolumn{2}{@{}l}{{\tt double max\_gro} (default: {\tt 1e+10})} \\
662 &Maximal growth of elements of factor $U$, \verb|max_gro| $\geq$ 1,
663 allowable on computing $LU$-factorization of the basis matrix. If on
664 some elimination step the ratio $u_{big}/b_{max}$ (where $u_{big}$ is
665 the largest magnitude of elements of factor $U$ appeared in its active
666 submatrix during all the factorization process, $b_{max}$ is the largest
667 magnitude of elements of the basis matrix to be factorized), the basis
668 matrix is considered as ill-conditioned.\\
669 \end{tabular}
670
671 \medskip
672
673 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
674 \multicolumn{2}{@{}l}{{\tt int nfs\_max} (default: {\tt 100})} \\
675 &Maximal number of additional row-like factors (entries of the eta
676 file), \verb|nfs_max| $\geq$ 1, which can be added to $LU$-factorization
677 of the basis matrix on updating it with the Forrest--Tomlin technique.
678 This parameter is used only once, before $LU$-factorization is computed
679 for the first time, to allocate working arrays. As a rule, each update
680 adds one new factor (however, some updates may need no addition), so
681 this parameter limits the number of updates between refactorizations.\\
682 \end{tabular}
683
684 \medskip
685
686 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
687 \multicolumn{2}{@{}l}{{\tt double upd\_tol} (default: {\tt 1e-6})} \\
688 &Update tolerance, 0 $<$ \verb|upd_tol| $<$ 1, used on updating
689 $LU$-factorization of the basis matrix with the Forrest--Tomlin
690 technique. If after updating the magnitude of some diagonal element
691 $u_{kk}$ of factor $U$ becomes less than
692 ${\tt upd\_tol}\cdot\max(|u_{k*}|, |u_{*k}|)$, the factorization is
693 considered as inaccurate.\\
694 \end{tabular}
695
696 \medskip
697
698 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
699 \multicolumn{2}{@{}l}{{\tt int nrs\_max} (default: {\tt 100})} \\
700 &Maximal number of additional rows and columns, \verb|nrs_max| $\geq$ 1,
701 which can be added to $LU$-factorization of the basis matrix on updating
702 it with the Schur complement technique. This parameter is used only
703 once, before $LU$-factorization is computed for the first time, to
704 allocate working arrays. As a rule, each update adds one new row and
705 column (however, some updates may need no addition), so this parameter
706 limits the number of updates between refactorizations.\\
707 \end{tabular}
708
709 \medskip
710
711 \noindent\begin{tabular}{@{}p{17pt}@{}p{120.5mm}@{}}
712 \multicolumn{2}{@{}l}{{\tt int rs\_size} (default: {\tt 0})} \\
713 &The initial size of the Sparse Vector Area, in non-zeros, used to
714 store non-zero elements of additional rows and columns introduced on
715 updating $LU$-factorization of the basis matrix with the Schur
716 complement technique. If this parameter is set to 0, the initial SVA
717 size is determined automatically.\\
718 \end{tabular}
719
720 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
721
722 \newpage
723
724 \subsection{glp\_get\_bhead---retrieve the basis header information}
725
726 \subsubsection*{Synopsis}
727
728 \begin{verbatim}
729 int glp_get_bhead(glp_prob *lp, int k);
730 \end{verbatim}
731
732 \subsubsection*{Description}
733
734 The routine \verb|glp_get_bhead| returns the basis header information
735 for the current basis associated with the specified problem object.
736
737 \subsubsection*{Returns}
738
739 If basic variable $(x_B)_k$, $1\leq k\leq m$, is $i$-th auxiliary
740 variable ($1\leq i\leq m$), the routine returns $i$. Otherwise, if
741 $(x_B)_k$ is $j$-th structural variable ($1\leq j\leq n$), the routine
742 returns $m+j$. Here $m$ is the number of rows and $n$ is the number of
743 columns in the problem object.
744
745 \subsubsection*{Comments}
746
747 Sometimes the application program may need to know which original
748 (auxiliary and structural) variable correspond to a given basic
749 variable, or, that is the same, which column of the augmented constraint
750 matrix $(I\ |-\!A)$ correspond to a given column of the basis matrix
751 $B$.
752
753 \def\arraystretch{1}
754
755 The correspondence is defined as follows:\footnote{For more details see
756 Subsection \ref{basbgd}, page \pageref{basbgd}.}
757 $$\left(\begin{array}{@{}c@{}}x_B\\x_N\\\end{array}\right)=
758 \Pi\left(\begin{array}{@{}c@{}}x_R\\x_S\\\end{array}\right)
759 \ \ \Leftrightarrow
760 \ \ \left(\begin{array}{@{}c@{}}x_R\\x_S\\\end{array}\right)=
761 \Pi^T\left(\begin{array}{@{}c@{}}x_B\\x_N\\\end{array}\right),$$
762 where $x_B$ is the vector of basic variables, $x_N$ is the vector of
763 non-basic variables, $x_R$ is the vector of auxiliary variables
764 following in their original order,\footnote{The original order of
765 auxiliary and structural variables is defined by the ordinal numbers
766 of corresponding rows and columns in the problem object.} $x_S$ is the
767 vector of structural variables following in their original order, $\Pi$
768 is a permutation matrix (which is a component of the basis
769 factorization).
770
771 Thus, if $(x_B)_k=(x_R)_i$ is $i$-th auxiliary variable, the routine
772 returns $i$, and if $(x_B)_k=(x_S)_j$ is $j$-th structural variable,
773 the routine returns $m+j$, where $m$ is the number of rows in the
774 problem object.
775
776 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
777
778 \newpage
779
780 \subsection{glp\_get\_row\_bind---retrieve row index in the basis\\
781 header}
782
783 \subsubsection*{Synopsis}
784
785 \begin{verbatim}
786 int glp_get_row_bind(glp_prob *lp, int i);
787 \end{verbatim}
788
789 \subsubsection*{Returns}
790
791 The routine \verb|glp_get_row_bind| returns the index $k$ of basic
792 variable $(x_B)_k$, $1\leq k\leq m$, which is $i$-th auxiliary variable
793 (that is, the auxiliary variable corresponding to $i$-th row),
794 $1\leq i\leq m$, in the current basis associated with the specified
795 problem object, where $m$ is the number of rows. However, if $i$-th
796 auxiliary variable is non-basic, the routine returns zero.
797
798 \subsubsection*{Comments}
799
800 The routine \verb|glp_get_row_bind| is an inverse to the routine
801 \verb|glp_get_bhead|: if \verb|glp_get_bhead|$(lp,k)$ returns $i$,
802 \verb|glp_get_row_bind|$(lp,i)$ returns $k$, and vice versa.
803
804 \subsection{glp\_get\_col\_bind---retrieve column index in the basis
805 header}
806
807 \subsubsection*{Synopsis}
808
809 \begin{verbatim}
810 int glp_get_col_bind(glp_prob *lp, int j);
811 \end{verbatim}
812
813 \subsubsection*{Returns}
814
815 The routine \verb|glp_get_col_bind| returns the index $k$ of basic
816 variable $(x_B)_k$, $1\leq k\leq m$, which is $j$-th structural
817 variable (that is, the structural variable corresponding to $j$-th
818 column), $1\leq j\leq n$, in the current basis associated with the
819 specified problem object, where $m$ is the number of rows, $n$ is the
820 number of columns. However, if $j$-th structural variable is non-basic,
821 the routine returns zero.
822
823 \subsubsection*{Comments}
824
825 The routine \verb|glp_get_col_bind| is an inverse to the routine
826 \verb|glp_get_bhead|: if \verb|glp_get_bhead|$(lp,k)$ returns $m+j$,
827 \verb|glp_get_col_bind|$(lp,j)$ returns $k$, and vice versa.
828
829 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
830
831 \newpage
832
833 \subsection{glp\_ftran---perform forward transformation}
834
835 \subsubsection*{Synopsis}
836
837 \begin{verbatim}
838 void glp_ftran(glp_prob *lp, double x[]);
839 \end{verbatim}
840
841 \subsubsection*{Description}
842
843 The routine \verb|glp_ftran| performs forward transformation (FTRAN),
844 i.e. it solves the system $Bx=b$, where $B$ is the basis matrix
845 associated with the specified problem object, $x$ is the vector of
846 unknowns to be computed, $b$ is the vector of right-hand sides.
847
848 On entry to the routine elements of the vector $b$ should be stored in
849 locations \verb|x[1]|, \dots, \verb|x[m]|, where $m$ is the number of
850 rows. On exit the routine stores elements of the vector $x$ in the same
851 locations.
852
853 \subsection{glp\_btran---perform backward transformation}
854
855 \subsubsection*{Synopsis}
856
857 \begin{verbatim}
858 void glp_btran(glp_prob *lp, double x[]);
859 \end{verbatim}
860
861 \subsubsection*{Description}
862
863 The routine \verb|glp_btran| performs backward transformation (BTRAN),
864 i.e. it solves the system $B^Tx=b$, where $B^T$ is a matrix transposed
865 to the basis matrix $B$ associated with the specified problem object,
866 $x$ is the vector of unknowns to be computed, $b$ is the vector of
867 right-hand sides.
868
869 On entry to the routine elements of the vector $b$ should be stored in
870 locations \verb|x[1]|, \dots, \verb|x[m]|, where $m$ is the number of
871 rows. On exit the routine stores elements of the vector $x$ in the same
872 locations.
873
874 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
875
876 \newpage
877
878 \subsection{glp\_warm\_up---``warm up'' LP basis}
879
880 \subsubsection*{Synopsis}
881
882 \begin{verbatim}
883 int glp_warm_up(glp_prob *P);
884 \end{verbatim}
885
886 \subsubsection*{Description}
887
888 The routine \verb|glp_warm_up| ``warms up'' the LP basis for the
889 specified problem object using current statuses assigned to rows and
890 columns (that is, to auxiliary and structural variables).
891
892 This operation includes computing factorization of the basis matrix
893 (if it does not exist), computing primal and dual components of basic
894 solution, and determining the solution status.
895
896 \subsubsection*{Returns}
897
898 \begin{tabular}{@{}p{25mm}p{97.3mm}@{}}
899 0 & The operation has been successfully performed.\\
900 \verb|GLP_EBADB| & The basis matrix is invalid, because the number of
901 basic (auxiliary and structural) variables is not the same as the number
902 of rows in the problem object.\\
903 \verb|GLP_ESING| & The basis matrix is singular within the working
904 precision.\\
905 \verb|GLP_ECOND| & The basis matrix is ill-conditioned, i.e. its
906 condition number is too large.\\
907 \end{tabular}
908
909 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
910
911 \newpage
912
913 \section{Simplex tableau routines}
914
915 \subsection{glp\_eval\_tab\_row---compute row of the tableau}
916
917 \subsubsection*{Synopsis}
918
919 \begin{verbatim}
920 int glp_eval_tab_row(glp_prob *lp, int k, int ind[],
921 double val[]);
922 \end{verbatim}
923
924 \subsubsection*{Description}
925
926 The routine \verb|glp_eval_tab_row| computes a row of the current
927 simplex tableau (see Subsection 3.1.1, formula (3.12)), which (row)
928 corresponds to some basic variable specified by the parameter $k$ as
929 follows: if $1\leq k\leq m$, the basic variable is $k$-th auxiliary
930 variable, and if $m+1\leq k\leq m+n$, the basic variable is $(k-m)$-th
931 structural variable, where $m$ is the number of rows and $n$ is the
932 number of columns in the specified problem object. The basis
933 factorization must exist.
934
935 The computed row shows how the specified basic variable depends on
936 non-basic variables:
937 $$x_k=(x_B)_i=\xi_{i1}(x_N)_1+\xi_{i2}(x_N)_2+\dots+\xi_{in}(x_N)_n,$$
938 where $\xi_{i1}$, $\xi_{i2}$, \dots, $\xi_{in}$ are elements of the
939 simplex table row, $(x_N)_1$, $(x_N)_2$, \dots, $(x_N)_n$ are non-basic
940 (auxiliary and structural) variables.
941
942 The routine stores column indices and corresponding numeric values of
943 non-zero elements of the computed row in unordered sparse format in
944 locations \verb|ind[1]|, \dots, \verb|ind[len]| and \verb|val[1]|,
945 \dots, \verb|val[len]|, respectively, where $0\leq{\tt len}\leq n$ is
946 the number of non-zero elements in the row returned on exit.
947
948 Element indices stored in the array \verb|ind| have the same sense as
949 index $k$, i.e. indices 1 to $m$ denote auxiliary variables while
950 indices $m+1$ to $m+n$ denote structural variables (all these variables
951 are obviously non-basic by definition).
952
953 \subsubsection*{Returns}
954
955 The routine \verb|glp_eval_tab_row| returns \verb|len|, which is the
956 number of non-zero elements in the simplex table row stored in the
957 arrays \verb|ind| and \verb|val|.
958
959 \subsubsection*{Comments}
960
961 A row of the simplex table is computed as follows. At first, the
962 routine checks that the specified variable $x_k$ is basic and uses the
963 permutation matrix $\Pi$ (3.7) to determine index $i$ of basic variable
964 $(x_B)_i$, which corresponds to $x_k$.
965
966 The row to be computed is $i$-th row of the matrix $\Xi$ (3.12),
967 therefore:
968 $$\xi_i=e_i^T\Xi=-e_i^TB^{-1}N=-(B^{-T}e_i)^TN,$$
969 where $e_i$ is $i$-th unity vector. So the routine performs BTRAN to
970 obtain $i$-th row of the inverse $B^{-1}$:
971 $$\varrho_i=B^{-T}e_i,$$
972 and then computes elements of the simplex table row as inner products:
973 $$\xi_{ij}=-\varrho_i^TN_j,\ \ j=1,2,\dots,n,$$
974 where $N_j$ is $j$-th column of matrix $N$ (3.9), which (column)
975 corresponds to non-basic variable $(x_N)_j$. The permutation matrix
976 $\Pi$ is used again to convert indices $j$ of non-basic columns to
977 original ordinal numbers of auxiliary and structural variables.
978
979 \subsection{glp\_eval\_tab\_col---compute column of the tableau}
980
981 \subsubsection*{Synopsis}
982
983 \begin{verbatim}
984 int glp_eval_tab_col(glp_prob *lp, int k, int ind[],
985 double val[]);
986 \end{verbatim}
987
988 \subsubsection*{Description}
989
990 The routine \verb|glp_eval_tab_col| computes a column of the current
991 simplex tableau (see Subsection 3.1.1, formula (3.12)), which (column)
992 corresponds to some non-basic variable specified by the parameter $k$:
993 if $1\leq k\leq m$, the non-basic variable is $k$-th auxiliary variable,
994 and if $m+1\leq k\leq m+n$, the non-basic variable is $(k-m)$-th
995 structural variable, where $m$ is the number of rows and $n$ is the
996 number of columns in the specified problem object. The basis
997 factorization must exist.
998
999 The computed column shows how basic variables depends on the specified
1000 non-basic variable $x_k=(x_N)_j$:
1001 $$
1002 \begin{array}{r@{\ }c@{\ }l@{\ }l}
1003 (x_B)_1&=&\dots+\xi_{1j}(x_N)_j&+\dots\\
1004 (x_B)_2&=&\dots+\xi_{2j}(x_N)_j&+\dots\\
1005 .\ \ .&.&.\ \ .\ \ .\ \ .\ \ .\ \ .\ \ .\\
1006 (x_B)_m&=&\dots+\xi_{mj}(x_N)_j&+\dots\\
1007 \end{array}
1008 $$
1009 where $\xi_{1j}$, $\xi_{2j}$, \dots, $\xi_{mj}$ are elements of the
1010 simplex table column, $(x_B)_1$, $(x_B)_2$, \dots, $(x_B)_m$ are basic
1011 (auxiliary and structural) variables.
1012
1013 The routine stores row indices and corresponding numeric values of
1014 non-zero elements of the computed column in unordered sparse format in
1015 locations \verb|ind[1]|, \dots, \verb|ind[len]| and \verb|val[1]|,
1016 \dots, \verb|val[len]|, respectively, where $0\leq{\tt len}\leq m$ is
1017 the number of non-zero elements in the column returned on exit.
1018
1019 Element indices stored in the array \verb|ind| have the same sense as
1020 index $k$, i.e. indices 1 to $m$ denote auxiliary variables while
1021 indices $m+1$ to $m+n$ denote structural variables (all these variables
1022 are obviously basic by definition).
1023
1024 \subsubsection*{Returns}
1025
1026 The routine \verb|glp_eval_tab_col| returns \verb|len|, which is the
1027 number of non-zero elements in the simplex table column stored in the
1028 arrays \verb|ind| and \verb|val|.
1029
1030 \subsubsection*{Comments}
1031
1032 A column of the simplex table is computed as follows. At first, the
1033 routine checks that the specified variable $x_k$ is non-basic and uses
1034 the permutation matrix $\Pi$ (3.7) to determine index $j$ of non-basic
1035 variable $(x_N)_j$, which corresponds to $x_k$.
1036
1037 The column to be computed is $j$-th column of the matrix $\Xi$ (3.12),
1038 therefore:
1039 $$\Xi_j=\Xi e_j=-B^{-1}Ne_j=-B^{-1}N_j,$$
1040 where $e_j$ is $j$-th unity vector, $N_j$ is $j$-th column of matrix
1041 $N$ (3.9). So the routine performs FTRAN to transform $N_j$ to the
1042 simplex table column $\Xi_j=(\xi_{ij})$ and uses the permutation matrix
1043 $\Pi$ to convert row indices $i$ to original ordinal numbers of
1044 auxiliary and structural variables.
1045
1046 \newpage
1047
1048 \subsection{glp\_transform\_row---transform explicitly specified\\
1049 row}
1050
1051 \subsubsection*{Synopsis}
1052
1053 \begin{verbatim}
1054 int glp_transform_row(glp_prob *P, int len, int ind[],
1055 double val[]);
1056 \end{verbatim}
1057
1058 \subsubsection*{Description}
1059
1060 The routine \verb|glp_transform_row| performs the same operation as the
1061 routine \verb|glp_eval_tab_row| with exception that the row to be
1062 transformed is specified explicitly as a sparse vector.
1063
1064 The explicitly specified row may be thought as a linear form:
1065 $$x=a_1x_{m+1}+a_2x_{m+2}+\dots+a_nx_{m+n},$$
1066 where $x$ is an auxiliary variable for this row, $a_j$ are coefficients
1067 of the linear form, $x_{m+j}$ are structural variables.
1068
1069 On entry column indices and numerical values of non-zero coefficients
1070 $a_j$ of the specified row should be placed in locations \verb|ind[1]|,
1071 \dots, \verb|ind[len]| and \verb|val[1]|, \dots, \verb|val[len]|, where
1072 \verb|len| is number of non-zero coefficients.
1073
1074 This routine uses the system of equality constraints and the current
1075 basis in order to express the auxiliary variable $x$ through the current
1076 non-basic variables (as if the transformed row were added to the problem
1077 object and the auxiliary variable $x$ were basic), i.e. the resultant
1078 row has the form:
1079 $$x=\xi_1(x_N)_1+\xi_2(x_N)_2+\dots+\xi_n(x_N)_n,$$
1080 where $\xi_j$ are influence coefficients, $(x_N)_j$ are non-basic
1081 (auxiliary and structural) variables, $n$ is the number of columns in
1082 the problem object.
1083
1084 On exit the routine stores indices and numerical values of non-zero
1085 coefficients $\xi_j$ of the resultant row in locations \verb|ind[1]|,
1086 \dots, \verb|ind[len']| and \verb|val[1]|, \dots, \verb|val[len']|,
1087 where $0\leq{\tt len'}\leq n$ is the number of non-zero coefficients in
1088 the resultant row returned by the routine. Note that indices of
1089 non-basic variables stored in the array \verb|ind| correspond to
1090 original ordinal numbers of variables: indices 1 to $m$ mean auxiliary
1091 variables and indices $m+1$ to $m+n$ mean structural ones.
1092
1093 \subsubsection*{Returns}
1094
1095 The routine \verb|glp_transform_row| returns \verb|len'|, the number of
1096 non-zero coefficients in the resultant row stored in the arrays
1097 \verb|ind| and \verb|val|.
1098
1099 \subsection{glp\_transform\_col---transform explicitly specified\\
1100 column}
1101
1102 \subsubsection*{Synopsis}
1103
1104 \begin{verbatim}
1105 int glp_transform_col(glp_prob *P, int len, int ind[],
1106 double val[]);
1107 \end{verbatim}
1108
1109 \subsubsection*{Description}
1110
1111 The routine \verb|glp_transform_col| performs the same operation as the
1112 routine \verb|glp_eval_tab_col| with exception that the column to be
1113 transformed is specified explicitly as a sparse vector.
1114
1115 The explicitly specified column may be thought as it were added to
1116 the original system of equality constraints:
1117 $$
1118 \begin{array}{l@{\ }c@{\ }r@{\ }c@{\ }r@{\ }c@{\ }r}
1119 x_1&=&a_{11}x_{m+1}&+\dots+&a_{1n}x_{m+n}&+&a_1x \\
1120 x_2&=&a_{21}x_{m+1}&+\dots+&a_{2n}x_{m+n}&+&a_2x \\
1121 \multicolumn{7}{c}
1122 {.\ \ .\ \ .\ \ .\ \ .\ \ .\ \ .\ \ .\ \ .\ \ .\ \ .\ \ .\ \ .\ \ .}\\
1123 x_m&=&a_{m1}x_{m+1}&+\dots+&a_{mn}x_{m+n}&+&a_mx \\
1124 \end{array}
1125 $$
1126 where $x_i$ are auxiliary variables, $x_{m+j}$ are structural variables
1127 (presented in the problem object), $x$ is a structural variable for the
1128 explicitly specified column, $a_i$ are constraint coefficients at $x$.
1129
1130 On entry row indices and numerical values of non-zero coefficients
1131 $a_i$ of the specified column should be placed in locations
1132 \verb|ind[1]|, \dots, \verb|ind[len]| and \verb|val[1]|, \dots,
1133 \verb|val[len]|, where \verb|len| is number of non-zero coefficients.
1134
1135 This routine uses the system of equality constraints and the current
1136 basis in order to express the current basic variables through the
1137 structural variable $x$ (as if the transformed column were added to the
1138 problem object and the variable $x$ were non-basic):
1139 $$
1140 \begin{array}{l@{\ }c@{\ }r}
1141 (x_B)_1&=\dots+&\xi_{1}x\\
1142 (x_B)_2&=\dots+&\xi_{2}x\\
1143 \multicolumn{3}{c}{.\ \ .\ \ .\ \ .\ \ .\ \ .}\\
1144 (x_B)_m&=\dots+&\xi_{m}x\\
1145 \end{array}
1146 $$
1147 where $\xi_i$ are influence coefficients, $x_B$ are basic (auxiliary
1148 and structural) variables, $m$ is the number of rows in the problem
1149 object.
1150
1151 On exit the routine stores indices and numerical values of non-zero
1152 coefficients $\xi_i$ of the resultant column in locations \verb|ind[1]|,
1153 \dots, \verb|ind[len']| and \verb|val[1]|, \dots, \verb|val[len']|,
1154 where $0\leq{\tt len'}\leq m$ is the number of non-zero coefficients in
1155 the resultant column returned by the routine. Note that indices of basic
1156 variables stored in the array \verb|ind| correspond to original ordinal
1157 numbers of variables, i.e. indices 1 to $m$ mean auxiliary variables,
1158 indices $m+1$ to $m+n$ mean structural ones.
1159
1160 \subsubsection*{Returns}
1161
1162 The routine \verb|glp_transform_col| returns \verb|len'|, the number of
1163 non-zero coefficients in the resultant column stored in the arrays
1164 \verb|ind| and \verb|val|.
1165
1166 \subsection{glp\_prim\_rtest---perform primal ratio test}
1167
1168 \subsubsection*{Synopsis}
1169
1170 \begin{verbatim}
1171 int glp_prim_rtest(glp_prob *P, int len, const int ind[],
1172 const double val[], int dir, double eps);
1173 \end{verbatim}
1174
1175 \subsubsection*{Description}
1176
1177 The routine \verb|glp_prim_rtest| performs the primal ratio test using
1178 an explicitly specified column of the simplex table.
1179
1180 The current basic solution associated with the LP problem object must be
1181 primal feasible.
1182
1183 The explicitly specified column of the simplex table shows how the basic
1184 variables $x_B$ depend on some non-basic variable $x$ (which is not
1185 necessarily presented in the problem object):
1186 $$
1187 \begin{array}{l@{\ }c@{\ }r}
1188 (x_B)_1&=\dots+&\xi_{1}x\\
1189 (x_B)_2&=\dots+&\xi_{2}x\\
1190 \multicolumn{3}{c}{.\ \ .\ \ .\ \ .\ \ .\ \ .}\\
1191 (x_B)_m&=\dots+&\xi_{m}x\\
1192 \end{array}
1193 $$
1194
1195 The column is specifed on entry to the routine in sparse format. Ordinal
1196 numbers of basic variables $(x_B)_i$ should be placed in locations
1197 \verb|ind[1]|, \dots, \verb|ind[len]|, where ordinal number 1 to $m$
1198 denote auxiliary variables, and ordinal numbers $m+1$ to $m+n$ denote
1199 structural variables. The corresponding non-zero coefficients $\xi_i$
1200 should be placed in locations \verb|val[1]|, \dots, \verb|val[len]|. The
1201 arrays \verb|ind| and \verb|val| are not changed by the routine.
1202
1203 The parameter \verb|dir| specifies direction in which the variable $x$
1204 changes on entering the basis: $+1$ means increasing, $-1$ means
1205 decreasing.
1206
1207 The parameter \verb|eps| is an absolute tolerance (small positive
1208 number, say, $10^{-9}$) used by the routine to skip $\xi_i$'s whose
1209 magnitude is less than \verb|eps|.
1210
1211 The routine determines which basic variable (among those specified in
1212 \verb|ind[1]|, \dots, \verb|ind[len]|) reaches its (lower or upper)
1213 bound first before any other basic variables do, and which, therefore,
1214 should leave the basis in order to keep primal feasibility.
1215
1216 \subsubsection*{Returns}
1217
1218 The routine \verb|glp_prim_rtest| returns the index, \verb|piv|, in the
1219 arrays \verb|ind| and \verb|val| corresponding to the pivot element
1220 chosen, $1\leq$ \verb|piv| $\leq$ \verb|len|. If the adjacent basic
1221 solution is primal unbounded, and therefore the choice cannot be made,
1222 the routine returns zero.
1223
1224 \subsubsection*{Comments}
1225
1226 If the non-basic variable $x$ is presented in the LP problem object, the
1227 input column can be computed with the routine \verb|glp_eval_tab_col|;
1228 otherwise, it can be computed with the routine \verb|glp_transform_col|.
1229
1230 \subsection{glp\_dual\_rtest---perform dual ratio test}
1231
1232 \subsubsection*{Synopsis}
1233
1234 \begin{verbatim}
1235 int glp_dual_rtest(glp_prob *P, int len, const int ind[],
1236 const double val[], int dir, double eps);
1237 \end{verbatim}
1238
1239 \subsubsection*{Description}
1240
1241 The routine \verb|glp_dual_rtest| performs the dual ratio test using
1242 an explicitly specified row of the simplex table.
1243
1244 The current basic solution associated with the LP problem object must be
1245 dual feasible.
1246
1247 The explicitly specified row of the simplex table is a linear form
1248 that shows how some basic variable $x$ (which is not necessarily
1249 presented in the problem object) depends on non-basic variables $x_N$:
1250 $$x=\xi_1(x_N)_1+\xi_2(x_N)_2+\dots+\xi_n(x_N)_n.$$
1251
1252 The row is specified on entry to the routine in sparse format. Ordinal
1253 numbers of non-basic variables $(x_N)_j$ should be placed in locations
1254 \verb|ind[1]|, \dots, \verb|ind[len]|, where ordinal numbers 1 to $m$
1255 denote auxiliary variables, and ordinal numbers $m+1$ to $m+n$ denote
1256 structural variables. The corresponding non-zero coefficients $\xi_j$
1257 should be placed in locations \verb|val[1]|, \dots, \verb|val[len]|.
1258 The arrays \verb|ind| and \verb|val| are not changed by the routine.
1259
1260 The parameter \verb|dir| specifies direction in which the variable $x$
1261 changes on leaving the basis: $+1$ means that $x$ goes on its lower
1262 bound, so its reduced cost (dual variable) is increasing (minimization)
1263 or decreasing (maximization); $-1$ means that $x$ goes on its upper
1264 bound, so its reduced cost is decreasing (minimization) or increasing
1265 (maximization).
1266
1267 The parameter \verb|eps| is an absolute tolerance (small positive
1268 number, say, $10^{-9}$) used by the routine to skip $\xi_j$'s whose
1269 magnitude is less than \verb|eps|.
1270
1271 The routine determines which non-basic variable (among those specified
1272 in \verb|ind[1]|, \dots, \verb|ind[len]|) should enter the basis in
1273 order to keep dual feasibility, because its reduced cost reaches the
1274 (zero) bound first before this occurs for any other non-basic variables.
1275
1276 \subsubsection*{Returns}
1277
1278 The routine \verb|glp_dual_rtest| returns the index, \verb|piv|, in the
1279 arrays \verb|ind| and \verb|val| corresponding to the pivot element
1280 chosen, $1\leq$ \verb|piv| $\leq$ \verb|len|. If the adjacent basic
1281 solution is dual unbounded, and therefore the choice cannot be made,
1282 the routine returns zero.
1283
1284 \subsubsection*{Comments}
1285
1286 If the basic variable $x$ is presented in the LP problem object, the
1287 input row can be computed with the routine \verb|glp_eval_tab_row|;
1288 otherwise, it can be computed with the routine \verb|glp_transform_row|.
1289
1290 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1291
1292 \newpage
1293
1294 \section{Post-optimal analysis routines}
1295
1296 \subsection{glp\_analyze\_bound---analyze active bound of non-basic
1297 variable}
1298
1299 \subsubsection*{Synopsis}
1300
1301 \begin{verbatim}
1302 void glp_analyze_bound(glp_prob *P, int k, double *limit1,
1303 int *var1, double *limit2, int *var2);
1304 \end{verbatim}
1305
1306 \subsubsection*{Description}
1307
1308 The routine \verb|glp_analyze_bound| analyzes the effect of varying the
1309 active bound of specified non-basic variable.
1310
1311 The non-basic variable is specified by the parameter $k$, where
1312 $1\leq k\leq m$ means auxiliary variable of corresponding row, and
1313 $m+1\leq k\leq m+n$ means structural variable (column).
1314
1315 Note that the current basic solution must be optimal, and the basis
1316 factorization must exist.
1317
1318 Results of the analysis have the following meaning.
1319
1320 \verb|value1| is the minimal value of the active bound, at which the
1321 basis still remains primal feasible and thus optimal. \verb|-DBL_MAX|
1322 means that the active bound has no lower limit.
1323
1324 \verb|var1| is the ordinal number of an auxiliary (1 to $m$) or
1325 structural ($m+1$ to $m+n$) basic variable, which reaches its bound
1326 first and thereby limits further decreasing the active bound being
1327 analyzed. if \verb|value1| = \verb|-DBL_MAX|, \verb|var1| is set to 0.
1328
1329 \verb|value2| is the maximal value of the active bound, at which the
1330 basis still remains primal feasible and thus optimal. \verb|+DBL_MAX|
1331 means that the active bound has no upper limit.
1332
1333 \verb|var2| is the ordinal number of an auxiliary (1 to $m$) or
1334 structural ($m+1$ to $m+n$) basic variable, which reaches its bound
1335 first and thereby limits further increasing the active bound being
1336 analyzed. if \verb|value2| = \verb|+DBL_MAX|, \verb|var2| is set to 0.
1337
1338 The parameters \verb|value1|, \verb|var1|, \verb|value2|, \verb|var2|
1339 can be specified as \verb|NULL|, in which case corresponding information
1340 is not stored.
1341
1342 \newpage
1343
1344 \subsection{glp\_analyze\_coef---analyze objective coefficient at basic
1345 variable}
1346
1347 \subsubsection*{Synopsis}
1348
1349 \begin{verbatim}
1350 void glp_analyze_coef(glp_prob *P, int k, double *coef1,
1351 int *var1, double *value1, double *coef2, int *var2,
1352 double *value2);
1353 \end{verbatim}
1354
1355 \subsubsection*{Description}
1356
1357 The routine \verb|glp_analyze_coef| analyzes the effect of varying the
1358 objective coefficient at specified basic variable.
1359
1360 The basic variable is specified by the parameter $k$, where
1361 $1\leq k\leq m$ means auxiliary variable of corresponding row, and
1362 $m+1\leq k\leq m+n$ means structural variable (column).
1363
1364 Note that the current basic solution must be optimal, and the basis
1365 factorization must exist.
1366
1367 Results of the analysis have the following meaning.
1368
1369 \verb|coef1| is the minimal value of the objective coefficient, at
1370 which the basis still remains dual feasible and thus optimal.
1371 \verb|-DBL_MAX| means that the objective coefficient has no lower limit.
1372
1373 \verb|var1| is the ordinal number of an auxiliary (1 to $m$) or
1374 structural ($m+1$ to $m+n$) non-basic variable, whose reduced cost
1375 reaches its zero bound first and thereby limits further decreasing the
1376 objective coefficient being analyzed. If \verb|coef1| = \verb|-DBL_MAX|,
1377 \verb|var1| is set to 0.
1378
1379 \verb|value1| is value of the basic variable being analyzed in an
1380 adjacent basis, which is defined as follows. Let the objective
1381 coefficient reaches its minimal value (\verb|coef1|) and continues
1382 decreasing. Then the reduced cost of the limiting non-basic variable
1383 (\verb|var1|) becomes dual infeasible and the current basis becomes
1384 non-optimal that forces the limiting non-basic variable to enter the
1385 basis replacing there some basic variable that leaves the basis to keep
1386 primal feasibility. Should note that on determining the adjacent basis
1387 current bounds of the basic variable being analyzed are ignored as if
1388 it were free (unbounded) variable, so it cannot leave the basis. It may
1389 happen that no dual feasible adjacent basis exists, in which case
1390 \verb|value1| is set to \verb|-DBL_MAX| or \verb|+DBL_MAX|.
1391
1392 \verb|coef2| is the maximal value of the objective coefficient, at
1393 which the basis still remains dual feasible and thus optimal.
1394 \verb|+DBL_MAX| means that the objective coefficient has no upper limit.
1395
1396 \verb|var2| is the ordinal number of an auxiliary (1 to $m$) or
1397 structural ($m+1$ to $m+n$) non-basic variable, whose reduced cost
1398 reaches its zero bound first and thereby limits further increasing the
1399 objective coefficient being analyzed. If \verb|coef2| = \verb|+DBL_MAX|,
1400 \verb|var2| is set to 0.
1401
1402 \verb|value2| is value of the basic variable being analyzed in an
1403 adjacent basis, which is defined exactly in the same way as
1404 \verb|value1| above with exception that now the objective coefficient
1405 is increasing.
1406
1407 The parameters \verb|coef1|, \verb|var1|, \verb|value1|, \verb|coef2|,
1408 \verb|var2|, \verb|value2| can be specified as \verb|NULL|, in which
1409 case corresponding information is not stored.
1410
1411 %* eof *%