COIN-OR::LEMON - Graph Library

source: lemon-project-template-glpk/deps/glpk/src/glpnpp02.c @ 10:5545663ca997

subpack-glpk
Last change on this file since 10:5545663ca997 was 9:33de93886c88, checked in by Alpar Juttner <alpar@…>, 13 years ago

Import GLPK 4.47

File size: 43.2 KB
Line 
1/* glpnpp02.c */
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, 2011 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#include "glpnpp.h"
26
27/***********************************************************************
28*  NAME
29*
30*  npp_free_row - process free (unbounded) row
31*
32*  SYNOPSIS
33*
34*  #include "glpnpp.h"
35*  void npp_free_row(NPP *npp, NPPROW *p);
36*
37*  DESCRIPTION
38*
39*  The routine npp_free_row processes row p, which is free (i.e. has
40*  no finite bounds):
41*
42*     -inf < sum a[p,j] x[j] < +inf.                                 (1)
43*             j
44*
45*  PROBLEM TRANSFORMATION
46*
47*  Constraint (1) cannot be active, so it is redundant and can be
48*  removed from the original problem.
49*
50*  Removing row p leads to removing a column of multiplier pi[p] for
51*  this row in the dual system. Since row p has no bounds, pi[p] = 0,
52*  so removing the column does not affect the dual solution.
53*
54*  RECOVERING BASIC SOLUTION
55*
56*  In solution to the original problem row p is inactive constraint,
57*  so it is assigned status GLP_BS, and multiplier pi[p] is assigned
58*  zero value.
59*
60*  RECOVERING INTERIOR-POINT SOLUTION
61*
62*  In solution to the original problem row p is inactive constraint,
63*  so its multiplier pi[p] is assigned zero value.
64*
65*  RECOVERING MIP SOLUTION
66*
67*  None needed. */
68
69struct free_row
70{     /* free (unbounded) row */
71      int p;
72      /* row reference number */
73};
74
75static int rcv_free_row(NPP *npp, void *info);
76
77void npp_free_row(NPP *npp, NPPROW *p)
78{     /* process free (unbounded) row */
79      struct free_row *info;
80      /* the row must be free */
81      xassert(p->lb == -DBL_MAX && p->ub == +DBL_MAX);
82      /* create transformation stack entry */
83      info = npp_push_tse(npp,
84         rcv_free_row, sizeof(struct free_row));
85      info->p = p->i;
86      /* remove the row from the problem */
87      npp_del_row(npp, p);
88      return;
89}
90
91static int rcv_free_row(NPP *npp, void *_info)
92{     /* recover free (unbounded) row */
93      struct free_row *info = _info;
94      if (npp->sol == GLP_SOL)
95         npp->r_stat[info->p] = GLP_BS;
96      if (npp->sol != GLP_MIP)
97         npp->r_pi[info->p] = 0.0;
98      return 0;
99}
100
101/***********************************************************************
102*  NAME
103*
104*  npp_geq_row - process row of 'not less than' type
105*
106*  SYNOPSIS
107*
108*  #include "glpnpp.h"
109*  void npp_geq_row(NPP *npp, NPPROW *p);
110*
111*  DESCRIPTION
112*
113*  The routine npp_geq_row processes row p, which is 'not less than'
114*  inequality constraint:
115*
116*     L[p] <= sum a[p,j] x[j] (<= U[p]),                             (1)
117*              j
118*
119*  where L[p] < U[p], and upper bound may not exist (U[p] = +oo).
120*
121*  PROBLEM TRANSFORMATION
122*
123*  Constraint (1) can be replaced by equality constraint:
124*
125*     sum a[p,j] x[j] - s = L[p],                                    (2)
126*      j
127*
128*  where
129*
130*     0 <= s (<= U[p] - L[p])                                        (3)
131*
132*  is a non-negative surplus variable.
133*
134*  Since in the primal system there appears column s having the only
135*  non-zero coefficient in row p, in the dual system there appears a
136*  new row:
137*
138*     (-1) pi[p] + lambda = 0,                                       (4)
139*
140*  where (-1) is coefficient of column s in row p, pi[p] is multiplier
141*  of row p, lambda is multiplier of column q, 0 is coefficient of
142*  column s in the objective row.
143*
144*  RECOVERING BASIC SOLUTION
145*
146*  Status of row p in solution to the original problem is determined
147*  by its status and status of column q in solution to the transformed
148*  problem as follows:
149*
150*     +--------------------------------------+------------------+
151*     |         Transformed problem          | Original problem |
152*     +-----------------+--------------------+------------------+
153*     | Status of row p | Status of column s | Status of row p  |
154*     +-----------------+--------------------+------------------+
155*     |     GLP_BS      |       GLP_BS       |       N/A        |
156*     |     GLP_BS      |       GLP_NL       |      GLP_BS      |
157*     |     GLP_BS      |       GLP_NU       |      GLP_BS      |
158*     |     GLP_NS      |       GLP_BS       |      GLP_BS      |
159*     |     GLP_NS      |       GLP_NL       |      GLP_NL      |
160*     |     GLP_NS      |       GLP_NU       |      GLP_NU      |
161*     +-----------------+--------------------+------------------+
162*
163*  Value of row multiplier pi[p] in solution to the original problem
164*  is the same as in solution to the transformed problem.
165*
166*  1. In solution to the transformed problem row p and column q cannot
167*     be basic at the same time; otherwise the basis matrix would have
168*     two linear dependent columns: unity column of auxiliary variable
169*     of row p and unity column of variable s.
170*
171*  2. Though in the transformed problem row p is equality constraint,
172*     it may be basic due to primal degenerate solution.
173*
174*  RECOVERING INTERIOR-POINT SOLUTION
175*
176*  Value of row multiplier pi[p] in solution to the original problem
177*  is the same as in solution to the transformed problem.
178*
179*  RECOVERING MIP SOLUTION
180*
181*  None needed. */
182
183struct ineq_row
184{     /* inequality constraint row */
185      int p;
186      /* row reference number */
187      int s;
188      /* column reference number for slack/surplus variable */
189};
190
191static int rcv_geq_row(NPP *npp, void *info);
192
193void npp_geq_row(NPP *npp, NPPROW *p)
194{     /* process row of 'not less than' type */
195      struct ineq_row *info;
196      NPPCOL *s;
197      /* the row must have lower bound */
198      xassert(p->lb != -DBL_MAX);
199      xassert(p->lb < p->ub);
200      /* create column for surplus variable */
201      s = npp_add_col(npp);
202      s->lb = 0.0;
203      s->ub = (p->ub == +DBL_MAX ? +DBL_MAX : p->ub - p->lb);
204      /* and add it to the transformed problem */
205      npp_add_aij(npp, p, s, -1.0);
206      /* create transformation stack entry */
207      info = npp_push_tse(npp,
208         rcv_geq_row, sizeof(struct ineq_row));
209      info->p = p->i;
210      info->s = s->j;
211      /* replace the row by equality constraint */
212      p->ub = p->lb;
213      return;
214}
215
216static int rcv_geq_row(NPP *npp, void *_info)
217{     /* recover row of 'not less than' type */
218      struct ineq_row *info = _info;
219      if (npp->sol == GLP_SOL)
220      {  if (npp->r_stat[info->p] == GLP_BS)
221         {  if (npp->c_stat[info->s] == GLP_BS)
222            {  npp_error();
223               return 1;
224            }
225            else if (npp->c_stat[info->s] == GLP_NL ||
226                     npp->c_stat[info->s] == GLP_NU)
227               npp->r_stat[info->p] = GLP_BS;
228            else
229            {  npp_error();
230               return 1;
231            }
232         }
233         else if (npp->r_stat[info->p] == GLP_NS)
234         {  if (npp->c_stat[info->s] == GLP_BS)
235               npp->r_stat[info->p] = GLP_BS;
236            else if (npp->c_stat[info->s] == GLP_NL)
237               npp->r_stat[info->p] = GLP_NL;
238            else if (npp->c_stat[info->s] == GLP_NU)
239               npp->r_stat[info->p] = GLP_NU;
240            else
241            {  npp_error();
242               return 1;
243            }
244         }
245         else
246         {  npp_error();
247            return 1;
248         }
249      }
250      return 0;
251}
252
253/***********************************************************************
254*  NAME
255*
256*  npp_leq_row - process row of 'not greater than' type
257*
258*  SYNOPSIS
259*
260*  #include "glpnpp.h"
261*  void npp_leq_row(NPP *npp, NPPROW *p);
262*
263*  DESCRIPTION
264*
265*  The routine npp_leq_row processes row p, which is 'not greater than'
266*  inequality constraint:
267*
268*     (L[p] <=) sum a[p,j] x[j] <= U[p],                             (1)
269*                j
270*
271*  where L[p] < U[p], and lower bound may not exist (L[p] = +oo).
272*
273*  PROBLEM TRANSFORMATION
274*
275*  Constraint (1) can be replaced by equality constraint:
276*
277*     sum a[p,j] x[j] + s = L[p],                                    (2)
278*      j
279*
280*  where
281*
282*     0 <= s (<= U[p] - L[p])                                        (3)
283*
284*  is a non-negative slack variable.
285*
286*  Since in the primal system there appears column s having the only
287*  non-zero coefficient in row p, in the dual system there appears a
288*  new row:
289*
290*     (+1) pi[p] + lambda = 0,                                       (4)
291*
292*  where (+1) is coefficient of column s in row p, pi[p] is multiplier
293*  of row p, lambda is multiplier of column q, 0 is coefficient of
294*  column s in the objective row.
295*
296*  RECOVERING BASIC SOLUTION
297*
298*  Status of row p in solution to the original problem is determined
299*  by its status and status of column q in solution to the transformed
300*  problem as follows:
301*
302*     +--------------------------------------+------------------+
303*     |         Transformed problem          | Original problem |
304*     +-----------------+--------------------+------------------+
305*     | Status of row p | Status of column s | Status of row p  |
306*     +-----------------+--------------------+------------------+
307*     |     GLP_BS      |       GLP_BS       |       N/A        |
308*     |     GLP_BS      |       GLP_NL       |      GLP_BS      |
309*     |     GLP_BS      |       GLP_NU       |      GLP_BS      |
310*     |     GLP_NS      |       GLP_BS       |      GLP_BS      |
311*     |     GLP_NS      |       GLP_NL       |      GLP_NU      |
312*     |     GLP_NS      |       GLP_NU       |      GLP_NL      |
313*     +-----------------+--------------------+------------------+
314*
315*  Value of row multiplier pi[p] in solution to the original problem
316*  is the same as in solution to the transformed problem.
317*
318*  1. In solution to the transformed problem row p and column q cannot
319*     be basic at the same time; otherwise the basis matrix would have
320*     two linear dependent columns: unity column of auxiliary variable
321*     of row p and unity column of variable s.
322*
323*  2. Though in the transformed problem row p is equality constraint,
324*     it may be basic due to primal degeneracy.
325*
326*  RECOVERING INTERIOR-POINT SOLUTION
327*
328*  Value of row multiplier pi[p] in solution to the original problem
329*  is the same as in solution to the transformed problem.
330*
331*  RECOVERING MIP SOLUTION
332*
333*  None needed. */
334
335static int rcv_leq_row(NPP *npp, void *info);
336
337void npp_leq_row(NPP *npp, NPPROW *p)
338{     /* process row of 'not greater than' type */
339      struct ineq_row *info;
340      NPPCOL *s;
341      /* the row must have upper bound */
342      xassert(p->ub != +DBL_MAX);
343      xassert(p->lb < p->ub);
344      /* create column for slack variable */
345      s = npp_add_col(npp);
346      s->lb = 0.0;
347      s->ub = (p->lb == -DBL_MAX ? +DBL_MAX : p->ub - p->lb);
348      /* and add it to the transformed problem */
349      npp_add_aij(npp, p, s, +1.0);
350      /* create transformation stack entry */
351      info = npp_push_tse(npp,
352         rcv_leq_row, sizeof(struct ineq_row));
353      info->p = p->i;
354      info->s = s->j;
355      /* replace the row by equality constraint */
356      p->lb = p->ub;
357      return;
358}
359
360static int rcv_leq_row(NPP *npp, void *_info)
361{     /* recover row of 'not greater than' type */
362      struct ineq_row *info = _info;
363      if (npp->sol == GLP_SOL)
364      {  if (npp->r_stat[info->p] == GLP_BS)
365         {  if (npp->c_stat[info->s] == GLP_BS)
366            {  npp_error();
367               return 1;
368            }
369            else if (npp->c_stat[info->s] == GLP_NL ||
370                     npp->c_stat[info->s] == GLP_NU)
371               npp->r_stat[info->p] = GLP_BS;
372            else
373            {  npp_error();
374               return 1;
375            }
376         }
377         else if (npp->r_stat[info->p] == GLP_NS)
378         {  if (npp->c_stat[info->s] == GLP_BS)
379               npp->r_stat[info->p] = GLP_BS;
380            else if (npp->c_stat[info->s] == GLP_NL)
381               npp->r_stat[info->p] = GLP_NU;
382            else if (npp->c_stat[info->s] == GLP_NU)
383               npp->r_stat[info->p] = GLP_NL;
384            else
385            {  npp_error();
386               return 1;
387            }
388         }
389         else
390         {  npp_error();
391            return 1;
392         }
393      }
394      return 0;
395}
396
397/***********************************************************************
398*  NAME
399*
400*  npp_free_col - process free (unbounded) column
401*
402*  SYNOPSIS
403*
404*  #include "glpnpp.h"
405*  void npp_free_col(NPP *npp, NPPCOL *q);
406*
407*  DESCRIPTION
408*
409*  The routine npp_free_col processes column q, which is free (i.e. has
410*  no finite bounds):
411*
412*     -oo < x[q] < +oo.                                              (1)
413*
414*  PROBLEM TRANSFORMATION
415*
416*  Free (unbounded) variable can be replaced by the difference of two
417*  non-negative variables:
418*
419*     x[q] = s' - s'',   s', s'' >= 0.                               (2)
420*
421*  Assuming that in the transformed problem x[q] becomes s',
422*  transformation (2) causes new column s'' to appear, which differs
423*  from column s' only in the sign of coefficients in constraint and
424*  objective rows. Thus, if in the dual system the following row
425*  corresponds to column s':
426*
427*     sum a[i,q] pi[i] + lambda' = c[q],                             (3)
428*      i
429*
430*  the row which corresponds to column s'' is the following:
431*
432*     sum (-a[i,q]) pi[i] + lambda'' = -c[q].                        (4)
433*      i
434*
435*  Then from (3) and (4) it follows that:
436*
437*     lambda' + lambda'' = 0   =>   lambda' = lmabda'' = 0,          (5)
438*
439*  where lambda' and lambda'' are multipliers for columns s' and s'',
440*  resp.
441*
442*  RECOVERING BASIC SOLUTION
443*
444*  With respect to (5) status of column q in solution to the original
445*  problem is determined by statuses of columns s' and s'' in solution
446*  to the transformed problem as follows:
447*
448*     +--------------------------------------+------------------+
449*     |         Transformed problem          | Original problem |
450*     +------------------+-------------------+------------------+
451*     | Status of col s' | Status of col s'' | Status of col q  |
452*     +------------------+-------------------+------------------+
453*     |      GLP_BS      |      GLP_BS       |       N/A        |
454*     |      GLP_BS      |      GLP_NL       |      GLP_BS      |
455*     |      GLP_NL      |      GLP_BS       |      GLP_BS      |
456*     |      GLP_NL      |      GLP_NL       |      GLP_NF      |
457*     +------------------+-------------------+------------------+
458*
459*  Value of column q is computed with formula (2).
460*
461*  1. In solution to the transformed problem columns s' and s'' cannot
462*     be basic at the same time, because they differ only in the sign,
463*     hence, are linear dependent.
464*
465*  2. Though column q is free, it can be non-basic due to dual
466*     degeneracy.
467*
468*  3. If column q is integral, columns s' and s'' are also integral.
469*
470*  RECOVERING INTERIOR-POINT SOLUTION
471*
472*  Value of column q is computed with formula (2).
473*
474*  RECOVERING MIP SOLUTION
475*
476*  Value of column q is computed with formula (2). */
477
478struct free_col
479{     /* free (unbounded) column */
480      int q;
481      /* column reference number for variables x[q] and s' */
482      int s;
483      /* column reference number for variable s'' */
484};
485
486static int rcv_free_col(NPP *npp, void *info);
487
488void npp_free_col(NPP *npp, NPPCOL *q)
489{     /* process free (unbounded) column */
490      struct free_col *info;
491      NPPCOL *s;
492      NPPAIJ *aij;
493      /* the column must be free */
494      xassert(q->lb == -DBL_MAX && q->ub == +DBL_MAX);
495      /* variable x[q] becomes s' */
496      q->lb = 0.0, q->ub = +DBL_MAX;
497      /* create variable s'' */
498      s = npp_add_col(npp);
499      s->is_int = q->is_int;
500      s->lb = 0.0, s->ub = +DBL_MAX;
501      /* duplicate objective coefficient */
502      s->coef = -q->coef;
503      /* duplicate column of the constraint matrix */
504      for (aij = q->ptr; aij != NULL; aij = aij->c_next)
505         npp_add_aij(npp, aij->row, s, -aij->val);
506      /* create transformation stack entry */
507      info = npp_push_tse(npp,
508         rcv_free_col, sizeof(struct free_col));
509      info->q = q->j;
510      info->s = s->j;
511      return;
512}
513
514static int rcv_free_col(NPP *npp, void *_info)
515{     /* recover free (unbounded) column */
516      struct free_col *info = _info;
517      if (npp->sol == GLP_SOL)
518      {  if (npp->c_stat[info->q] == GLP_BS)
519         {  if (npp->c_stat[info->s] == GLP_BS)
520            {  npp_error();
521               return 1;
522            }
523            else if (npp->c_stat[info->s] == GLP_NL)
524               npp->c_stat[info->q] = GLP_BS;
525            else
526            {  npp_error();
527               return -1;
528            }
529         }
530         else if (npp->c_stat[info->q] == GLP_NL)
531         {  if (npp->c_stat[info->s] == GLP_BS)
532               npp->c_stat[info->q] = GLP_BS;
533            else if (npp->c_stat[info->s] == GLP_NL)
534               npp->c_stat[info->q] = GLP_NF;
535            else
536            {  npp_error();
537               return -1;
538            }
539         }
540         else
541         {  npp_error();
542            return -1;
543         }
544      }
545      /* compute value of x[q] with formula (2) */
546      npp->c_value[info->q] -= npp->c_value[info->s];
547      return 0;
548}
549
550/***********************************************************************
551*  NAME
552*
553*  npp_lbnd_col - process column with (non-zero) lower bound
554*
555*  SYNOPSIS
556*
557*  #include "glpnpp.h"
558*  void npp_lbnd_col(NPP *npp, NPPCOL *q);
559*
560*  DESCRIPTION
561*
562*  The routine npp_lbnd_col processes column q, which has (non-zero)
563*  lower bound:
564*
565*     l[q] <= x[q] (<= u[q]),                                        (1)
566*
567*  where l[q] < u[q], and upper bound may not exist (u[q] = +oo).
568*
569*  PROBLEM TRANSFORMATION
570*
571*  Column q can be replaced as follows:
572*
573*     x[q] = l[q] + s,                                               (2)
574*
575*  where
576*
577*     0 <= s (<= u[q] - l[q])                                        (3)
578*
579*  is a non-negative variable.
580*
581*  Substituting x[q] from (2) into the objective row, we have:
582*
583*     z = sum c[j] x[j] + c0 =
584*          j
585*
586*       = sum c[j] x[j] + c[q] x[q] + c0 =
587*         j!=q
588*
589*       = sum c[j] x[j] + c[q] (l[q] + s) + c0 =
590*         j!=q
591*
592*       = sum c[j] x[j] + c[q] s + c~0,
593*
594*  where
595*
596*     c~0 = c0 + c[q] l[q]                                           (4)
597*
598*  is the constant term of the objective in the transformed problem.
599*  Similarly, substituting x[q] into constraint row i, we have:
600*
601*     L[i] <= sum a[i,j] x[j] <= U[i]  ==>
602*              j
603*
604*     L[i] <= sum a[i,j] x[j] + a[i,q] x[q] <= U[i]  ==>
605*             j!=q
606*
607*     L[i] <= sum a[i,j] x[j] + a[i,q] (l[q] + s) <= U[i]  ==>
608*             j!=q
609*
610*     L~[i] <= sum a[i,j] x[j] + a[i,q] s <= U~[i],
611*              j!=q
612*
613*  where
614*
615*     L~[i] = L[i] - a[i,q] l[q],  U~[i] = U[i] - a[i,q] l[q]        (5)
616*
617*  are lower and upper bounds of row i in the transformed problem,
618*  resp.
619*
620*  Transformation (2) does not affect the dual system.
621*
622*  RECOVERING BASIC SOLUTION
623*
624*  Status of column q in solution to the original problem is the same
625*  as in solution to the transformed problem (GLP_BS, GLP_NL or GLP_NU).
626*  Value of column q is computed with formula (2).
627*
628*  RECOVERING INTERIOR-POINT SOLUTION
629*
630*  Value of column q is computed with formula (2).
631*
632*  RECOVERING MIP SOLUTION
633*
634*  Value of column q is computed with formula (2). */
635
636struct bnd_col
637{     /* bounded column */
638      int q;
639      /* column reference number for variables x[q] and s */
640      double bnd;
641      /* lower/upper bound l[q] or u[q] */
642};
643
644static int rcv_lbnd_col(NPP *npp, void *info);
645
646void npp_lbnd_col(NPP *npp, NPPCOL *q)
647{     /* process column with (non-zero) lower bound */
648      struct bnd_col *info;
649      NPPROW *i;
650      NPPAIJ *aij;
651      /* the column must have non-zero lower bound */
652      xassert(q->lb != 0.0);
653      xassert(q->lb != -DBL_MAX);
654      xassert(q->lb < q->ub);
655      /* create transformation stack entry */
656      info = npp_push_tse(npp,
657         rcv_lbnd_col, sizeof(struct bnd_col));
658      info->q = q->j;
659      info->bnd = q->lb;
660      /* substitute x[q] into objective row */
661      npp->c0 += q->coef * q->lb;
662      /* substitute x[q] into constraint rows */
663      for (aij = q->ptr; aij != NULL; aij = aij->c_next)
664      {  i = aij->row;
665         if (i->lb == i->ub)
666            i->ub = (i->lb -= aij->val * q->lb);
667         else
668         {  if (i->lb != -DBL_MAX)
669               i->lb -= aij->val * q->lb;
670            if (i->ub != +DBL_MAX)
671               i->ub -= aij->val * q->lb;
672         }
673      }
674      /* column x[q] becomes column s */
675      if (q->ub != +DBL_MAX)
676         q->ub -= q->lb;
677      q->lb = 0.0;
678      return;
679}
680
681static int rcv_lbnd_col(NPP *npp, void *_info)
682{     /* recover column with (non-zero) lower bound */
683      struct bnd_col *info = _info;
684      if (npp->sol == GLP_SOL)
685      {  if (npp->c_stat[info->q] == GLP_BS ||
686             npp->c_stat[info->q] == GLP_NL ||
687             npp->c_stat[info->q] == GLP_NU)
688            npp->c_stat[info->q] = npp->c_stat[info->q];
689         else
690         {  npp_error();
691            return 1;
692         }
693      }
694      /* compute value of x[q] with formula (2) */
695      npp->c_value[info->q] = info->bnd + npp->c_value[info->q];
696      return 0;
697}
698
699/***********************************************************************
700*  NAME
701*
702*  npp_ubnd_col - process column with upper bound
703*
704*  SYNOPSIS
705*
706*  #include "glpnpp.h"
707*  void npp_ubnd_col(NPP *npp, NPPCOL *q);
708*
709*  DESCRIPTION
710*
711*  The routine npp_ubnd_col processes column q, which has upper bound:
712*
713*     (l[q] <=) x[q] <= u[q],                                        (1)
714*
715*  where l[q] < u[q], and lower bound may not exist (l[q] = -oo).
716*
717*  PROBLEM TRANSFORMATION
718*
719*  Column q can be replaced as follows:
720*
721*     x[q] = u[q] - s,                                               (2)
722*
723*  where
724*
725*     0 <= s (<= u[q] - l[q])                                        (3)
726*
727*  is a non-negative variable.
728*
729*  Substituting x[q] from (2) into the objective row, we have:
730*
731*     z = sum c[j] x[j] + c0 =
732*          j
733*
734*       = sum c[j] x[j] + c[q] x[q] + c0 =
735*         j!=q
736*
737*       = sum c[j] x[j] + c[q] (u[q] - s) + c0 =
738*         j!=q
739*
740*       = sum c[j] x[j] - c[q] s + c~0,
741*
742*  where
743*
744*     c~0 = c0 + c[q] u[q]                                           (4)
745*
746*  is the constant term of the objective in the transformed problem.
747*  Similarly, substituting x[q] into constraint row i, we have:
748*
749*     L[i] <= sum a[i,j] x[j] <= U[i]  ==>
750*              j
751*
752*     L[i] <= sum a[i,j] x[j] + a[i,q] x[q] <= U[i]  ==>
753*             j!=q
754*
755*     L[i] <= sum a[i,j] x[j] + a[i,q] (u[q] - s) <= U[i]  ==>
756*             j!=q
757*
758*     L~[i] <= sum a[i,j] x[j] - a[i,q] s <= U~[i],
759*              j!=q
760*
761*  where
762*
763*     L~[i] = L[i] - a[i,q] u[q],  U~[i] = U[i] - a[i,q] u[q]        (5)
764*
765*  are lower and upper bounds of row i in the transformed problem,
766*  resp.
767*
768*  Note that in the transformed problem coefficients c[q] and a[i,q]
769*  change their sign. Thus, the row of the dual system corresponding to
770*  column q:
771*
772*     sum a[i,q] pi[i] + lambda[q] = c[q]                            (6)
773*      i
774*
775*  in the transformed problem becomes the following:
776*
777*     sum (-a[i,q]) pi[i] + lambda[s] = -c[q].                       (7)
778*      i
779*
780*  Therefore:
781*
782*     lambda[q] = - lambda[s],                                       (8)
783*
784*  where lambda[q] is multiplier for column q, lambda[s] is multiplier
785*  for column s.
786*
787*  RECOVERING BASIC SOLUTION
788*
789*  With respect to (8) status of column q in solution to the original
790*  problem is determined by status of column s in solution to the
791*  transformed problem as follows:
792*
793*     +-----------------------+--------------------+
794*     |  Status of column s   | Status of column q |
795*     | (transformed problem) | (original problem) |
796*     +-----------------------+--------------------+
797*     |        GLP_BS         |       GLP_BS       |
798*     |        GLP_NL         |       GLP_NU       |
799*     |        GLP_NU         |       GLP_NL       |
800*     +-----------------------+--------------------+
801*
802*  Value of column q is computed with formula (2).
803*
804*  RECOVERING INTERIOR-POINT SOLUTION
805*
806*  Value of column q is computed with formula (2).
807*
808*  RECOVERING MIP SOLUTION
809*
810*  Value of column q is computed with formula (2). */
811
812static int rcv_ubnd_col(NPP *npp, void *info);
813
814void npp_ubnd_col(NPP *npp, NPPCOL *q)
815{     /* process column with upper bound */
816      struct bnd_col *info;
817      NPPROW *i;
818      NPPAIJ *aij;
819      /* the column must have upper bound */
820      xassert(q->ub != +DBL_MAX);
821      xassert(q->lb < q->ub);
822      /* create transformation stack entry */
823      info = npp_push_tse(npp,
824         rcv_ubnd_col, sizeof(struct bnd_col));
825      info->q = q->j;
826      info->bnd = q->ub;
827      /* substitute x[q] into objective row */
828      npp->c0 += q->coef * q->ub;
829      q->coef = -q->coef;
830      /* substitute x[q] into constraint rows */
831      for (aij = q->ptr; aij != NULL; aij = aij->c_next)
832      {  i = aij->row;
833         if (i->lb == i->ub)
834            i->ub = (i->lb -= aij->val * q->ub);
835         else
836         {  if (i->lb != -DBL_MAX)
837               i->lb -= aij->val * q->ub;
838            if (i->ub != +DBL_MAX)
839               i->ub -= aij->val * q->ub;
840         }
841         aij->val = -aij->val;
842      }
843      /* column x[q] becomes column s */
844      if (q->lb != -DBL_MAX)
845         q->ub -= q->lb;
846      else
847         q->ub = +DBL_MAX;
848      q->lb = 0.0;
849      return;
850}
851
852static int rcv_ubnd_col(NPP *npp, void *_info)
853{     /* recover column with upper bound */
854      struct bnd_col *info = _info;
855      if (npp->sol == GLP_BS)
856      {  if (npp->c_stat[info->q] == GLP_BS)
857            npp->c_stat[info->q] = GLP_BS;
858         else if (npp->c_stat[info->q] == GLP_NL)
859            npp->c_stat[info->q] = GLP_NU;
860         else if (npp->c_stat[info->q] == GLP_NU)
861            npp->c_stat[info->q] = GLP_NL;
862         else
863         {  npp_error();
864            return 1;
865         }
866      }
867      /* compute value of x[q] with formula (2) */
868      npp->c_value[info->q] = info->bnd - npp->c_value[info->q];
869      return 0;
870}
871
872/***********************************************************************
873*  NAME
874*
875*  npp_dbnd_col - process non-negative column with upper bound
876*
877*  SYNOPSIS
878*
879*  #include "glpnpp.h"
880*  void npp_dbnd_col(NPP *npp, NPPCOL *q);
881*
882*  DESCRIPTION
883*
884*  The routine npp_dbnd_col processes column q, which is non-negative
885*  and has upper bound:
886*
887*     0 <= x[q] <= u[q],                                             (1)
888*
889*  where u[q] > 0.
890*
891*  PROBLEM TRANSFORMATION
892*
893*  Upper bound of column q can be replaced by the following equality
894*  constraint:
895*
896*     x[q] + s = u[q],                                               (2)
897*
898*  where s >= 0 is a non-negative complement variable.
899*
900*  Since in the primal system along with new row (2) there appears a
901*  new column s having the only non-zero coefficient in this row, in
902*  the dual system there appears a new row:
903*
904*     (+1)pi + lambda[s] = 0,                                        (3)
905*
906*  where (+1) is coefficient at column s in row (2), pi is multiplier
907*  for row (2), lambda[s] is multiplier for column s, 0 is coefficient
908*  at column s in the objective row.
909*
910*  RECOVERING BASIC SOLUTION
911*
912*  Status of column q in solution to the original problem is determined
913*  by its status and status of column s in solution to the transformed
914*  problem as follows:
915*
916*     +-----------------------------------+------------------+
917*     |         Transformed problem       | Original problem |
918*     +-----------------+-----------------+------------------+
919*     | Status of col q | Status of col s | Status of col q  |
920*     +-----------------+-----------------+------------------+
921*     |     GLP_BS      |     GLP_BS      |      GLP_BS      |
922*     |     GLP_BS      |     GLP_NL      |      GLP_NU      |
923*     |     GLP_NL      |     GLP_BS      |      GLP_NL      |
924*     |     GLP_NL      |     GLP_NL      |      GLP_NL (*)  |
925*     +-----------------+-----------------+------------------+
926*
927*  Value of column q in solution to the original problem is the same as
928*  in solution to the transformed problem.
929*
930*  1. Formally, in solution to the transformed problem columns q and s
931*     cannot be non-basic at the same time, since the constraint (2)
932*     would be violated. However, if u[q] is close to zero, violation
933*     may be less than a working precision even if both columns q and s
934*     are non-basic. In this degenerate case row (2) can be only basic,
935*     i.e. non-active constraint (otherwise corresponding row of the
936*     basis matrix would be zero). This allows to pivot out auxiliary
937*     variable and pivot in column s, in which case the row becomes
938*     active while column s becomes basic.
939*
940*  2. If column q is integral, column s is also integral.
941*
942*  RECOVERING INTERIOR-POINT SOLUTION
943*
944*  Value of column q in solution to the original problem is the same as
945*  in solution to the transformed problem.
946*
947*  RECOVERING MIP SOLUTION
948*
949*  Value of column q in solution to the original problem is the same as
950*  in solution to the transformed problem. */
951
952struct dbnd_col
953{     /* double-bounded column */
954      int q;
955      /* column reference number for variable x[q] */
956      int s;
957      /* column reference number for complement variable s */
958};
959
960static int rcv_dbnd_col(NPP *npp, void *info);
961
962void npp_dbnd_col(NPP *npp, NPPCOL *q)
963{     /* process non-negative column with upper bound */
964      struct dbnd_col *info;
965      NPPROW *p;
966      NPPCOL *s;
967      /* the column must be non-negative with upper bound */
968      xassert(q->lb == 0.0);
969      xassert(q->ub > 0.0);
970      xassert(q->ub != +DBL_MAX);
971      /* create variable s */
972      s = npp_add_col(npp);
973      s->is_int = q->is_int;
974      s->lb = 0.0, s->ub = +DBL_MAX;
975      /* create equality constraint (2) */
976      p = npp_add_row(npp);
977      p->lb = p->ub = q->ub;
978      npp_add_aij(npp, p, q, +1.0);
979      npp_add_aij(npp, p, s, +1.0);
980      /* create transformation stack entry */
981      info = npp_push_tse(npp,
982         rcv_dbnd_col, sizeof(struct dbnd_col));
983      info->q = q->j;
984      info->s = s->j;
985      /* remove upper bound of x[q] */
986      q->ub = +DBL_MAX;
987      return;
988}
989
990static int rcv_dbnd_col(NPP *npp, void *_info)
991{     /* recover non-negative column with upper bound */
992      struct dbnd_col *info = _info;
993      if (npp->sol == GLP_BS)
994      {  if (npp->c_stat[info->q] == GLP_BS)
995         {  if (npp->c_stat[info->s] == GLP_BS)
996               npp->c_stat[info->q] = GLP_BS;
997            else if (npp->c_stat[info->s] == GLP_NL)
998               npp->c_stat[info->q] = GLP_NU;
999            else
1000            {  npp_error();
1001               return 1;
1002            }
1003         }
1004         else if (npp->c_stat[info->q] == GLP_NL)
1005         {  if (npp->c_stat[info->s] == GLP_BS ||
1006                npp->c_stat[info->s] == GLP_NL)
1007               npp->c_stat[info->q] = GLP_NL;
1008            else
1009            {  npp_error();
1010               return 1;
1011            }
1012         }
1013         else
1014         {  npp_error();
1015            return 1;
1016         }
1017      }
1018      return 0;
1019}
1020
1021/***********************************************************************
1022*  NAME
1023*
1024*  npp_fixed_col - process fixed column
1025*
1026*  SYNOPSIS
1027*
1028*  #include "glpnpp.h"
1029*  void npp_fixed_col(NPP *npp, NPPCOL *q);
1030*
1031*  DESCRIPTION
1032*
1033*  The routine npp_fixed_col processes column q, which is fixed:
1034*
1035*     x[q] = s[q],                                                   (1)
1036*
1037*  where s[q] is a fixed column value.
1038*
1039*  PROBLEM TRANSFORMATION
1040*
1041*  The value of a fixed column can be substituted into the objective
1042*  and constraint rows that allows removing the column from the problem.
1043*
1044*  Substituting x[q] = s[q] into the objective row, we have:
1045*
1046*     z = sum c[j] x[j] + c0 =
1047*          j
1048*
1049*       = sum c[j] x[j] + c[q] x[q] + c0 =
1050*         j!=q
1051*
1052*       = sum c[j] x[j] + c[q] s[q] + c0 =
1053*         j!=q
1054*
1055*       = sum c[j] x[j] + c~0,
1056*         j!=q
1057*
1058*  where
1059*
1060*     c~0 = c0 + c[q] s[q]                                           (2)
1061*
1062*  is the constant term of the objective in the transformed problem.
1063*  Similarly, substituting x[q] = s[q] into constraint row i, we have:
1064*
1065*     L[i] <= sum a[i,j] x[j] <= U[i]  ==>
1066*              j
1067*
1068*     L[i] <= sum a[i,j] x[j] + a[i,q] x[q] <= U[i]  ==>
1069*             j!=q
1070*
1071*     L[i] <= sum a[i,j] x[j] + a[i,q] s[q] <= U[i]  ==>
1072*             j!=q
1073*
1074*     L~[i] <= sum a[i,j] x[j] + a[i,q] s <= U~[i],
1075*              j!=q
1076*
1077*  where
1078*
1079*     L~[i] = L[i] - a[i,q] s[q],  U~[i] = U[i] - a[i,q] s[q]        (3)
1080*
1081*  are lower and upper bounds of row i in the transformed problem,
1082*  resp.
1083*
1084*  RECOVERING BASIC SOLUTION
1085*
1086*  Column q is assigned status GLP_NS and its value is assigned s[q].
1087*
1088*  RECOVERING INTERIOR-POINT SOLUTION
1089*
1090*  Value of column q is assigned s[q].
1091*
1092*  RECOVERING MIP SOLUTION
1093*
1094*  Value of column q is assigned s[q]. */
1095
1096struct fixed_col
1097{     /* fixed column */
1098      int q;
1099      /* column reference number for variable x[q] */
1100      double s;
1101      /* value, at which x[q] is fixed */
1102};
1103
1104static int rcv_fixed_col(NPP *npp, void *info);
1105
1106void npp_fixed_col(NPP *npp, NPPCOL *q)
1107{     /* process fixed column */
1108      struct fixed_col *info;
1109      NPPROW *i;
1110      NPPAIJ *aij;
1111      /* the column must be fixed */
1112      xassert(q->lb == q->ub);
1113      /* create transformation stack entry */
1114      info = npp_push_tse(npp,
1115         rcv_fixed_col, sizeof(struct fixed_col));
1116      info->q = q->j;
1117      info->s = q->lb;
1118      /* substitute x[q] = s[q] into objective row */
1119      npp->c0 += q->coef * q->lb;
1120      /* substitute x[q] = s[q] into constraint rows */
1121      for (aij = q->ptr; aij != NULL; aij = aij->c_next)
1122      {  i = aij->row;
1123         if (i->lb == i->ub)
1124            i->ub = (i->lb -= aij->val * q->lb);
1125         else
1126         {  if (i->lb != -DBL_MAX)
1127               i->lb -= aij->val * q->lb;
1128            if (i->ub != +DBL_MAX)
1129               i->ub -= aij->val * q->lb;
1130         }
1131      }
1132      /* remove the column from the problem */
1133      npp_del_col(npp, q);
1134      return;
1135}
1136
1137static int rcv_fixed_col(NPP *npp, void *_info)
1138{     /* recover fixed column */
1139      struct fixed_col *info = _info;
1140      if (npp->sol == GLP_SOL)
1141         npp->c_stat[info->q] = GLP_NS;
1142      npp->c_value[info->q] = info->s;
1143      return 0;
1144}
1145
1146/***********************************************************************
1147*  NAME
1148*
1149*  npp_make_equality - process row with almost identical bounds
1150*
1151*  SYNOPSIS
1152*
1153*  #include "glpnpp.h"
1154*  int npp_make_equality(NPP *npp, NPPROW *p);
1155*
1156*  DESCRIPTION
1157*
1158*  The routine npp_make_equality processes row p:
1159*
1160*     L[p] <= sum a[p,j] x[j] <= U[p],                               (1)
1161*              j
1162*
1163*  where -oo < L[p] < U[p] < +oo, i.e. which is double-sided inequality
1164*  constraint.
1165*
1166*  RETURNS
1167*
1168*  0 - row bounds have not been changed;
1169*
1170*  1 - row has been replaced by equality constraint.
1171*
1172*  PROBLEM TRANSFORMATION
1173*
1174*  If bounds of row (1) are very close to each other:
1175*
1176*     U[p] - L[p] <= eps,                                            (2)
1177*
1178*  where eps is an absolute tolerance for row value, the row can be
1179*  replaced by the following almost equivalent equiality constraint:
1180*
1181*     sum a[p,j] x[j] = b,                                           (3)
1182*      j
1183*
1184*  where b = (L[p] + U[p]) / 2. If the right-hand side in (3) happens
1185*  to be very close to its nearest integer:
1186*
1187*     |b - floor(b + 0.5)| <= eps,                                   (4)
1188*
1189*  it is reasonable to use this nearest integer as the right-hand side.
1190*
1191*  RECOVERING BASIC SOLUTION
1192*
1193*  Status of row p in solution to the original problem is determined
1194*  by its status and the sign of its multiplier pi[p] in solution to
1195*  the transformed problem as follows:
1196*
1197*     +-----------------------+---------+--------------------+
1198*     |    Status of row p    | Sign of |  Status of row p   |
1199*     | (transformed problem) |  pi[p]  | (original problem) |
1200*     +-----------------------+---------+--------------------+
1201*     |        GLP_BS         |  + / -  |       GLP_BS       |
1202*     |        GLP_NS         |    +    |       GLP_NL       |
1203*     |        GLP_NS         |    -    |       GLP_NU       |
1204*     +-----------------------+---------+--------------------+
1205*
1206*  Value of row multiplier pi[p] in solution to the original problem is
1207*  the same as in solution to the transformed problem.
1208*
1209*  RECOVERING INTERIOR POINT SOLUTION
1210*
1211*  Value of row multiplier pi[p] in solution to the original problem is
1212*  the same as in solution to the transformed problem.
1213*
1214*  RECOVERING MIP SOLUTION
1215*
1216*  None needed. */
1217
1218struct make_equality
1219{     /* row with almost identical bounds */
1220      int p;
1221      /* row reference number */
1222};
1223
1224static int rcv_make_equality(NPP *npp, void *info);
1225
1226int npp_make_equality(NPP *npp, NPPROW *p)
1227{     /* process row with almost identical bounds */
1228      struct make_equality *info;
1229      double b, eps, nint;
1230      /* the row must be double-sided inequality */
1231      xassert(p->lb != -DBL_MAX);
1232      xassert(p->ub != +DBL_MAX);
1233      xassert(p->lb < p->ub);
1234      /* check row bounds */
1235      eps = 1e-9 + 1e-12 * fabs(p->lb);
1236      if (p->ub - p->lb > eps) return 0;
1237      /* row bounds are very close to each other */
1238      /* create transformation stack entry */
1239      info = npp_push_tse(npp,
1240         rcv_make_equality, sizeof(struct make_equality));
1241      info->p = p->i;
1242      /* compute right-hand side */
1243      b = 0.5 * (p->ub + p->lb);
1244      nint = floor(b + 0.5);
1245      if (fabs(b - nint) <= eps) b = nint;
1246      /* replace row p by almost equivalent equality constraint */
1247      p->lb = p->ub = b;
1248      return 1;
1249}
1250
1251int rcv_make_equality(NPP *npp, void *_info)
1252{     /* recover row with almost identical bounds */
1253      struct make_equality *info = _info;
1254      if (npp->sol == GLP_SOL)
1255      {  if (npp->r_stat[info->p] == GLP_BS)
1256            npp->r_stat[info->p] = GLP_BS;
1257         else if (npp->r_stat[info->p] == GLP_NS)
1258         {  if (npp->r_pi[info->p] >= 0.0)
1259               npp->r_stat[info->p] = GLP_NL;
1260            else
1261               npp->r_stat[info->p] = GLP_NU;
1262         }
1263         else
1264         {  npp_error();
1265            return 1;
1266         }
1267      }
1268      return 0;
1269}
1270
1271/***********************************************************************
1272*  NAME
1273*
1274*  npp_make_fixed - process column with almost identical bounds
1275*
1276*  SYNOPSIS
1277*
1278*  #include "glpnpp.h"
1279*  int npp_make_fixed(NPP *npp, NPPCOL *q);
1280*
1281*  DESCRIPTION
1282*
1283*  The routine npp_make_fixed processes column q:
1284*
1285*     l[q] <= x[q] <= u[q],                                          (1)
1286*
1287*  where -oo < l[q] < u[q] < +oo, i.e. which has both lower and upper
1288*  bounds.
1289*
1290*  RETURNS
1291*
1292*  0 - column bounds have not been changed;
1293*
1294*  1 - column has been fixed.
1295*
1296*  PROBLEM TRANSFORMATION
1297*
1298*  If bounds of column (1) are very close to each other:
1299*
1300*     u[q] - l[q] <= eps,                                            (2)
1301*
1302*  where eps is an absolute tolerance for column value, the column can
1303*  be fixed:
1304*
1305*     x[q] = s[q],                                                   (3)
1306*
1307*  where s[q] = (l[q] + u[q]) / 2. And if the fixed column value s[q]
1308*  happens to be very close to its nearest integer:
1309*
1310*     |s[q] - floor(s[q] + 0.5)| <= eps,                             (4)
1311*
1312*  it is reasonable to use this nearest integer as the fixed value.
1313*
1314*  RECOVERING BASIC SOLUTION
1315*
1316*  In the dual system of the original (as well as transformed) problem
1317*  column q corresponds to the following row:
1318*
1319*     sum a[i,q] pi[i] + lambda[q] = c[q].                           (5)
1320*      i
1321*
1322*  Since multipliers pi[i] are known for all rows from solution to the
1323*  transformed problem, formula (5) allows computing value of multiplier
1324*  (reduced cost) for column q:
1325*
1326*     lambda[q] = c[q] - sum a[i,q] pi[i].                           (6)
1327*                         i
1328*
1329*  Status of column q in solution to the original problem is determined
1330*  by its status and the sign of its multiplier lambda[q] in solution to
1331*  the transformed problem as follows:
1332*
1333*     +-----------------------+-----------+--------------------+
1334*     |  Status of column q   |  Sign of  | Status of column q |
1335*     | (transformed problem) | lambda[q] | (original problem) |
1336*     +-----------------------+-----------+--------------------+
1337*     |        GLP_BS         |   + / -   |       GLP_BS       |
1338*     |        GLP_NS         |     +     |       GLP_NL       |
1339*     |        GLP_NS         |     -     |       GLP_NU       |
1340*     +-----------------------+-----------+--------------------+
1341*
1342*  Value of column q in solution to the original problem is the same as
1343*  in solution to the transformed problem.
1344*
1345*  RECOVERING INTERIOR POINT SOLUTION
1346*
1347*  Value of column q in solution to the original problem is the same as
1348*  in solution to the transformed problem.
1349*
1350*  RECOVERING MIP SOLUTION
1351*
1352*  None needed. */
1353
1354struct make_fixed
1355{     /* column with almost identical bounds */
1356      int q;
1357      /* column reference number */
1358      double c;
1359      /* objective coefficient at x[q] */
1360      NPPLFE *ptr;
1361      /* list of non-zero coefficients a[i,q] */
1362};
1363
1364static int rcv_make_fixed(NPP *npp, void *info);
1365
1366int npp_make_fixed(NPP *npp, NPPCOL *q)
1367{     /* process column with almost identical bounds */
1368      struct make_fixed *info;
1369      NPPAIJ *aij;
1370      NPPLFE *lfe;
1371      double s, eps, nint;
1372      /* the column must be double-bounded */
1373      xassert(q->lb != -DBL_MAX);
1374      xassert(q->ub != +DBL_MAX);
1375      xassert(q->lb < q->ub);
1376      /* check column bounds */
1377      eps = 1e-9 + 1e-12 * fabs(q->lb);
1378      if (q->ub - q->lb > eps) return 0;
1379      /* column bounds are very close to each other */
1380      /* create transformation stack entry */
1381      info = npp_push_tse(npp,
1382         rcv_make_fixed, sizeof(struct make_fixed));
1383      info->q = q->j;
1384      info->c = q->coef;
1385      info->ptr = NULL;
1386      /* save column coefficients a[i,q] (needed for basic solution
1387         only) */
1388      if (npp->sol == GLP_SOL)
1389      {  for (aij = q->ptr; aij != NULL; aij = aij->c_next)
1390         {  lfe = dmp_get_atom(npp->stack, sizeof(NPPLFE));
1391            lfe->ref = aij->row->i;
1392            lfe->val = aij->val;
1393            lfe->next = info->ptr;
1394            info->ptr = lfe;
1395         }
1396      }
1397      /* compute column fixed value */
1398      s = 0.5 * (q->ub + q->lb);
1399      nint = floor(s + 0.5);
1400      if (fabs(s - nint) <= eps) s = nint;
1401      /* make column q fixed */
1402      q->lb = q->ub = s;
1403      return 1;
1404}
1405
1406static int rcv_make_fixed(NPP *npp, void *_info)
1407{     /* recover column with almost identical bounds */
1408      struct make_fixed *info = _info;
1409      NPPLFE *lfe;
1410      double lambda;
1411      if (npp->sol == GLP_SOL)
1412      {  if (npp->c_stat[info->q] == GLP_BS)
1413            npp->c_stat[info->q] = GLP_BS;
1414         else if (npp->c_stat[info->q] == GLP_NS)
1415         {  /* compute multiplier for column q with formula (6) */
1416            lambda = info->c;
1417            for (lfe = info->ptr; lfe != NULL; lfe = lfe->next)
1418               lambda -= lfe->val * npp->r_pi[lfe->ref];
1419            /* assign status to non-basic column */
1420            if (lambda >= 0.0)
1421               npp->c_stat[info->q] = GLP_NL;
1422            else
1423               npp->c_stat[info->q] = GLP_NU;
1424         }
1425         else
1426         {  npp_error();
1427            return 1;
1428         }
1429      }
1430      return 0;
1431}
1432
1433/* eof */
Note: See TracBrowser for help on using the repository browser.