source: branches/CglWorking101215/src/CglProbing/CglProbingProbe.cpp @ 960

Last change on this file since 960 was 960, checked in by lou, 11 years ago

(nonworking) Code labelled MOVE_SINGLETON now a separate routine
singletonRows. Code that strengthened coefficients now a separate routine
strengthenCoeff. Nontrivial code motion in the code blocks that finish
off down feasible and down feasible, up feasible.

File size: 116.2 KB
Line 
1
2/*
3  Copyright (C) 2010, International Business Machines Corporation and others.
4  All Rights Reserved.
5
6  This code is licensed under the terms of the Eclipse Public License (EPL).
7
8  $Id$
9*/
10
11#if defined(_MSC_VER)
12// Turn off compiler warning about long names
13#  pragma warning(disable:4786)
14#endif
15
16#define PROBING100 0
17
18// #define PRINT_DEBUG
19
20#include "CoinPackedMatrix.hpp"
21#include "CglProbing.hpp"
22#include "CglProbingRowCut.hpp"
23
24// To enable debugging, set CGL_DEBUG in CglProbingDebug.hpp
25#include "CglProbingDebug.hpp"
26
27namespace {
28
29/*
30  Allocating one big block for data arrays is efficient but makes it
31  difficult for debugging tools to detect accesses outside the boundary of an
32  array. The assert checks that a double is larger than an int and that the
33  ratio is a power of two.
34*/
35#if CGL_DEBUG == 0
36
37const unsigned int DIratio = sizeof(double)/sizeof(int) ;
38
39#define ONE_ARRAY
40
41#endif
42
43/*
44  Bit constants used to specify bound information.
45
46   0x0: no bounds tightened
47   0x1: u<j> tightened
48   0x2: l<j> tightened
49   0x4: l<j> < -1e10 (-infty)
50   0x8: u<j> >  1e10 (+infty)
51*/
52const unsigned int tightenUpper = 0x1 ;
53const unsigned int tightenLower = 0x2 ;
54const unsigned int infLower = 0x4 ;
55const unsigned int infUpper = 0x8 ;
56
57/*
58  Integer and bit constants used to specify probing.
59
60  The second down pass is used to restore the down probe state when down
61  is feasible but up is not.
62
63  The first set (downProbe, upProbe, downProbe2, oneProbeTooMany) *must*
64  match the iteration number for the corresponding pass in the PROBE
65  loop. The symbolic names make the code a bit easier to read.
66*/
67const unsigned int downProbe = 0 ;
68const unsigned int upProbe = 1 ;
69const unsigned int downProbe2 = 2 ;
70const unsigned int oneProbeTooMany = 3 ;
71
72const unsigned int probeDown = tightenUpper ;
73const unsigned int probeUp = tightenLower ;
74
75const unsigned int downProbeFeas = 0x1 ;
76const unsigned int upProbeFeas = 0x2 ;
77const unsigned int downProbe2Feas = 0x4 ;
78
79}  // end file-local namespace
80
81// ===============================================
82
83
84/*
85  Run through stackC and create disaggregation cuts. See the typeset
86  documentation for the full derivation. Suppose the probe variable is x<p>
87  and we're trying to generate a disaggregation cut for target x<t>.
88
89  For a down probe, we have orig_u<p> reduced to u<p> = down = colUpper[p].
90  It may be that this causes orig_u<t> to be reduced to u<t>. The upper
91  disaggregation cut is
92    -(orig_u<t>-u<t>)x<p> + x<t> <= u<t> - u<p>(orig_u<t>-u<t>)
93
94  It may be that this causes orig_l<t> to be increased to l<t>. The lower
95  disaggregation cut is
96    -(orig_l<t>-l<t>)x<p> + x<t> >= l<t> - u<p>(orig_l<t>-l<t>)
97
98  For an up probe, we have orig_l<p> increased to l<p> = up = colLower[p].
99  It may be that this causes orig_u<t> to be reduced to u<t>. The upper
100  disaggregation cut is
101     (orig_u<t>-u<t>)x<p> + x<t> <= u<t> + l<p>(orig_u<t>-u<t>)
102
103  It may be that this causes orig_l<t> to be increased to l<t>. The lower
104  disaggregation cut is
105     (orig_l<t>-l<t>)x<p> + x<t> >= l<t> + l<p>(orig_l<t>-l<t>)
106
107  Note that stackC[0] = pndx, the index of the probed variable, so we
108  don't want to process stackC[0] in the main loop.
109
110  It must be true that goingToTrueBound == 2 and justReplace == false when
111  this method is called.
112
113  TODO: If you look at the math, the critical thing to get the cut form
114        used here is that the old and new bounds on x<p> be separated by
115        one. E.g., for a down probe, (orig_u<p>-u<p>) = 1.  It looks to
116        me like goingToTrueBound might have been an attempt to guarantee this
117        condition, but the test is wrong, hence it's patched up with tweak
118        that essentially cuts this cut back to binary probing variables.
119        -- lh, 101214 --
120
121  probeDir should be coded using probeDown, probeUp.
122*/
123void disaggCuts (int nstackC, unsigned int probeDir,
124                 double primalTolerance_, double disaggEffectiveness,
125                 const OsiSolverInterface &si,
126                 CglProbingRowCut &rowCut, const int *const stackC,
127                 const double *const colsol,
128                 const double *const colUpper, const double *const colLower,
129                 const double *const saveU, const double *const saveL,
130                 int *const index, double *const element)
131{ 
132  int pndx = stackC[0] ;
133
134# if CGL_DEBUG > 0
135  std::cout
136    << "Entering disaggCuts, probed variable "
137    << si.getColName(pndx) << "(" << pndx << ")." << std::endl ;
138  const OsiRowCutDebugger *debugger = si.getRowCutDebugger() ;
139  // down or up are the only possibilities
140  assert((probeDir&(~(probeDown|probeUp))) == 0) ;
141# endif
142
143  int plusMinus = 0 ;
144  double lu_p = 0.0 ;
145  double origsol_p = colsol[pndx] ;
146  double origdelta_p = 0.0 ;
147/*
148  Set up the quantities that vary depending on whether this is an up or down
149  probe:
150  * lu_p is the new upper (down probe) or lower (up probe) bound on x<p>
151  * plusMinus will be -1 (down probe) or +1 (up probe).
152  * origdelta_p will be x*<p>-u<p> (down probe) or l<p>-x*<p> (up probe)
153*/
154  if (probeDir == probeDown) {
155    lu_p = colUpper[pndx] ;
156    plusMinus = -1 ;
157    origdelta_p = origsol_p-lu_p ;
158  } else {
159    lu_p = colLower[pndx] ;
160    plusMinus = 1 ;
161    origdelta_p = lu_p-origsol_p ;
162  }
163
164  for (int istackC = nstackC-1 ; istackC > 0 ; istackC--) {
165    int icol = stackC[istackC] ;
166    double origu_t = saveU[istackC] ;
167    double origl_t = saveL[istackC] ;
168    double u_t = colUpper[icol] ;
169    double sol_t = colsol[icol] ;
170    double l_t = colLower[icol] ;
171    double boundChange = origu_t-u_t ;
172    double coeff = plusMinus*boundChange ;
173/*
174  Generate the upper disaggregation cut, if it's violated and the coefficients
175  will be reasonable. The effectiveness is the improvement (change) in
176  x<pndx>.
177*/
178    if (boundChange > 0.0 && origu_t < 1.0e10 &&
179        (sol_t > u_t+boundChange*origdelta_p+primalTolerance_)) {
180      OsiRowCut rc ;
181      rc.setLb(-DBL_MAX) ;
182      rc.setUb(u_t+coeff*lu_p) ;
183      index[0] = icol ;
184      element[0] = 1.0 ;
185      index[1] = pndx ;
186      element[1] = coeff ;
187/*
188  Effectiveness is how far x<p> moves. We can save effort by transforming
189  sol_p = lu_p + ((u_t-sol_t)/coeff) to delta_p = ((u_t-sol_t)/coeff).
190*/
191      double delta_p = ((u_t-sol_t)/coeff) ;
192      assert(delta_p > origdelta_p) ;
193      rc.setEffectiveness(delta_p-origdelta_p) ;
194      if (rc.effectiveness() > disaggEffectiveness) {
195        rc.setRow(2,index,element,false) ;
196#       if CGL_DEBUG > 0
197        if (debugger) assert(!debugger->invalidCut(rc)); 
198#       endif
199        rowCut.addCutIfNotDuplicate(rc) ;
200      }
201    }
202/*
203  Generate the upper disaggregation cut.
204*/
205    boundChange = origl_t-l_t ;
206    coeff = plusMinus*boundChange ;
207    if (boundChange < 0.0 && origl_t > -1.0e10 &&
208        (sol_t < l_t+boundChange*origdelta_p-primalTolerance_)) {
209      OsiRowCut rc ;
210      rc.setLb(l_t+coeff*lu_p) ;
211      rc.setUb(DBL_MAX) ;
212      index[0] = icol ;
213      element[0] = 1.0 ;
214      index[1] = pndx ;
215      element[1] = coeff ;
216      double delta_p = ((sol_t-l_t)/coeff) ;
217      assert(delta_p > origdelta_p) ;
218      rc.setEffectiveness(delta_p-origdelta_p) ;
219      if (rc.effectiveness() > disaggEffectiveness) {
220        rc.setRow(2,index,element,false) ;
221#       if CGL_DEBUG > 0
222        if (debugger) assert(!debugger->invalidCut(rc)); 
223#       endif
224        rowCut.addCutIfNotDuplicate(rc) ;
225#if 0
226        printf("%d original bounds %g, %g new Lo %g sol= %g int %d sol= %g\n",icol,origl_t,origu_t,colLower[icol],colsol[icol], pndx, colsol[pndx]) ;
227        printf("-1.0 * x(%d) + %g * y(%d) <= %g\n",
228               icol,boundChange,pndx,rc.ub()) ;
229#endif
230      }
231    }
232  }
233
234# if CGL_DEBUG > 0
235  std::cout
236    << "Leaving disaggCuts, "
237    << rowCut.numberCuts() << " cuts." << std::endl ;
238# endif
239  return ;
240}
241
242
243
244
245/*
246  Walk the column for this variable and propagate a bound change on x<j> to
247  the row bounds for rows referencing x<j>. Recall that the row bounds are
248
249    L<i> = SUM{P}(l<j>) + SUM{M}(u<j>)
250    U<i> = SUM{P}(u<j>) + SUM{M}(l<j>)
251
252  for negative coefficients M and positive coefficients P.
253  The update cases are:
254
255    column bound    a<ij>    row bound
256        l<j>         >0        L<i>
257        l<j>         <0        U<i>
258        u<j>         >0        U<i>
259        u<j>         <0        L<i>
260 
261  Movement is used to determine the column bound that's being tightened. It
262  should be > 0 to tighten l<j>, movement < 0 to tighten u<j>.  Once we've
263  selected the row bound to update, the update is always += value*movement.
264*/
265
266bool updateRowBounds (int j, double movement,
267                      const int *const columnStart,
268                      const int *const columnLength,
269                      const int *const row,
270                      const double *const columnElements, 
271                      const double *const rowUpper,
272                      const double *const rowLower,
273                      int &nstackR, int *const stackR, int *const markR,
274                      double *const minR, double *const maxR,
275                      double *const saveMin, double *const saveMax)
276{
277  bool feasible = true ;
278
279  for (int k = columnStart[j] ; k < columnStart[j]+columnLength[j] ; k++) {
280    int irow = row[k] ;
281    double value = columnElements[k] ;
282    double delta = value*movement ;
283/*
284  markR is initialised by tighten2. A code of -1 indicates this row should be
285  processed; -2 says ignore. Other codes should not happen.
286
287  TODO: This assert, and the disabled continue, go back to the change in
288        tighten2 where rows with infinite bounds are handed arbitrary
289        bounds of 1e60 and labelled suitable for further processing. The
290        assert should be removed and the continue clause reinstated.
291        -- lh, 101127 --
292
293  TODO: Replace the original assert (!= -2) with a more effective one.
294        But don't remove the rest of the code structure. I still think I'll
295        end up reintroducing the `ignore' code.  -- lh, 101203 --
296*/
297#   if CGL_DEBUG > 0
298    if (markR[irow] < -1 || markR[irow]%10000 >= nstackR) {
299      std::cout
300        << "Row " << irow << " marked " << markR[irow]
301        << "; expected between -1 and " << nstackR << "." << std::endl ;
302     
303      dump_row(irow,rowLower[irow],rowUpper[irow],minR[irow],maxR[irow],
304               0,false,false,0,0,0,0,0,0,0) ;
305      std::cout << std::endl ;
306    }
307    assert (markR[irow] >= -1 && markR[irow]%10000 < nstackR) ;
308#   endif
309    if (markR[irow] == -1) {
310      stackR[nstackR] = irow ;
311      markR[irow] = nstackR ;
312      saveMin[nstackR] = minR[irow] ;
313      saveMax[nstackR] = maxR[irow] ;
314      nstackR++ ;
315#   if 0
316    } else if (markR[irow] == -2) {
317      continue ;
318#   endif
319    }
320/*
321  LOU_DEBUG: clear count as row bound will change.
322*/
323    markR[irow] = markR[irow]%10000 ;
324/*
325  Case analysis to tighten the appropriate bound.
326*/
327    if (movement > 0.0) {
328      if (value > 0.0) {
329        if (minR[irow] > -1.0e10)
330          minR[irow] += delta ;
331        if (minR[irow] > rowUpper[irow]+1.0e-5) {
332          feasible = false ;
333          break ;
334        }
335      } else {
336        if (maxR[irow] < 1.0e10)
337          maxR[irow] += delta ;
338        if (maxR[irow] < rowLower[irow]-1.0e-5) {
339          feasible = false ;
340          break ;
341        }
342      }
343    } else {
344      if (value > 0.0) {
345        if (maxR[irow] < 1.0e10)
346          maxR[irow] += delta ;
347        if (maxR[irow] < rowLower[irow]-1.0e-5) {
348          feasible = false ;
349          break ;
350        }
351      } else {
352        if (minR[irow] > -1.0e10)
353          minR[irow] += delta ;
354        if (minR[irow] > rowUpper[irow]+1.0e-5) {
355          feasible = false ;
356          break ;
357        }
358      }
359    }
360  }
361  return (feasible) ;
362}
363
364
365/*
366  jjf: clean up djs and solution
367
368  Calculate domains (u<j>-l<j>) for variables, and potential swings in rows. Do
369  some cleanup of reduced costs.
370
371  If CGL_DEBUG > 0, we'll also do some matrix consistency checks.
372*/
373void groomSoln (double direction, double primalTolerance, double *const djs,
374                const double *const colLower, double *const colsol,
375                const double *const colUpper, double *const columnGap,
376                const CoinPackedMatrix *const rowCopy,
377                const int *const rowStartPos,
378                double *const largestNegativeInRow,
379                double *const largestPositiveInRow)
380{
381  int nRows = rowCopy->getNumRows() ;
382  int nCols = rowCopy->getNumCols() ;
383  const int *rowStart = rowCopy->getVectorStarts() ;
384  const int *column = rowCopy->getIndices() ;
385  const double *rowElements = rowCopy->getElements() ;
386
387# ifndef NDEBUG
388  // Note that rowLength is only used in the assert
389  const int *rowLength = rowCopy->getVectorLengths() ;
390# endif
391
392/*
393  Find the largest potential postive and negative swing in a row, where
394  swing is defined as a<ij>(u<j>-l<j>). Check correctness of rowStartPos if
395  we're debugging.
396*/
397  for (int i = 0 ; i < nRows ; i++) {
398
399    assert (rowStart[i]+rowLength[i] == rowStart[i+1]) ;
400
401    int kk ;
402    double value ;
403
404#if CGL_DEBUG > 0
405    // Check correctness of rowStartPos
406    for (kk = rowStart[i] ; kk < rowStart[i+1] ; kk++) {
407      value = rowElements[kk] ;
408      if (value > 0.0) 
409        break ;
410    }
411    assert (rowStartPos[i] == kk) ;
412#endif
413
414    value = 0.0 ;
415    for (kk = rowStart[i] ; kk < rowStartPos[i] ; kk++) {
416      int iColumn = column[kk] ;
417      double gap = CoinMin(1.0e100,colUpper[iColumn]-colLower[iColumn]) ;
418      value = CoinMin(value,gap*rowElements[kk]) ;
419    }
420    largestNegativeInRow[i] = value ;
421    value = 0.0 ;
422    for ( ; kk < rowStart[i+1] ; kk++) {
423      int iColumn = column[kk] ;
424      double gap = CoinMin(1.0e100,colUpper[iColumn]-colLower[iColumn]) ;
425      value = CoinMax(value,gap*rowElements[kk]) ;
426    }
427    largestPositiveInRow[i] = value ;
428  }
429/*
430  Convert reduced costs to minimisation and clean up any that are on the wrong
431  side of zero.
432
433  TODO: Seems like we could move this ahead of the previous loops and then
434        use columnGap[i] to calculate largest[Negative,Positive]InRow.
435        -- lh, 101125 --
436*/
437  for (int i = 0 ; i < nCols ; ++i) {
438    double djValue = djs[i]*direction ;
439    double gap = colUpper[i]-colLower[i] ;
440    if (gap > 1.0e-8) {
441      if (colsol[i] < colLower[i]+primalTolerance) {
442        colsol[i] = colLower[i] ;
443        djs[i] = CoinMax(0.0,djValue) ;
444      } else if (colsol[i] > colUpper[i]-primalTolerance) {
445        colsol[i] = colUpper[i] ;
446        djs[i] = CoinMin(0.0,djValue) ;
447      } else {
448        djs[i] = 0.0 ;
449      }
450    }
451    columnGap[i] = gap-primalTolerance ;
452  }
453
454  return ;
455}
456
457/*
458  Given that the variable being probed has been discovered to be monotone,
459  this method will save the resulting bound changes for integer variables as
460  column cuts.
461
462  One could do the same for continuous variables, but at present we don't.
463  Presumably we need to test against origColLower, origColUpper (the
464  bounds held in the si) to decide if we have real cuts on the original
465  problem. One could argue that we could test against saveL and saveU and
466  save the trouble of repeatedly recording the same column cut.
467
468  Index and elements are scratch arrays.
469
470  Returns true if any of the bound improvements are actually cuts, false
471  otherwise.
472*/
473bool monotoneActions (double primalTolerance_,
474                      const OsiSolverInterface &si,
475                      OsiCuts &cs,
476                      int nstackC, const int *const stackC,
477                      const char *const intVar,
478                      const double *const colLower,
479                      const double *const colsol,
480                      const double *const colUpper,
481                      int *const index, double *const element)
482{
483  OsiColCut cc ;
484  bool anyColumnCuts = false ;
485  int nTot = 0 ;
486  int nFix = 0 ;
487  const double *origColLower = si.getColLower() ;
488  const double *origColUpper = si.getColUpper() ;
489# if CGL_DEBUG > 0
490  int probedVar = stackC[0] ;
491  std::cout
492    << "    " << si.getColName(probedVar) << " (" << probedVar
493    << ") monotone [" << colLower[probedVar] << ", "
494    << colUpper[probedVar] << "] from " << colsol[probedVar]
495    << "." << std::endl ;
496# endif
497/*
498  Stash the improved column bounds, u<j> first, then l<j>. If any of them
499  actually cut off the current solution, indicate that we have cuts.
500*/
501  bool ifCut = false ;
502  for (int istackC = 0 ; istackC < nstackC ; istackC++) {
503    int icol = stackC[istackC] ;
504    if (intVar[icol]) {
505      if (colUpper[icol] < origColUpper[icol]-1.0e-4) {
506        element[nFix] = colUpper[icol] ;
507        index[nFix++] = icol ;
508        if (colsol[icol] > colUpper[icol]+primalTolerance_) {
509          ifCut = true ;
510          anyColumnCuts = true ;
511        }
512      }
513    }
514  }
515  if (nFix) {
516    nTot = nFix ;
517    cc.setUbs(nFix,index,element) ;
518    nFix = 0 ;
519  }
520  for (int istackC = 0 ; istackC < nstackC ; istackC++) {
521    int icol = stackC[istackC] ;
522    if (intVar[icol]) {
523      if (colLower[icol] > origColLower[icol]+1.0e-4) {
524        element[nFix] = colLower[icol] ;
525        index[nFix++] = icol ;
526        if (colsol[icol] < colLower[icol]-primalTolerance_) {
527          ifCut = true ;
528          anyColumnCuts = true ;
529        }
530      }
531    }
532  }
533  if (nFix) {
534    nTot += nFix ;
535    cc.setLbs(nFix,index,element) ;
536  }
537  if (nTot) {
538    if (ifCut) {
539      cc.setEffectiveness(100.0) ;
540    } else {
541      cc.setEffectiveness(1.0e-5) ;
542    }
543#if CGL_DEBUG > 0
544    checkBounds(si,cc) ;
545#endif
546    cs.insert(cc) ;
547  }
548
549  return (anyColumnCuts) ;
550}
551
552void clearStacks (double primalTolerance_,
553                  int nstackC, int *const stackC, int *const markC,
554                  const double *const colLower, const double *const colUpper,
555                  int nstackR, int *const stackR, int *const markR)
556{
557/*
558  Clear off the bound change markers, except for fixed variables. We'll be
559  moving on to another variable.
560*/
561  for (int istackC = 0 ; istackC < nstackC ; istackC++) {
562    int icol = stackC[istackC] ;
563    if (colUpper[icol]-colLower[icol] > primalTolerance_) {
564      markC[icol] &= ~(tightenLower|tightenUpper) ;
565    } else {
566      markC[icol] = tightenLower|tightenUpper ;
567    }
568  }
569  for (int istackR = 0 ; istackR < nstackR ; istackR++) {
570    int irow = stackR[istackR] ;
571    markR[irow] = -1 ;
572  }
573}
574
575
576
577/*
578  This method examines the rows on stackR looking for redundant rows that
579  consist solely of singleton variables (i.e., variables that appear in just
580  this constraint). It then looks to see if any of the variables in the row
581  can be fixed at bound.
582
583  Consider the case where U<i> < b<i> and all unfixed variables in the row
584  are singletons (occur in no other constraint). Given x<j> is a singleton,
585  the reduced cost is c<j> - ya<j> = c<j> - ye<i> = c<j> - y<i>.  But if
586  U<i> < b<i>, the constraint is loose by definition, hence y<i> = 0 and
587  the reduced cost is c<j>. If a<ij> > 0, x<j> should be u<j> in U<i>. If
588  c<j> < 0, minimising will tend to drive x<j> to u<j>. This cannot violate
589  constraint i (U<i> < b<i>) and (since x<j> is a singleton) will have no
590  effect elsewhere. Hence we can fix x<j> at u<j>.
591
592  Do the case analysis, and what you find is that against U<i> we want
593    a<ij> > 0  ==>  c<j> < 0  ==>  push x<j> to u<j>
594    a<ij> < 0  ==>  c<j> > 0  ==>  push x<j> to l<j>
595  and against L<i> we want
596    a<ij> > 0  ==>  c<j> > 0  ==>  push x<j> to l<j>
597    a<ij> < 0  ==>  c<j> < 0  ==>  push x<j> to u<j>
598
599  Extend it one more time by taking the objective direction (max/min) into
600  account (dir = 1.0 for min, -1.0 for max) and you have
601    against U<i> ==> a<ij>*c<j>*dir < 0
602    against L<i> ==> a<ij>*c<j>*dir > 0
603
604  John's original comment for this code was
605    // also see if singletons can go to good objective
606    // Taken out as should be found elsewhere
607    // and has to be original column length
608  but he reinstated it. Apparently there were cases where fixing the probing
609  variable was required to satisfy the condition that all unfixed variables be
610  singletons. Enough of them to justify reinstating this code.
611*/
612
613void singletonRows (int jProbe, double primalTolerance_,
614                    const OsiSolverInterface &si,
615                    const CoinPackedMatrix *rowCopy,
616                    int *const markC,
617                    int &nstackC, int *const stackC,
618                    double *const saveL, double *const saveU,
619                    double *const colUpper, double *const colLower,
620                    double *const colGap,
621                    const int nstackR, const int *const stackR,
622                    const double *const rowUpper,
623                    const double *const rowLower,
624                    const double *const maxR, const double *const minR)
625{
626/*
627  Unpack a few vectors from the row-major matrix.
628*/
629  const CoinBigIndex *rowStart = rowCopy->getVectorStarts() ;
630  const int *column = rowCopy->getIndices() ;
631  const double *rowElements = rowCopy->getElements() ;
632  const int nCols = rowCopy->getNumCols() ;
633/*
634  `Singleton' must be based on the column length in the original system.
635*/
636  const double *objective = si.getObjCoefficients() ;
637  const int *columnLengths = si.getMatrixByCol()->getVectorLengths() ;
638  const double objSense = si.getObjSense() ;
639/*
640  Open a loop to work through the rows on stackR.
641*/
642  for (int istackR = 0 ; istackR < nstackR ; istackR++) {
643    int i = stackR[istackR] ;
644/*
645  Check the gaps. If the constraint is potentially tight in both directions,
646  there's nothing more to do here.
647*/
648    const double uGap = rowUpper[i]-maxR[i] ;
649    const double lGap = minR[i]-rowLower[i] ;
650    if (uGap < primalTolerance_ && lGap < primalTolerance_) continue ;
651/*
652  Scan the row to see if it meets the `all singletons' condition. Again,
653  if this fails, there's nothing more to be done.
654
655  Note that the original code didn't check the probing variable x<p>,
656  because if you're probing a binary variable it's fixed. But for general
657  integer variables a probe does not in general fix the variable.  So we
658  check all variables.
659
660  We should not be executing this method if we have prima facie infeasibility.
661*/
662    bool canFix = true ;
663    for (int kk = rowStart[i] ; kk < rowStart[i+1] ; kk++) {
664      int j = column[kk] ;
665      assert(colUpper[j]-colLower[j] > -primalTolerance_) ;
666      if (colUpper[j] > colLower[j] && columnLengths[j] != 1) {
667        canFix = false ;
668        break ;
669      }
670    }
671    if (!canFix) continue ;
672/*
673  If we've passed the tests, we can look for variables suitable to drive to
674  bound. Work against U<i> first. We're looking for variables with a<ij> > 0
675  that will be naturally driven to u<j>, or variables with a<ij> < 0 that will
676  be naturally driven to l<j>.
677 
678  Don't overwrite the saved bounds if we've
679  tightened this variable already!
680*/
681    if (uGap > primalTolerance_) {
682      for (int kk = rowStart[i] ; kk < rowStart[i+1] ; kk++) {
683        int j = column[kk] ;
684        const double lj = colLower[j] ;
685        const double uj = colUpper[j] ;
686        if (uj > lj) {
687          double value = rowElements[kk] ;
688          if (objSense*objective[j]*value < 0.0)
689          { assert(jProbe != j) ;
690            if (!(markC[j]&(tightenLower|tightenUpper))) {
691              stackC[nstackC] = j ;
692              saveL[nstackC] = lj ;
693              saveU[nstackC] = uj ;
694            }
695            assert(saveU[nstackC] > saveL[nstackC]) ;
696            if (value > 0.0) {
697              colLower[j] = uj ;
698            } else {
699              colUpper[j] = lj ;
700            }
701            markC[j] |= tightenLower|tightenUpper ;
702            colGap[j] = -primalTolerance_ ;
703            assert(nstackC < nCols) ;
704            nstackC++ ;
705          }
706        }
707      }
708    }
709/*
710  And now the identical code, except that we're working against L<i>, hence
711  the sense of the elibigility test is reversed and we want variables with
712  a<ij> > 0 that will be naturally driven to l<j>, or variables with
713  a<ij> < 0 that will be naturally driven to u<j>.
714*/
715    if (lGap > primalTolerance_) {
716      for (int kk = rowStart[i] ; kk < rowStart[i+1] ; kk++) {
717        int j = column[kk] ;
718        const double lj = colLower[j] ;
719        const double uj = colUpper[j] ;
720        if (uj > lj) {
721          double value = rowElements[kk] ;
722          if (objSense*objective[j]*value > 0.0)
723          { assert(jProbe != j) ;
724            if (!(markC[j]&(tightenLower|tightenUpper))) {
725              stackC[nstackC] = j ;
726              saveL[nstackC] = lj ;
727              saveU[nstackC] = uj ;
728            }
729            assert(saveU[nstackC] > saveL[nstackC]) ;
730            if (value < 0.0) {
731              colLower[j] = uj ;
732            } else {
733              colUpper[j] = lj ;
734            }
735            markC[j] |= tightenLower|tightenUpper ;
736            colGap[j] = -primalTolerance_ ;
737            assert(nstackC < nCols) ;
738            nstackC++ ;
739          }
740        }
741      }
742    }
743  }
744  return ;
745}
746
747
748/*
749  Generate coefficient strengthening cuts.
750 
751  Assume we have a binary probe variable x<p> with a<ip> > 0. This will
752  be a down probe, so assume that U<i> > b<i> before the probe, but now
753  U'<i> < b<i> (i.e., the constraint is redundant against b<i> after the
754  down probe forces x<p> to 0 and reduces U<i> to U'<i>). We would like to
755  have the following:
756    * When x<p> = 0, b'<i> = U'<i>  (after all, the lhs can't be any larger)
757    * When x<p> = 1, b'<i> = b<i>   (the original value)
758  Define delta<i> = b<i> - U'<i>. Let b'<i> = b<i>-delta<i> and let
759  a'<ip> = a<ip>-delta<i>. Then when x<p> = 1, the delta<i> on each side
760  cancels and we're left with the original constraint. When x<p> = 0, the rhs
761  becomes b'<i> = b<i>+delta<i> = U'<i>.
762
763  The usual case analysis applies; see the typeset documentation. It's not
764  necessarily true that U'<i> = U<i> - a<ip>u<p>; additional propagation could
765  have reduced it even further. That doesn't alter the reasoning above.
766
767  TODO: Figure out why goingToTrueBound is relevant here, other than it's
768        set to zero if we're infeasible.  -- lh, 101127 --
769
770        I think that might be all there is to it.   -- lh, 101128 --
771
772        Well, maybe not. Check the math for the contained cuts to see if
773        it's valid only for binary variables; when within distance 1 of
774        the upper or lower bound of a general integer; or more generally.
775        GoingToTrueBound indicates the first of these two. -- lh, 110105 --
776
777        After working the math a few times, it looks to me like the important
778        aspect is that goingToTrueBound = 2 indicates binary variables. I
779        can't get the math to agree with the code otherwise. See the problem
780        noted below. I'm going to restrict this method to binary variables
781        until I can resolve questions about general integers.  -- lh, 110113 --
782*/
783
784#define STRENGTHEN_PRINT
785void strengthenCoeff (
786                      int jProbe,
787                      double primalTolerance_,
788                      bool justReplace, bool canReplace,
789                      double needEffectiveness,
790                      const OsiSolverInterface &si,
791                      CglProbingRowCut &rowCut,
792                      const CoinPackedMatrix *rowCopy,
793                      double *const colUpper, double *const colLower,
794                      const double *const colsol,
795                      const int nstackR, const int *const stackR,
796                      const double *const rowUpper,
797                      const double *const rowLower,
798                      const double *const maxR,
799                      const double *const minR,
800                      const int *const realRows,
801                      double *const element,
802                      int *const index,
803                      CglTreeInfo *const info
804                     )
805{
806# if CGL_DEBUG > 0
807  std::cout
808    << "Entering strengthenCoeff, probed variable "
809    << si.getColName(jProbe) << "(" << jProbe << ")." << std::endl ;
810  const OsiRowCutDebugger *debugger = si.getRowCutDebugger() ;
811# endif
812/*
813  Unpack a few vectors from the row-major matrix.
814*/
815  const CoinBigIndex *rowStart = rowCopy->getVectorStarts() ;
816  const int *column = rowCopy->getIndices() ;
817  const double *rowElements = rowCopy->getElements() ;
818/*
819  Open up an outer loop to walk stackR and look for interesting constraints.
820*/
821  for (int istackR = 0 ; istackR < nstackR ; istackR++) {
822    int irow = stackR[istackR] ;
823/*
824  We can't get anywhere unless probing has made one end or the other of the
825  constraint redundant.
826
827  Really, we shouldn't be here if b<i> or blow<i> are not finite. Check with
828  an assert.
829*/
830    double uGap = rowUpper[irow]-maxR[irow] ;
831    assert(uGap < 1.0e8) ;
832    double lGap = minR[irow]-rowLower[irow] ;
833    assert(lGap < 1.0e8) ;
834    if (uGap < primalTolerance_ && lGap < primalTolerance_) continue ;
835/*
836  We'll need the lhs activity to evaluate the effectiveness of the cut. Do a
837  base calculation which we'll correct later.
838 
839  TODO: Sort out why anyColumnCuts is an obstacle. Theoretical or practical?
840        Could we use the new bound instead?  -- lh, 101128 --
841
842        After working the math, looks to me like the only effect will be where
843        a column cut cuts off the current solution, in which case the lhs
844        value (sum) may be incorrect, affecting the calculation of
845        effectiveness. We could correct for this, if it was worth the effort
846        (a max or min when calculating the sum below). -- lh, 110113 --
847*/
848    double sum = 0.0 ;
849    for (int kk = rowStart[irow] ; kk < rowStart[irow+1] ; kk++) {
850      sum += rowElements[kk]*colsol[column[kk]] ;
851    }
852/*
853  Start out by working against the row upper bound U(i).  We've calculated
854  sum using the current a<ik>, so reduce it by (-a<ip>+a'<ip>)x*<p> =
855  (-uGap)(x*<p>) to do the check.
856
857  We're not willing to generate a cut if it doesn't cut off the current
858  solution, but we are willing to strengthen the existing constraint in place.
859  Can't hurt, eh? Excluding range constraints in this case avoids the issue
860  of conflicting a'<ij> if it turns out we can strengthen against b<i>
861  and blow<i> in the same constraint.
862*/
863    if (uGap > primalTolerance_ &&
864        (sum-uGap*colsol[jProbe] > maxR[irow]+primalTolerance_ ||
865         (info->strengthenRow && rowLower[irow] < -1.0e20))) {
866/*
867  Generate the coefficients. For all except the probing variable, we just copy
868  the coefficient. The probing variable becomes a'<ij> = (a<ij> - uGap). Allow
869  for the possibility that a<ij> starts or ends as 0.
870
871  TODO: The value of sum2 calculated here should be precisely
872        sum-uGap*colsol[j], unless my math is way wrong.  -- lh, 110107 --
873*/
874      int n = 0 ;
875      bool saw_aip = false ;
876      double sum2 = 0.0 ;
877      for (int kk = rowStart[irow] ; kk < rowStart[irow+1] ; kk++) {
878        int kColumn = column[kk] ;
879        double el = rowElements[kk] ;
880        if (kColumn != jProbe) {
881          index[n] = kColumn ;
882          element[n++] = el ;
883        } else {
884          el = el-uGap ;
885          if (fabs(el) > 1.0e-12) {
886            index[n] = kColumn ;
887            element[n++] = el ;
888          }
889          saw_aip = true ;
890        }
891        sum2 += colsol[kColumn]*el ;
892      }
893      if (!saw_aip) {
894        index[n] = jProbe ;
895        element[n++] = -uGap ;
896        sum2 -= colsol[jProbe]*uGap ;
897      }
898      assert(sum == sum2) ;
899/*
900  Check effectiveness. Add the cut only if it's sufficiently effective. (If
901  we're strengthening existing constraints, needEffectiveness is quite low.)
902
903  TODO: (*) I can't justify colLower in uGap*(l<p>+1.0). When I work
904        the math, b'<i> = b<i> - (b<i> - U'<i>) = b<i> - uGap. I'll keep
905        trying, on the theory that it's some attempt at general integers, but
906        I'm thinking it's wrong. Note that for binary variables, l<p> = 0 and
907        there's no problem.   -- lh, 110113 --
908
909  TODO: I don't understand what's going on with the first two calculations
910        for effectiveness. First off, sum2 and sum-uGap*colsol[j] should be
911        exactly the same value (see notes with previous calculations). Then
912        there's the question `Why divide by (b<i>-U<i>)?' in the second
913        calculation. This inflates the effectiveness for a tiny uGap. A
914        normalisation, perhaps?  -- lh, 110107 --
915
916  TODO: And why go through the effort of setting the row bounds and
917        effectiveness in the OsiRowCut, before we've decided that the cut will
918        be used?  -- lh, 110107 --
919*/
920      OsiRowCut rc ;
921      rc.setLb(-DBL_MAX) ;
922      // (*) double ub = rowUpper[irow]-uGap*(colLower[jProbe]+1.0) ;
923      double ub = rowUpper[irow]-uGap ;
924      rc.setUb(ub) ;
925      double effectiveness = sum2-ub ;
926      effectiveness =
927          CoinMax(effectiveness,
928                  (sum-uGap*colsol[jProbe]-maxR[irow])/uGap) ;
929      if (!saw_aip)
930        effectiveness = CoinMax(1.0e-7,effectiveness) ;
931      rc.setEffectiveness(effectiveness) ;
932      // rc.setEffectiveness((sum-uGap*colsol[jProbe]-maxR[irow])/uGap) ;
933      if (rc.effectiveness() > needEffectiveness) {
934        rc.setRow(n,index,element,false) ;
935#       if CGL_DEBUG > 0
936        if (debugger) assert(!debugger->invalidCut(rc)); 
937#       endif
938#       ifdef STRENGTHEN_PRINT
939        if (canReplace && rowLower[irow] < -1.0e20) {
940          printf("Strengthen Cut (1):\n") ;
941          dump_row(irow,rc.lb(),rc.ub(),
942                   nan(""),nan(""),&si,true,false,4,
943                   n,index,element,
944                   1.0e-10,colLower,colUpper) ;
945          printf("Original Row:\n") ;
946          int k = rowStart[irow] ;
947          dump_row(irow,rowLower[irow],rowUpper[irow],
948                   minR[irow],maxR[irow],&si,true,false,4,
949                   rowStart[irow+1]-k,&column[k],&rowElements[k],
950                   1.0e-10,colLower,colUpper) ;
951        }
952#       endif
953/*
954  If we're in preprocessing, we might try to simply replace the existing
955  constraint (justReplace = true; preprocessing is the only place this will
956  happen). Otherwise, drop the cut into the cut set.
957
958  realRows comes in as a parameter. This is the translation array created if
959  we modified the constraint system during preprocessing in gutsOfGenerateCuts.
960
961  TODO: Seems a bit disingenuous to fail at replacement now, given that
962        effectiveness is artificially low. Notice again the inconsistent use
963        of finite infinities on the row lb.  -- lh, 101128 --
964
965  TODO: To elaborate on the above, it seems to me that we can get here with
966        justReplace = true and canReplace = false, and this condition is
967        constant over an iteration of the way loop. In other words, we've done
968        all the work to generate this cut and we knew before we started that
969        we would discard it here.  -- lh, 110107 --
970*/
971        int realRow = (canReplace && rowLower[irow] < -1.0e20)?(irow):(-1) ;
972        if (realRows && realRow >= 0)
973          realRow = realRows[realRow] ;
974        if (!justReplace) {
975          rowCut.addCutIfNotDuplicate(rc,realRow) ;
976        } else if (realRow >= 0) {
977          double effectiveness = 0.0 ;
978          for (int i = 0 ; i < n ; i++)
979            effectiveness += fabs(element[i]) ;
980          if (!info->strengthenRow[realRow] ||
981              info->strengthenRow[realRow]->effectiveness() > effectiveness) {
982            delete info->strengthenRow[realRow] ;
983            rc.setEffectiveness(effectiveness) ;
984            info->strengthenRow[realRow]=rc.clone() ;
985          }
986        }
987      }
988    }
989/*
990  And repeat working against the lower bound L(i). As in so many other places
991  in this code, it's the identical functionality, with some subtle edits that
992  distinguish the L(i) math from the U(i) math.
993*/
994    if (lGap > primalTolerance_ &&
995        (sum+lGap*colsol[jProbe] < minR[irow]-primalTolerance_ ||
996         (info->strengthenRow && rowUpper[irow] > -1.0e20))) {
997      int n = 0 ;
998      bool saw_aip = false ;
999      double sum2 = 0.0 ;
1000      for (int kk = rowStart[irow] ; kk < rowStart[irow+1] ; kk++) {
1001        int kColumn = column[kk] ;
1002        double el = rowElements[kk] ;
1003        if (kColumn != jProbe) {
1004          index[n] = kColumn ;
1005          element[n++] = el ;
1006        } else {
1007          el = el+lGap ;
1008          if (fabs(el) > 1.0e-12) {
1009            index[n] = kColumn ;
1010            element[n++] = el ;
1011          }
1012          saw_aip = true ;
1013        }
1014        sum2 += colsol[kColumn]*el ;
1015      }
1016      if (!saw_aip) {
1017        index[n] = jProbe ;
1018        element[n++] = lGap ;
1019        sum2 += colsol[jProbe]*lGap ;
1020      }
1021      OsiRowCut rc ;
1022      double lb = rowLower[irow]+lGap*(colLower[jProbe]+1.0) ;
1023      rc.setLb(lb) ;
1024      rc.setUb(DBL_MAX) ;
1025      // effectiveness
1026      double effectiveness = lb-sum2 ;
1027      effectiveness = CoinMax(effectiveness,
1028                              (minR[irow]-sum-lGap*colsol[jProbe])/lGap) ;
1029      if (!saw_aip)
1030        effectiveness = CoinMax(1.0e-7,effectiveness) ;
1031      rc.setEffectiveness(effectiveness) ;
1032      if (rc.effectiveness() > needEffectiveness) {
1033        rc.setRow(n,index,element,false) ;
1034#       if CGL_DEBUG > 0
1035        if (debugger) assert(!debugger->invalidCut(rc)); 
1036#       endif
1037#       ifdef STRENGTHEN_PRINT
1038        if (canReplace && rowUpper[irow] > 1.0e20) {
1039          printf("Strengthen Cut (2):\n") ;
1040          dump_row(irow,rc.lb(),rc.ub(),
1041                   nan(""),nan(""),&si,true,false,4,
1042                   n,index,element,
1043                   1.0e-10,colLower,colUpper) ;
1044          printf("Original Row:\n") ;
1045          int k = rowStart[irow] ;
1046          dump_row(irow,rowLower[irow],rowUpper[irow],
1047                   minR[irow],maxR[irow],&si,true,false,4,
1048                   rowStart[irow+1]-k,&column[k],&rowElements[k],
1049                   1.0e-10,colLower,colUpper) ;
1050        }
1051#       endif
1052
1053        int realRow = (canReplace&&rowUpper[irow]>1.0e20) ? irow : -1 ;
1054        if (realRows && realRow >= 0)
1055          realRow = realRows[realRow] ;
1056        if (!justReplace) {
1057          rowCut.addCutIfNotDuplicate(rc,realRow) ;
1058        } else if (realRow >= 0) {
1059          double effectiveness = 0.0 ;
1060          for (int i = 0 ; i < n ; i++)
1061            effectiveness += fabs(element[i]) ;
1062          if (!info->strengthenRow[realRow] ||
1063              info->strengthenRow[realRow]->effectiveness() > effectiveness) {
1064            delete info->strengthenRow[realRow] ;
1065            rc.setEffectiveness(effectiveness) ;
1066            info->strengthenRow[realRow] = rc.clone() ;
1067          }
1068        }
1069      }
1070    }
1071  }
1072
1073# if CGL_DEBUG > 0
1074  std::cout
1075    << "Leaving strengthenCoeff, "
1076    << rowCut.numberCuts() << " cuts." << std::endl ;
1077# endif
1078  return ;
1079}
1080
1081
1082
1083
1084// =========================================================
1085
1086
1087
1088/*
1089  jjf: Does probing and adding cuts
1090
1091  Note that this method is heavily commented and instrumented in CbcAnn.
1092
1093  It would appear that probe() has received more attention that probeClique or
1094  probeSlack. Neither of them has the ONE_ARRAY option.
1095*/
1096
1097/*
1098  We're going to probe integer variables, and we're interested in propagating
1099  the effect of each probe (force a variable up or down).
1100
1101  The bulk of the method consists of three nested loops:
1102    * PASS LOOP: makes multiple passes, probing all variables in lookedAt_
1103    * LOOKEDAT LOOP: iterates over the variables in lookedAt_
1104    * PROBE LOOP: probes down/up/down for each variable.
1105
1106  The body of the probe loop runs a bit over 2600 lines and propagates
1107  the effect of forcing the probed variable down or up.
1108
1109    * The start of the body updates the row bounds of affected rows, then
1110      initialises stackC with the probed variable.
1111
1112    * The middle is a 1,000 line loop that does the bulk of propagation.
1113      At the start there's some dead code (PROBING100) associated with
1114      CglTreeProbingInfo. The purpose remains to be determined. This is
1115      followed by six specialised replications of the same functionality:
1116      walk the column of a variable from stackC, and for each affected
1117      row, try to tighten the bounds of other variables in the row. If
1118      we successfully tighten bounds on a variable, walk the column of
1119      that variable, tighten the bounds of affected rows, and place the
1120      variable on stackC.
1121
1122    * The end of the body deals with the result of the probe. At the end
1123      of each iteration, there's a decision: have we proven infeasibility
1124      (INFEASIBILITY) or monotonicity (MONOTONE), or do we need another
1125      iteration (ITERATE). Monotonicity and iteration are large blocks
1126      in themselves.
1127
1128  There was a fair bit of preamble and postamble around the PASS loop.
1129  The preamble code is protected by PROBING_EXTRA_STUFF and has been disabled
1130  since 2006. It appears to find disaggregation cuts.  The postamble code
1131  purports to find bigM constraints. It's clearly a work in progress. I've
1132  moved both out to tmp files.  -- lh, 101203 --
1133
1134  Comments on the data structures:
1135
1136  markC  0: not fixed
1137         1: tightened upper bound
1138         2: tightened lower bound
1139         3: fixed
1140
1141    or perhaps the above should be interpreted as bit codes, along with the
1142    bits used to indicate infinite bounds:
1143
1144       0x0: no bounds tightened
1145       0x1: u<j> tightened
1146       0x2: l<j> tightened
1147       0x4: l<j> < -1e10 (-infty)
1148       0x8: u<j> >  1e10 (+infty)
1149
1150    Note that we'll only tighten a bound once --- this is enforced during
1151    bounds propagation.
1152
1153  stackC      Records columns to be processed. This record is started anew
1154              each time we probe, pushing a variable in a given direction.
1155              The variable being probed is always entry 0 on stackC.
1156              saveL and saveU are correlated with stackC.
1157
1158              This isn't a stack --- it's a record of every variable
1159              processed. nstackC is the tail of the queue, and istackC is
1160              the head. At the end of a probing pass it's scanned to recover
1161              the actual implications that were generated. It's not clear
1162              to me why it's allocated at 2*ncols. Perhaps, at one time,
1163              there were separate entries for upper and lower bound changes?
1164
1165              No, just as I'm thinking I'm almost done, I notice there's some
1166              fancy footwork in the portion of the code that handles the case
1167              where both the up and down probes were feasible. Seems likely
1168              that we're allowing for each direction to process all variables
1169              at most once.
1170
1171              And as I'm working through the code that handles the case
1172              where down and up probes were both feasible, there in the
1173              move singleton code is an assert that nstackC < nCols.
1174
1175              And by the time I finally reach the end, I'm pretty well
1176              convinced that it's an error if we exceed nCols on stackC.
1177              The only thing to worry about is whether singleton processing
1178              could add a variable that's already been added due to
1179              propagation, and a check of the code says the answer is 'no'.
1180
1181       TODO:  When I rewrite this, allocate stackC as nCols and put in some
1182              asserts to catch any violations. There are already asserts
1183              in place at the point where a variable is placed on stackC
1184              during propagation. (Look for `Is there any change to
1185              propagate'.)   -- lh, 101128 --
1186
1187
1188  stackR      Stack rows to be processed. This stack is started anew each
1189              time we probe pushing a variable in a given direction.
1190
1191  minR, maxR, markR  Initialised externally by tighten2. For row i, minR[i] =
1192                     LB(i), maxR[i] = UB(i) (row lhs lower and upper bounds,
1193                     respectively).  markR is set to -1 if there's at least
1194                     one finite lhs bound, -2 if we shouldn't process this
1195                     row (no finite bounds, or too many coefficients).
1196
1197  Then, as we work, markR[i] will be set to point to the entry in stackR for
1198  row i. saveMin and saveMax are correlated with stackR, and hold the
1199  original values for LB(i) and UB(i) (from minR and maxR) while we're
1200  working. As it turns out, however, the the back pointer from markR is never
1201  used.
1202
1203  largestPositiveInRow  for row i, max{j in P} a<ij>(u<i>-l<j>), where P is
1204                        the set of positive coefficients a<ij>
1205  largestNegativeInRow  for row i, min{j in M} a<ij>(u<i>-l<j>), where M is
1206                        the set of negative coefficients a<ij>
1207
1208  element and index are scratch arrays used to construct cuts.
1209
1210  realRows is the translation array created if we modified the constraint
1211  system during preprocessing in gutsOfGenerateCuts.
1212
1213*/
1214int CglProbing::probe( const OsiSolverInterface &si, 
1215                       OsiCuts &cs, 
1216                       double *colLower, double *colUpper, 
1217                       CoinPackedMatrix *rowCopy,
1218                       CoinPackedMatrix *columnCopy,
1219                       const CoinBigIndex *rowStartPos, const int *realRows, 
1220                       const double *rowLower, const double *rowUpper,
1221                       const char *intVar, double *minR, double *maxR, 
1222                       int *markR, 
1223                       CglTreeInfo *info) const
1224
1225{
1226# if CGL_DEBUG > 0
1227  std::cout << "Entering CglProbing::probe." << std::endl ;
1228# endif
1229/*
1230  PREPARATION
1231
1232  All the code down through `PASS LOOP: HEAD' is preparation. Do all of the
1233  setup that will not change over the nested loops that do the work of probing.
1234  Note that rowCopy may have considerably fewer rows than the original system.
1235*/
1236
1237  int nRows = rowCopy->getNumRows() ;
1238  int nCols = rowCopy->getNumCols() ;
1239  int maxStack = info->inTree ? maxStack_ : maxStackRoot_ ;
1240
1241/*
1242  Fiddle maxStack. Odd sort of function:
1243     1 < maxStack <= 25 => 2*maxStack
1244    26 < maxStack <= 50 => 50
1245    51 < maxStack       => maxStack
1246
1247  TODO: Grepping through the code says there's no way that totalTimesCalled_
1248        can become negative. Perhaps in some derived class? More likely a
1249        magic number for something, or else a quick way to disable this
1250        code.  -- lh, 101021 --
1251*/
1252  if ((totalTimesCalled_%10) == -1) {
1253    int newMax = CoinMin(2*maxStack,50) ;
1254    maxStack = CoinMax(newMax,maxStack) ;
1255  }
1256  totalTimesCalled_++ ;
1257
1258# ifdef ONE_ARRAY
1259  assert((DIratio != 0) && ((DIratio&(DIratio-1)) == 0)) ;
1260  int nSpace = 8*nCols+4*nRows+2*maxStack ;
1261  nSpace += (4*nCols+nRows+maxStack+DIratio-1)>>(DIratio-1) ;
1262  double * colsol = new double[nSpace] ;
1263  double * djs = colsol + nCols ;
1264  double * columnGap = djs + nCols ;
1265  double * saveL = columnGap + nCols ;
1266  double * saveU = saveL + 2*nCols ;
1267  double * saveMin = saveU + 2*nCols ;
1268  double * saveMax = saveMin + nRows ;
1269  double * largestPositiveInRow = saveMax + nRows ;
1270  double * largestNegativeInRow = largestPositiveInRow + nRows ;
1271  double * element = largestNegativeInRow + nRows ;
1272  double * lo0 = element + nCols ;
1273  double * up0 = lo0 + maxStack ;
1274  int * markC = reinterpret_cast<int *> (up0+maxStack) ;
1275  int * stackC = markC + nCols ;
1276  int * stackR = stackC + 2*nCols ;
1277  int * index = stackR + nRows ;
1278  int * stackC0 = index + nCols ;
1279# else
1280  double * colsol = new double[nCols] ;
1281  double * djs = new double[nCols] ;
1282  double * columnGap = new double [nCols] ;
1283  double * saveL = new double [2*nCols] ;
1284  double * saveU = new double [2*nCols] ;
1285  double * saveMin = new double [nRows] ;
1286  double * saveMax = new double [nRows] ;
1287  double * largestPositiveInRow = new double [nRows] ;
1288  double * largestNegativeInRow = new double [nRows] ;
1289  double * element = new double[nCols] ;
1290  double * lo0 = new double[maxStack] ;
1291  double * up0 = new double[maxStack] ;
1292  int * markC = new int [nCols] ;
1293  int * stackC = new int [2*nCols] ;
1294  int * stackR = new int [nRows] ;
1295  int * index = new int[nCols] ;
1296  int * stackC0 = new int[maxStack] ;
1297# endif
1298
1299/*
1300  Create a local container to hold the cuts we generate. At the end we'll
1301  transfer these to the container passed as a parameter.
1302
1303  Typically, rowCopy (nRows) will be smaller than the original system in the
1304  si, but it could be one larger (all original constraints plus the objective).
1305  At the root, allow lots of room; even more if this is the first call at the
1306  root. In the tree, much less.
1307
1308  If all we're doing is strengthening rows in place (option 0x40 indicates
1309  we're in preprocessing, a good place to do this), we can only work with
1310  what we have in rowCopy (and we need a translation array).
1311
1312  TODO: There's a problem with realRows here --- if we didn't delete any
1313        rows during preprocessing in gutsOfGenerateCuts, realRows is not
1314        allocated. So we can't strengthen in place in that case. Not likely
1315        the intent.  -- lh, 101209 --
1316*/
1317  bool justReplace = ((info->options&0x40) != 0) && (realRows != NULL) ;
1318  int nRowsFake = 0 ;
1319  if (justReplace) {
1320    nRowsFake = nRows ;
1321  } else {
1322    int nRowsSafe = CoinMin(nRows,si.getNumRows()) ;
1323    nRowsFake = info->inTree ? nRowsSafe/3 : nRowsSafe*10 ;
1324    if (!info->inTree && info->pass == 0) nRowsFake *= 10 ;
1325  }
1326  CglProbingRowCut rowCut(nRowsFake,!info->inTree) ;
1327
1328
1329  // Unpack matrices
1330  const int * column = rowCopy->getIndices() ;
1331  const CoinBigIndex * rowStart = rowCopy->getVectorStarts() ;
1332  const double * rowElements = rowCopy->getElements() ;
1333
1334  const int * row = columnCopy->getIndices() ;
1335  const CoinBigIndex * columnStart = columnCopy->getVectorStarts() ;
1336  const int * columnLength = columnCopy->getVectorLengths(); 
1337  const double * columnElements = columnCopy->getElements() ;
1338
1339/*
1340  Grab the column solution and reduced costs and groom them.
1341*/
1342  CoinMemcpyN(si.getReducedCost(),nCols,djs) ;
1343  CoinMemcpyN(si.getColSolution(),nCols,colsol) ;
1344  double direction = si.getObjSense() ;
1345  groomSoln(direction,primalTolerance_,djs,colLower,colsol,colUpper,columnGap,
1346            rowCopy,rowStartPos,largestNegativeInRow,largestPositiveInRow) ;
1347/*
1348  Determine objective cutoff (minimisation convention).
1349
1350  TODO: Seems to me that this bit of code is a guaranteed fail for a
1351        maximisation problem with no cutoff. It also repeats in a number
1352        of places.  -- lh, 101125 --
1353*/
1354  double cutoff ;
1355  bool cutoff_available = si.getDblParam(OsiDualObjectiveLimit,cutoff) ;
1356  if (!cutoff_available || usingObjective_ < 0) {
1357    cutoff = si.getInfinity() ;
1358  }
1359  cutoff *= direction ;
1360  if (fabs(cutoff) > 1.0e30)
1361    assert (cutoff > 1.0e30) ;
1362  double current = si.getObjValue() ;
1363  current *= direction ;
1364/*
1365  Scan the variables, noting the ones that are fixed or have a bound too large
1366  to be useful.
1367*/
1368  //int nFix = 0 ;
1369  for (int i = 0 ; i < nCols ; i++) {
1370    if (colUpper[i]-colLower[i] < 1.0e-8) {
1371      markC[i] = tightenLower|tightenUpper ;
1372      //nFix++ ;
1373    } else {
1374      markC[i] = 0 ;
1375      if (colUpper[i] > 1.0e10)
1376        markC[i] |= infUpper ;
1377      if (colLower[i] < -1.0e10)
1378        markC[i] |= infLower ;
1379    }
1380  }
1381  //printf("PROBE %d fixed out of %d\n",nFix,nCols) ;
1382
1383/*
1384  jjf: If we are going to replace coefficient then we don't need to be
1385       effective
1386
1387  Seems like this comment does not agree with CglProbingRowCut.addCuts(),
1388  which checks effectiveness before adding a cut to the OsiCuts collection or
1389  entering it into the strenghtenRow array.
1390
1391  But it does agree with the way coefficient strengthening cuts are handled
1392  down in the end-of-iteration processing for the down/up/down probe loop.
1393
1394  It looks like the idea of strengthenRow is that a cut that's simply
1395  strengthening an existing row i is entered in slot i of strengthenRow. It's
1396  left to the client to deal with it on return. I need to get a better idea of
1397  how justReplace is handled, here and in the client.  -- lh, 101125 --
1398*/
1399  //double needEffectiveness = info->strengthenRow ? -1.0e10 : 1.0e-3 ;
1400  double needEffectiveness = info->strengthenRow ? 1.0e-8 : 1.0e-3 ;
1401  if (justReplace && ((info->pass&1) != 0))
1402    needEffectiveness = -1.0e10 ;
1403
1404  double tolerance = 1.0e1*primalTolerance_ ;
1405
1406  /* for both way coding */
1407  int nstackC0 = -1 ;
1408  int nstackR,nstackC ;
1409/*
1410  So, let's assume that the real business starts here, and the previous block
1411  is irrelevant.
1412*/
1413/*
1414  All the initialisation that occurs here is specific to CglTreeProbingInfo.
1415  Sort of begs the question `Why is this handled as a derived class of
1416  CglTreeInfo?' And why (what?) is CglTreeInfo, in the grand scheme of things?
1417
1418  Some grepping about in Cbc says that we might see a CglTreeProbingInfo
1419  object during processing of the root node. My tentative hypothesis is that
1420  it's a vehicle to pass implications back to cbc, where they will be stashed
1421  in a CglImplication object for further use. It's worth mentioning that the
1422  existence of a CglTreeProbingInfo object depends on an independent symbol
1423  (CLIQUE_ANALYSIS) being defined in CbcModel::branchAndBound. It's also worth
1424  mentioning that the code here will core dump if the info object is not a
1425  CglTreeProbingInfo object.
1426
1427  And finally, it's worth mentioning that this code is disabled --- PROBING100
1428  is defined to 0 at the head of CglProbing --- since 080722.
1429
1430  The trivial implementation CglTreeInfo::initializeFixing() always returns 0,
1431  so in effect this block ensures saveFixingInfo is false.
1432
1433  -- lh, 101126 --
1434*/
1435  bool saveFixingInfo = false ;
1436#if PROBING100
1437  CglTreeProbingInfo *probingInfo = dynamic_cast<CglTreeProbingInfo *> (info) ;
1438  const int *backward = NULL ;
1439  const int *integerVariable = NULL ;
1440  const int *toZero = NULL ;
1441  const int *toOne = NULL ;
1442  const fixEntry *fixEntries = NULL ;
1443#endif
1444  if (info->inTree) {
1445#if PROBING100
1446    backward = probingInfo->backward() ;
1447    integerVariable = probingInfo->integerVariable() ;
1448    toZero = probingInfo->toZero() ;
1449    toOne = probingInfo->toOne() ;
1450    fixEntries = probingInfo->fixEntries() ;
1451#endif
1452  } else {
1453    saveFixingInfo = (info->initializeFixing(&si) > 0) ;
1454  }
1455/*
1456  PASS LOOP: HEAD
1457
1458  Major block #2: Main pass loop.
1459
1460  From CglProbingAnn: anyColumnCuts is set only in the case that we've
1461  fixed a variable by probing (i.e., one of the up or down probe resulted
1462  in infeasibility) and that probe entailed column cuts. Once set, it is
1463  never rescinded. In the reworked code, it's set as the return value of
1464  monotoneActions().
1465*/
1466  bool anyColumnCuts = false ;
1467  int ninfeas = 0 ;
1468  int rowCuts = rowCuts_ ;
1469  double disaggEffectiveness = 1.0e-3 ;
1470  int ipass = 0 ;
1471  int nfixed = -1 ;
1472  int maxPass = info->inTree ? maxPass_ : maxPassRoot_ ;
1473
1474  while (ipass < maxPass && nfixed) {
1475    int iLook ;
1476    ipass++ ;
1477    //printf("pass %d\n",ipass) ;
1478    nfixed = 0 ;
1479/*
1480  We went through a fair bit of trouble in gutsOfGenerateCut to determine
1481  the set of variables to be probed and loaded the indices into lookedAt_;
1482  numberThisTime_ reflects that.
1483
1484  The net effect here is that the root gets special treatment
1485  (maxProbeRoot_) and the first pass at the root gets extra special treatment
1486  (numberThisTime_).
1487
1488  See CbcCutGenerator::generateCuts for the origin of magic number 123.
1489  Comments there indicate it's intended to take effect deep in the tree, but
1490  the code here will also work at the root.
1491
1492  Special cases aside, the net effect of the loops is to walk lookedAt_,
1493  promoting every nth variable (n = cutDown) to the front of lookedAt_. When
1494  the loop finishes, the rest of the variables (which were copied off to
1495  stackC) are appended to lookedAt_. If XXXXXX is not defined, we have an
1496  expensive noop.  -- lh, 101126 --
1497*/
1498    int justFix = (!info->inTree && !info->pass) ? -1 : 0 ;
1499    int maxProbe = info->inTree ? maxProbe_ : maxProbeRoot_ ;
1500    if (justFix < 0)
1501      maxProbe = numberThisTime_ ;
1502    if (maxProbe == 123) {
1503      maxProbe = 0 ;
1504      if (!info->inTree) {
1505        if (!info->pass || numberThisTime_ < -100) {
1506          maxProbe = numberThisTime_ ;
1507        } else {
1508          int cutDown = 4 ;
1509          int offset = info->pass%cutDown ;
1510          int i ;
1511          int k = 0 ;
1512          int kk = offset ;
1513          for (i = 0 ; i < numberThisTime_ ; i++) {
1514            if (!kk) {
1515#define XXXXXX
1516#ifdef XXXXXX
1517              lookedAt_[maxProbe] = lookedAt_[i] ;
1518#endif
1519              maxProbe++ ;
1520              kk = cutDown-1 ;
1521            } else {
1522              stackC[k++] = lookedAt_[i] ;
1523              kk-- ;
1524            }
1525          }
1526#ifdef XXXXXX
1527          memcpy(lookedAt_+maxProbe,stackC,k*sizeof(int)) ;
1528#endif
1529        }
1530      } else {
1531        if (numberThisTime_ < 200) {
1532          maxProbe = numberThisTime_ ;
1533        } else {
1534          int cutDown = CoinMax(numberThisTime_/100,4) ;
1535          int offset = info->pass%cutDown ;
1536          int i ;
1537          int k = 0 ;
1538          int kk = offset ;
1539          for (i = 0 ; i < numberThisTime_ ; i++) {
1540            if (!kk) {
1541#ifdef XXXXXX
1542              lookedAt_[maxProbe] = lookedAt_[i] ;
1543#endif
1544              maxProbe++ ;
1545              kk = cutDown-1 ;
1546            } else {
1547              stackC[k++] = lookedAt_[i] ;
1548              kk-- ;
1549            }
1550          }
1551#ifdef XXXXXX
1552          memcpy(lookedAt_+maxProbe,stackC,k*sizeof(int)) ;
1553#endif
1554        }
1555      }
1556    }  // end maxProbe == 123
1557
1558/*
1559  This looks to be an overall limit on probing. It's decremented every time a
1560  variable is popped off stackC for processing.
1561
1562  TODO: PROBING5 would disable this limit; currently inactive. -- lh, 101127 --
1563*/
1564    int leftTotalStack = maxStack*CoinMax(200,maxProbe) ;
1565#ifdef PROBING5
1566    if (!info->inTree&&!info->pass)
1567      leftTotalStack = 1234567890 ;
1568#endif
1569    //printf("maxStack %d maxPass %d numberThisTime %d info pass %d\n",
1570    //   maxStack,maxPass,numberThisTime_,info->pass) ;
1571/*
1572  LOOKEDAT LOOP: HEAD
1573
1574  Main loop to probe each variable in lookedAt_
1575
1576  We're finally ready to get down to the business of probing. Open a loop to
1577  probe each variable in lookedAt_.
1578*/
1579    for (iLook = 0 ; iLook < numberThisTime_ ; iLook++) {
1580      double solval ;
1581      double down ;
1582      double up ;
1583/*
1584  Too successful? Consider bailing out.
1585
1586  If we're generating row cuts, but haven't fixed any variables or we're in
1587  the tree, break from probing.
1588
1589  JustFix < 0 says this is the first pass at the root; for any other
1590  iteration of the pass loop, it'll be initialised to 0.  If it is the first
1591  pass at the root, turn off row cuts, keep on fixing variables, and stop
1592  at the end of this pass. Otherwise, if we haven't fixed any variables, break.
1593
1594  TODO: Otherwise keep going with no changes? That doesn't seem right. The
1595        logic here does not cover all situations.  -- lh, 101126 --
1596*/
1597      if (rowCut.outOfSpace() || leftTotalStack <= 0) {
1598        if (!justFix && (!nfixed || info->inTree)) {
1599#ifdef COIN_DEVELOP
1600          if (!info->inTree)
1601            printf("Exiting a on pass %d, maxProbe %d\n",
1602                   ipass,maxProbe) ;
1603#endif   
1604          break ;
1605        } else if (justFix <= 0) {
1606          if (!info->inTree) {
1607            rowCuts = 0 ;
1608            justFix = 1 ;
1609            disaggEffectiveness = COIN_DBL_MAX ;
1610            needEffectiveness = COIN_DBL_MAX ;
1611            //maxStack=10 ;
1612            maxPass = 1 ;
1613          } else if (!nfixed) {
1614#ifdef COIN_DEVELOP
1615            printf("Exiting b on pass %d, maxProbe %d\n",
1616                   ipass,maxProbe) ;
1617#endif   
1618            break ;
1619          }
1620        }
1621      }
1622/*
1623  awayFromBound isn't quite accurate --- we'll also set awayFromBound = 0
1624  for a general integer variable at an integer value strictly within bounds.
1625*/
1626      int awayFromBound = 1 ;
1627/*
1628  Have a look at the current variable. We're not interested in processing it
1629  if either or both bounds have been improved (0x02|0x01). If the variable has
1630  become fixed, claim both bounds have improved.
1631
1632  TODO: if x<j> is not an integer variable we have big problems. This should
1633        be an assert. -- lh, 101126 --
1634*/
1635      int j = lookedAt_[iLook] ;
1636      solval = colsol[j] ;
1637      down = floor(solval+tolerance) ;
1638      up = ceil(solval-tolerance) ;
1639      if (columnGap[j] < 0.0) markC[j] = tightenUpper|tightenLower ;
1640      if ((markC[j]&(tightenUpper|tightenLower)) != 0 || !intVar[j]) continue ;
1641/*
1642  `Normalize' variables that are near (or past) their bounds.
1643
1644  In normal circumstances, u<j> - l<j> is at least 1, while tol will be
1645  down around 10e-5 (given a primal tolerance of 10e-6, about as large as
1646  you'd normally see). Then
1647
1648  l<j> < l<j>+tol < l<j>+2tol << u<j>-2tol < u<j>-tol < u<j>
1649
1650  For x<j> > u<j>-tol, (x<j>,u<j>,u<j>) => (u<j>-.1, l<j>, u<j>)
1651  For x<j> < l<j>+tol, (x<j>,l<j>,l<j>) => (l<j>+.1, l<j>, u<j>)
1652  For up == down,      (x<j>,v,v)       => (v+.1,v,v+1)
1653
1654  So we're forcing a spread of 1 between up and down, and moving x<j>
1655  to something far enough (.1) from integer to avoid accuracy problems when
1656  calculating movement. A side effect is that primal infeasibility is
1657  corrected.
1658
1659  TODO: Why this structure? Another approach would be to test using up and
1660        down calculated just above. It works out pretty much the same if
1661        you stick with variables of type double. Moving to int would likely
1662        make it faster.  -- lh, 101127 --
1663
1664  TODO: Things will get seriously wierd if (u<j> - l<j>) <= 3tol. For
1665        u<j>-tol < x<j> < l<j>+2tol, values will change. But the symmetric
1666        case u<j>-2tol < x<j> < l<j>+tol will never change values because
1667        of the order of the if statements. Clearly the code is not designed
1668        to cope with this level of pathology.  -- lh, 101127 --
1669 
1670  TODO: It's worth noting that awayFromBound is never used.  -- lh, 101127 --
1671
1672  TODO: It's worth noting that the values set for solval are never used except
1673        in a debug print.  -- lh, 101213 --
1674
1675  TODO: The final case below corresponds to l<j>+1 <= down == up < u<j>-1,
1676        e.g., an interior integer value. Seems like we'd want to probe +1 and
1677        -1. Presumably the code isn't set up to do this, so we arbitrarily
1678        pick the up probe. Could the code be augmented to do the stronger
1679        probe?  -- lh, 101127 --
1680*/
1681      if (solval >= colUpper[j]-tolerance ||
1682          solval <= colLower[j]+tolerance || up == down) {
1683        awayFromBound = 0 ;
1684        if (solval <= colLower[j]+2.0*tolerance) {
1685          solval = colLower[j]+1.0e-1 ;
1686          down = colLower[j] ;
1687          up = down+1.0 ;
1688        } else if (solval >= colUpper[j]-2.0*tolerance) {
1689          solval = colUpper[j]-1.0e-1 ;
1690          up = colUpper[j] ;
1691          down = up-1 ;
1692        } else {
1693          up = down+1.0 ;
1694          solval = down+1.0e-1 ;
1695        }
1696      }
1697      assert (up <= colUpper[j]) ;
1698      assert (down >= colLower[j]) ;
1699      assert (up > down) ;
1700#     if CGL_DEBUG > 0
1701      const double lj = colLower[j] ;
1702      const double uj = colUpper[j] ;
1703      const bool downIsLower = (fabs(down-colLower[j]) < 1.0e-7) ;
1704      const bool upIsUpper = (fabs(up-colUpper[j]) < 1.0e-7) ;
1705#     endif
1706/*
1707  Set up to probe each variable down (1), up (2), down (1).
1708
1709  The notion is that we'll set a bit (1, 2, 4) in the result record if we're
1710  feasible for a particular way. Way is defined as:
1711    1: we're looking at movement between floor(x*<j>) and x*<j>, u<j>
1712    2: we're looking at movement between ceil(x*<j>) and x*<j>, l<j>
1713  As defined here, movement will be negative for way = 1.
1714
1715  Why do we have a second pass with way = 1? Glad you asked. About 1200
1716  lines farther down, it finally becomes apparent that we execute the
1717  final pass only when the initial down probe was feasible, but the up probe
1718  was not. In a nutshell, we're restoring the state produced by the down
1719  probe, which we'll then nail down in the section of code labelled `keep',
1720  a mere 600 lines farther on.
1721*/
1722      unsigned int iway ;
1723      unsigned int way[] = { probeDown, probeUp, probeDown } ;
1724      unsigned int feasValue[] =
1725          { downProbeFeas, upProbeFeas, downProbe2Feas } ;
1726      unsigned int feasRecord = 0 ;
1727      bool notFeasible = false ;
1728      int istackC = 0 ;
1729      int istackR = 0 ;
1730/*
1731  PROBE LOOP: HEAD
1732
1733  Open a loop to probe up and down for the current variable. Get things
1734  started by placing the variable on the propagation queue (stackC).
1735
1736  As with the previous loops (variables in lookedAt_, passes) this loop
1737  extends to the end of major block #2).
1738*/
1739      for (iway = downProbe ; iway < oneProbeTooMany ; iway++) {
1740        stackC[0] = j ;
1741        markC[j] = way[iway] ;
1742/*
1743  Calculate movement given the current direction. Given the setup above,
1744  movement should be at least +/-1.
1745
1746  TODO: Notice we reset up or down yet again. Really, this shouldn't be
1747        necessary. For that matter, why did we bother with awayFromBound
1748        immediately above? -- lh, 101127 --
1749
1750  TODO: The more I think about this, the less happy I am. Seems to me that
1751        neither colUpper or colLower can change during iterations of this
1752        loop. If we discover monotone, we save bounds and break. If we
1753        don't discover monotone, we restore. Let's test that.
1754        -- lh, 110105 --
1755*/
1756        double solMovement ;
1757        double movement ;
1758        int goingToTrueBound = 0 ;
1759
1760#       if CGL_DEBUG > 0
1761        assert(lj == colLower[j]) ;
1762        assert(uj == colUpper[j]) ;
1763        assert(downIsLower == (fabs(down-colLower[j]) < 1.0e-7)) ;
1764        assert(upIsUpper == (fabs(up-colUpper[j]) < 1.0e-7)) ;
1765#       endif
1766
1767        if (way[iway] == probeDown) {
1768          movement = down-colUpper[j] ;
1769          solMovement = down-colsol[j] ;
1770          assert(movement < -0.99999) ;
1771          if (fabs(down-colLower[j]) < 1.0e-7) {
1772            goingToTrueBound = 2 ;
1773            down = colLower[j] ;
1774          }
1775        } else {
1776          movement = up-colLower[j] ;
1777          solMovement = up-colsol[j] ;
1778          assert(movement > 0.99999) ;
1779          if (fabs(up-colUpper[j]) < 1.0e-7) {
1780            goingToTrueBound = 2 ;
1781            up = colUpper[j] ;
1782          }
1783        }
1784/*
1785  Coding for goingToTrueBound is:
1786    0: We're somewhere in the interior of a general integer.
1787    1: We're heading for one of the bounds on a general integer.
1788    2: We're processing a binary variable (u<j>-l<j> = 1 and l<j> = 0).
1789  You can view this next test as correcting for the fact that the previous
1790  code assumes binary variables are all there is.
1791
1792  TODO: Now that I've figured out the math for the constraints generated
1793        by disaggCuts (see typeset documentation), I see what's wrong
1794        (?) here. The disagg cut that's generated relies on the new probe
1795        bound being distance one from the original bound. But (for example)
1796        that's (colUpper-down) = 1 for a down probe.  Not quite what's
1797        being checked for here. But this next test ensures binary variables,
1798        and then the tests are equivalent. Before I willy-nilly change this, I
1799        need to sort out a couple of other places where goingToTrueBound
1800        controls behaviour.  -- lh, 101214 --
1801
1802  TODO: Apropos the above, for example, the test for coefficient strengthening
1803        cuts includes goingToTrueBound != 0  -- lh, 110105 --
1804*/
1805        if (goingToTrueBound && (colUpper[j]-colLower[j] > 1.5 || colLower[j]))
1806          goingToTrueBound = 1 ;
1807/*
1808  If the user doesn't want disaggregation cuts, pretend we're in the interior
1809  of a general integer variable.
1810
1811  TODO: Why is row strengthening limited to binary variables? If we can
1812        figure out what to do, the limitation seems artificial.
1813        -- lh, 101127 --
1814        The test here also says that we can't have coefficient strengthening
1815        without disaggregation. I don't see any reason for this dominance.
1816        -- lh, 110105 --
1817*/
1818        if ((rowCuts&1) == 0)
1819          goingToTrueBound = 0 ;
1820        bool canReplace = info->strengthenRow && (goingToTrueBound == 2) ;
1821
1822#ifdef PRINT_DEBUG
1823        // Alert for general integer with nontrivial range.
1824        // Also last use of solval
1825        if (fabs(movement)>1.01) {
1826          printf("big %d %g %g %g\n",j,colLower[j],solval,colUpper[j]) ;
1827        }
1828#endif
1829/*
1830  Recall that we adjusted the reduced costs (djs) and current objective
1831  (current) to the minimisation sign convention. objVal will accumulate
1832  objective change due to forced bound changes.
1833*/
1834        double objVal = current ;
1835        if (solMovement*djs[j] > 0.0)
1836          objVal += solMovement*djs[j] ;
1837        nstackC = 1 ;
1838        nstackR = 0 ;
1839        saveL[0] = colLower[j] ;
1840        saveU[0] = colUpper[j] ;
1841        assert (saveU[0] > saveL[0]) ;
1842        notFeasible = false ;
1843/*
1844  Recalculate the upper (down probe) or lower (up probe) bound. You can see
1845  from the debug print that arbitrary integers are an afterthought in this
1846  code.
1847
1848  TODO: And why not just set the bounds to down (colUpper) or up (colLower)?
1849        We set movement = down-colUpper just above, and we've taken a lot of
1850        care to groom up and down. This is pointless effort.  -- lh, 101214 --
1851*/
1852        if (movement < 0.0) {
1853          colUpper[j] += movement ;
1854          colUpper[j] = floor(colUpper[j]+0.5) ;
1855          columnGap[j] = colUpper[j]-colLower[j]-primalTolerance_ ;
1856#ifdef PRINT_DEBUG
1857          printf("** Trying %d down to 0\n",j) ;
1858#endif
1859        } else {
1860          colLower[j] += movement ;
1861          colLower[j] = floor(colLower[j]+0.5) ;
1862          columnGap[j] = colUpper[j]-colLower[j]-primalTolerance_ ;
1863#ifdef PRINT_DEBUG
1864          printf("** Trying %d up to 1\n",j) ;
1865#endif
1866        }
1867/*
1868  Is this variable now fixed?
1869*/
1870        if (fabs(colUpper[j]-colLower[j]) < 1.0e-6)
1871          markC[j] = tightenUpper|tightenLower ;
1872/*
1873  Clear the infinite bound bits (infUpper|infLower) and check again. Bounds
1874  may have improved.
1875*/
1876        markC[j] &= ~(infUpper|infLower) ;
1877        if (colUpper[j] > 1.0e10)
1878          markC[j] |= infUpper ;
1879        if (colLower[j] < -1.0e10)
1880          markC[j] |= infLower ;
1881        istackC = 0 ;
1882/*
1883  Update row bounds to reflect the change in variable bound.
1884
1885  TODO: Replacing code to update row bounds with updateRowBounds(). We'll
1886        see how it looks.  -- lh, 101203 --
1887*/
1888
1889  if (!updateRowBounds(j,movement,columnStart,columnLength,row,columnElements,
1890                       rowUpper,rowLower,nstackR,stackR,markR,
1891                       minR,maxR,saveMin,saveMax)) {
1892    notFeasible = true ;
1893    istackC = 1 ;
1894  }
1895/*
1896  PROBE LOOP: BEGIN PROPAGATION
1897
1898  Row bounds are now adjusted for all rows with a<ij> != 0. Time to consider
1899  the effects. nstackC is incremented each time we add a variable, istackC is
1900  incremented each time we finish processing it.
1901
1902  If we lost feasibility above, istackC = nstackC = 1 and this loop will not
1903  execute.
1904
1905  TODO: Is stackC really a stack? Or a queue? And what is the role of
1906        maxStack? stackC is allocated at 2*ncols, but maxStack defaults to
1907        50. -- lh, 101127 --
1908
1909  TODO: stackC is clearly a queue. Allocation strategy is still unclear.
1910        -- lh, 101128 --
1911
1912  jjf:  could be istackC<maxStack?
1913*/
1914        while (istackC < nstackC && nstackC < maxStack) {
1915          leftTotalStack-- ;
1916          int jway ;
1917          int jcol = stackC[istackC] ;
1918          jway = markC[jcol] ;
1919/*
1920  jjf: If not first and fixed then skip
1921
1922  It's clear that x<j>, the probing variable, can be marked as fixed at this
1923  point (happens just above, when the upper or lower bound is adjusted). And
1924  there's an earlier check that avoids probing a variable that's fixed when
1925  plucked from lookedAt_. But ... why leave the if with an empty body?
1926
1927  TODO: Disabled since r118 050225, but seems like the right thing to do. Is it
1928        somehow guaranteed that fixed variables are not stacked?
1929        -- lh, 101127 --
1930
1931  TODO: Based on what's happening down below, it looks like we can encounter
1932        variables here which have had both bounds tightened and are then
1933        pushed onto istackC to propagate the effect of that. This bit of code
1934        should simply go away.  -- lh, 101127 --
1935*/
1936          if (((jway&(tightenLower|tightenUpper)) ==
1937               (tightenLower|tightenUpper)) &&
1938              istackC) {
1939            //istackC++ ;
1940            //continue ;
1941            //printf("fixed %d on stack\n",jcol) ;
1942          }
1943#if PROBING100
1944/*
1945  This block of code has been disabled since r661 080722. I have no notes on
1946  it from my old annotated CglProbing. Given the use of backward, it's tied to
1947  CglTreeProbingInfo.  -- lh, 101127 --
1948*/
1949          if (backward) {
1950            int jColumn = backward[jcol] ;
1951            if (jColumn>=0) {
1952              int nFix=0 ;
1953              // 0-1 see what else could be fixed
1954              if (jway==tightenUpper) {
1955                // fixed to 0
1956                int j ;
1957                for (j = toZero_[jColumn];j<toOne_[jColumn];j++) {
1958                  int kColumn=fixEntry_[j].sequence ;
1959                  kColumn = integerVariable_[kColumn] ;
1960                  bool fixToOne = fixEntry_[j].oneFixed ;
1961                  if (fixToOne) {
1962                    if (colLower[kColumn]==0.0) {
1963                      if (colUpper[kColumn]==1.0) {
1964                        // See if on list
1965                        if (!(markC[kColumn]&(tightenLower|tightenUpper))) {
1966                          if(nStackC<nCols) {
1967                            stackC[nstackC]=kColumn ;
1968                            saveL[nstackC]=colLower[kColumn] ;
1969                            saveU[nstackC]=colUpper[kColumn] ;
1970                            assert (saveU[nstackC]>saveL[nstackC]) ;
1971                            assert (nstackC<nCols) ;
1972                            nstackC++ ;
1973                            markC[kColumn]|=tightenLower ;
1974                            nFix++ ;
1975                          }
1976                        } else if ((markC[kColumn]&(tightenLower|tightenUpper))==tightenUpper) {
1977                          notFeasible = true ;
1978                        }
1979                      } else {
1980                        // infeasible!
1981                        notFeasible = true ;
1982                      }
1983                    }
1984                  } else {
1985                    if (colUpper[kColumn]==1.0) {
1986                      if (colLower[kColumn]==0.0) {
1987                        // See if on list
1988                        if (!(markC[kColumn]&(tightenLower|tightenUpper))) {
1989                          if(nStackC<nCols) {
1990                            stackC[nstackC]=kColumn ;
1991                            saveL[nstackC]=colLower[kColumn] ;
1992                            saveU[nstackC]=colUpper[kColumn] ;
1993                            assert (saveU[nstackC]>saveL[nstackC]) ;
1994                            assert (nstackC<nCols) ;
1995                            nstackC++ ;
1996                            markC[kColumn]|=tightenUpper ;
1997                            nFix++ ;
1998                          }
1999                        } else if ((markC[kColumn]&(tightenLower|tightenUpper))==tightenLower) {
2000                          notFeasible = true ;
2001                        }
2002                      } else {
2003                        // infeasible!
2004                        notFeasible = true ;
2005                      }
2006                    }
2007                  }
2008                }
2009              } else if (jway==tightenLower) {
2010                int j ;
2011                for ( j=toOne_[jColumn];j<toZero_[jColumn+1];j++) {
2012                  int kColumn=fixEntry_[j].sequence ;
2013                  kColumn = integerVariable_[kColumn] ;
2014                  bool fixToOne = fixEntry_[j].oneFixed ;
2015                  if (fixToOne) {
2016                    if (colLower[kColumn]==0.0) {
2017                      if (colUpper[kColumn]==1.0) {
2018                        // See if on list
2019                        if (!(markC[kColumn]&(tightenLower|tightenUpper))) {
2020                          if(nStackC<nCols) {
2021                            stackC[nstackC]=kColumn ;
2022                            saveL[nstackC]=colLower[kColumn] ;
2023                            saveU[nstackC]=colUpper[kColumn] ;
2024                            assert (saveU[nstackC]>saveL[nstackC]) ;
2025                            assert (nstackC<nCols) ;
2026                            nstackC++ ;
2027                            markC[kColumn]|=tightenLower ;
2028                            nFix++ ;
2029                          }
2030                        } else if ((markC[kColumn]&(tightenLower|tightenUpper))==tightenUpper) {
2031                          notFeasible = true ;
2032                        }
2033                      } else {
2034                        // infeasible!
2035                        notFeasible = true ;
2036                      }
2037                    }
2038                  } else {
2039                    if (colUpper[kColumn]==1.0) {
2040                      if (colLower[kColumn]==0.0) {
2041                        // See if on list
2042                        if (!(markC[kColumn]&(tightenLower|tightenUpper))) {
2043                          if(nStackC<nCols) {
2044                            stackC[nstackC]=kColumn ;
2045                            saveL[nstackC]=colLower[kColumn] ;
2046                            saveU[nstackC]=colUpper[kColumn] ;
2047                            assert (saveU[nstackC]>saveL[nstackC]) ;
2048                            assert (nstackC<nCols) ;
2049                            nstackC++ ;
2050                            markC[kColumn]|=tightenUpper ;
2051                            nFix++ ;
2052                          }
2053                        } else if ((markC[kColumn]&(tightenLower|tightenUpper))==tightenLower) {
2054                          notFeasible = true ;
2055                        }
2056                      } else {
2057                        // infeasible!
2058                        notFeasible = true ;
2059                      }
2060                    }
2061                  }
2062                }
2063              }
2064            }
2065          }
2066#endif  // PROBING100 (disabled)
2067/*
2068  PROBE LOOP: WALK COLUMN
2069
2070  Loop to walk the column of a variable popped off stackC.
2071
2072  We've pulled x<jcol> off stackC. We're going to walk the column and process
2073  each row where a<i,jcol> != 0. Note that we've already updated the row
2074  bounds to reflect the change in bound for x<jcol>. In effect, we've stacked
2075  this variable as a surrogate for stacking the rows whose bounds were
2076  changed by the change to bounds on this variable. Now we're going to examine
2077  those rows and see if any look promising for further propagation.
2078
2079*/
2080          for (int k = columnStart[jcol] ;
2081               k < columnStart[jcol]+columnLength[jcol] ; k++) {
2082            if (notFeasible)
2083              break ;
2084            int irow = row[k] ;
2085            // Row already processed with these bounds.
2086            if (markR[irow]/10000 > 0) continue ;
2087
2088            int rStart = rowStart[irow] ;
2089            int rEnd = rowStartPos[irow] ;
2090            double rowUp = rowUpper[irow] ;
2091            double rowUp2 = 0.0 ;
2092            bool doRowUpN ;
2093            bool doRowUpP ;
2094/*
2095  So, is this row promising?
2096
2097  If our current L<i> (minR) is larger than the original b<i> (rowUp), we're
2098  infeasible.
2099
2100  Otherwise, a bit of linear algebra brings the conclusion that if the
2101  largest negative gap (min{j in N} a<ij>(u<j>-l<j>)) in the row is less than
2102  -(b<i>-L<i>), we aren't going to be able to make any progress shrinking the
2103  domains of variables with negative coefficients.  Similarly, if the largest
2104  positive gap (max{j in P} a<ij>(u<j>-l<j>) is greater than b<i>-L<i>, we
2105  can't shrink the domain of any variable with a positive coefficient.
2106
2107  There are analogous conditions using U<i> and blow<i>.
2108
2109  Note that we don't update the largest negative and positive gaps, so as
2110  propagation grinds on, this simple test becomes less accurate. On the other
2111  hand, the gaps (b<i>-L<i>) and (U<i>-blow<i>) are also shrinking. Still, it's
2112  another justification for the strict limits on propagation.
2113
2114  Summary:
2115
2116  minR = SUM{P}(a<ij>l<j>) + SUM{N}(a<ij>u<j>)
2117  maxR = SUM{P}(a<ij>u<j>) + SUM{N}(a<ij>l<j>)
2118
2119  doRowUpN: a<ij> < 0, minR can violate rowUp => increase l<j> (more negative)
2120  doRowLoN: a<ij> < 0, maxR can violate rowLo => decrease u<j> (more positive)
2121  doRowUpP: a<ij> > 0, minR can violate rowUp => decrease u<j> (less positive)
2122  doRowLoP: a<ij> > 0, maxR can violate rowLo => increase l<j> (less negative)
2123
2124  Note that the bound (l<j>, u<j>) that can be tightened in any given
2125  situation is *not* the bound involved in minR or maxR.
2126
2127  For binary variables, of course, `shrink the domain' corresponds to fixing
2128  the variable.
2129*/
2130
2131            if (rowUp < 1.0e10) {
2132              doRowUpN = true ;
2133              doRowUpP = true ;
2134              rowUp2 = rowUp-minR[irow] ;
2135              if (rowUp2 < -primalTolerance_) {
2136                notFeasible = true ;
2137                break ;
2138              } else {
2139                if (rowUp2+largestNegativeInRow[irow] > 0)
2140                  doRowUpN = false ;
2141                if (rowUp2-largestPositiveInRow[irow] > 0)
2142                  doRowUpP = false ;
2143              }
2144            } else {
2145              doRowUpN = false ;
2146              doRowUpP = false ;
2147              rowUp2 = COIN_DBL_MAX ;
2148            }
2149            double rowLo = rowLower[irow] ;
2150            double rowLo2 = 0.0 ;
2151            bool doRowLoN ;
2152            bool doRowLoP ;
2153            if (rowLo > -1.0e10) {
2154              doRowLoN = true ;
2155              doRowLoP = true ;
2156              rowLo2 = rowLo-maxR[irow] ;
2157              if (rowLo2 > primalTolerance_) {
2158                notFeasible = true ;
2159                break ;
2160              } else {
2161                if (rowLo2-largestNegativeInRow[irow] < 0)
2162                  doRowLoN = false ;
2163                if (rowLo2+largestPositiveInRow[irow] < 0)
2164                  doRowLoP = false ;
2165              }
2166            } else {
2167              doRowLoN = false ;
2168              doRowLoP = false ;
2169              rowLo2 = -COIN_DBL_MAX ;
2170            }
2171            markR[irow] += 10000 ;
2172/*
2173  LOU_DEBUG: see if we're processing again without a bound change.
2174            if (markR[irow]/10000 > 1)
2175              std::cout
2176                << "REDUNDANT: processing " << irow
2177                << " for the " << markR[irow]/10000 << " time." << std::endl ;
2178*/
2179/*
2180  PROBE LOOP: WALK ROW
2181
2182  Given the above analysis, work on the columns with negative coefficients.
2183
2184  doRowUpN: a<ij> < 0, minR can violate rowUp => increase l<j> (more negative)
2185  doRowLoN: a<ij> < 0, maxR can violate rowLo => decrease u<j> (more positive)
2186
2187  TODO: The asserts here should be compiled out unless we're debugging.
2188        -- lh, 101210 --
2189*/
2190            if (doRowUpN && doRowLoN) {
2191              //doRowUpN=doRowLoN=false ;
2192              // Start neg values loop
2193              for (int kk = rStart ; kk < rEnd ; kk++) {
2194                int kcol = column[kk] ;
2195                int markIt = markC[kcol] ;
2196                // Skip columns with fixed variables.
2197                if ((markIt&(tightenLower|tightenUpper)) != (tightenLower|tightenUpper)) {
2198                  double value2 = rowElements[kk] ;
2199                  if (colUpper[kcol] <= 1e10)
2200                    assert ((markIt&infUpper) == 0) ;
2201                  else
2202                    assert ((markIt&infUpper) != 0) ;
2203                  if (colLower[kcol] >= -1e10)
2204                    assert ((markIt&infLower) == 0) ;
2205                  else
2206                    assert ((markIt&infLower) != 0) ;
2207                  assert (value2 < 0.0) ;
2208/*
2209  Not every column will be productive. Can we do anything with this one?
2210*/
2211                  double gap = columnGap[kcol]*value2 ;
2212                  bool doUp = (rowUp2+gap < 0.0) ;
2213                  bool doDown = (rowLo2-gap > 0.0) ;
2214                  if (doUp || doDown) {
2215                    double moveUp = 0.0 ;
2216                    double moveDown = 0.0 ;
2217                    double newUpper = -1.0 ;
2218                    double newLower = 1.0 ;
2219/*
2220  doUp attempts to increase the lower bound. The math looks like this:
2221    a<irow,kcol>x<kcol> + (minR[irow]-a<irow,kcol>u<kcol>) <= b<irow>
2222    value2*x<kcol> - value2*u<kcol> <= (b<irow>-minR[irow]) = rowUp2
2223    x<kcol> - u<kcol> >= rowUp2/value2
2224    l<kcol> = u<kcol> + rowUp2/value2
2225  Hence we need u<kcol> to be finite (0x8 not set). We also refuse to tighten
2226  l<kcol> a second time (0x2 not set).
2227
2228  TODO: So, why don't we allow tightening a second time?  Clearly we might
2229        process a row on some other iteration that imposes a tighter bound.
2230        Is this an optimisation for binary variables, or is there some
2231        deeper issue at work in terms of correctness?  -- lh, 101127 --
2232
2233  TODO: Also, it's clear that this bit of code would like to handle
2234        continuous variables, but it's also clear that it does it incorrectly.
2235        When the new lower bound straddles the existing upper bound, that
2236        looks to me to be sufficient justification to fix the variable. But
2237        the code here does just the opposite: merely tightening the bound is
2238        flagged as fixing the variable, with no consideration that we may have
2239        gone infeasible. But if we prove we straddle the upper bound, all we
2240        claim is we've tightened the upper bound!  -- lh, 101127 --
2241*/
2242                    if (doUp && ((markIt&(tightenLower|infUpper)) == 0)) {
2243                      double dbound = colUpper[kcol]+rowUp2/value2 ;
2244                      if (intVar[kcol]) {
2245                        markIt |= tightenLower ;
2246                        newLower = ceil(dbound-primalTolerance_) ;
2247                      } else {
2248                        newLower = dbound ;
2249                        if (newLower+primalTolerance_ > colUpper[kcol] &&
2250                            newLower-primalTolerance_ <= colUpper[kcol]) {
2251                          newLower = colUpper[kcol] ;
2252                          markIt |= tightenLower ;
2253                          //markIt = (tightenUpper|tightenLower) ;
2254                        } else {
2255                          // avoid problems - fix later ?
2256                          markIt = (tightenUpper|tightenLower) ;
2257                        }
2258                      }
2259                      moveUp = newLower-colLower[kcol] ;
2260                    }
2261/*
2262  The same, but attempt to decrease the upper bound by comparison with the
2263  blow for the row.
2264*/
2265                    if (doDown && ((markIt&(tightenUpper|infLower)) == 0)) {
2266                      double dbound = colLower[kcol] + rowLo2/value2 ;
2267                      if (intVar[kcol]) {
2268                        markIt |= tightenUpper ;
2269                        newUpper = floor(dbound+primalTolerance_) ;
2270                      } else {
2271                        newUpper = dbound ;
2272                        if (newUpper-primalTolerance_ < colLower[kcol] &&
2273                            newUpper+primalTolerance_ >= colLower[kcol]) {
2274                          newUpper = colLower[kcol] ;
2275                          markIt |= tightenUpper ;
2276                          //markIt = (tightenUpper|tightenLower) ;
2277                        } else {
2278                          // avoid problems - fix later ?
2279                          markIt = (tightenUpper|tightenLower) ;
2280                        }
2281                      }
2282                      moveDown = newUpper-colUpper[kcol] ;
2283                    }
2284/*
2285  Is there any change to propagate? Assuming we haven't exceeded the limit
2286  for propagating, queue the variable. (If the variable's already on the list,
2287  that's not a problem.)
2288
2289  Note that markIt is initially loaded from markC and whatever bounds are
2290  changed are then or'd in. So changes accumulate.
2291*/
2292                    if (!moveUp && !moveDown)
2293                      continue ;
2294                    bool onList = ((markC[kcol]&(tightenLower|tightenUpper)) != 0) ;
2295                    if (nstackC < 2*maxStack) {
2296                      markC[kcol] = markIt ;
2297                    }
2298                    if (moveUp && (nstackC < 2*maxStack)) {
2299                      if (!onList) {
2300                        stackC[nstackC] = kcol ;
2301                        saveL[nstackC] = colLower[kcol] ;
2302                        saveU[nstackC] = colUpper[kcol] ;
2303                        assert (saveU[nstackC] > saveL[nstackC]) ;
2304                        assert (nstackC < nCols) ;
2305                        nstackC++ ;
2306                        onList = true ;
2307                      }
2308/*
2309  The first assert here says `If we're minimising, and d<j> < 0, then
2310  x<j> must be NBUB, so if the new l<j> > x*<j> = u<j>, we're infeasible.'
2311
2312  TODO: Why haven't we detected infeasibility? Put another way, why isn't
2313        there an assert(notFeasible)?  -- lh, 101127 --
2314
2315  TODO: Seems like the change in obj should be ((new l<j>) - x*<j>)*d<j>, but
2316        moveUp is (newLower-colLower). Was there some guarantee that x*<j>
2317        equals l<j>?  -- lh, 101127 --
2318  TODO: I think this is another symptom of incomplete conversion for general
2319        integers. For binary variables, it's not an issue.  -- lh, 101203 --
2320*/
2321                      if (newLower > colsol[kcol]) {
2322                        if (djs[kcol] < 0.0) {
2323                          assert (newLower > colUpper[kcol]+primalTolerance_) ;
2324                        } else {
2325                          objVal += moveUp*djs[kcol] ;
2326                        }
2327                      }
2328/*
2329  TODO: Up above we've already set newLower = ceil(dbound-primalTolerance_).
2330        It's highly unlikely that ceil(newLower-1.0e4) will be different.
2331        -- lh, 101127 --
2332  TODO: My first inclination is to claim that newLower should never be worse
2333        than the existing bound. But I'd better extend the benefit of the
2334        doubt for a bit. Some rows will produce stronger bounds than others.
2335        -- lh, 101127 --
2336*/
2337                      if (intVar[kcol]) 
2338                        newLower = CoinMax(colLower[kcol],
2339                                           ceil(newLower-1.0e-4)) ;
2340                      colLower[kcol] = newLower ;
2341                      columnGap[kcol] = colUpper[kcol]-newLower-primalTolerance_ ;
2342/*
2343  Propagate the column bound change.
2344
2345  TODO: Notice that we're not setting a variable here when we discover
2346        infeasibility; rather, we're setting the column bound to a ridiculous
2347        value which will be discovered farther down.  -- lh, 101127 --
2348  TODO: Replacing code to update row bounds with updateRowBounds(). We'll
2349        see how it looks.  -- lh, 101203 --
2350*/
2351                      if (fabs(colUpper[kcol]-colLower[kcol]) < 1.0e-6) {
2352                        markC[kcol] = (tightenLower|tightenUpper) ;
2353                      }
2354                      markC[kcol] &= ~(infLower|infUpper) ;
2355                      if (colUpper[kcol] > 1.0e10)
2356                        markC[kcol] |= infUpper ;
2357                      if (colLower[kcol] < -1.0e10)
2358                        markC[kcol] |= infLower ;
2359
2360                      if (!updateRowBounds(kcol,moveUp,
2361                               columnStart,columnLength,row,columnElements,
2362                               rowUpper,rowLower,nstackR,stackR,markR,
2363                               minR,maxR,saveMin,saveMax)) {
2364                        colLower[kcol] = 1.0e50 ;
2365                      }
2366                    }
2367/*
2368  Repeat the whole business, in the down direction. Stack the variable, if
2369  it's not already there, do the infeasibility check / objective update, and
2370  update minR / maxR for affected constraints.
2371*/
2372                    if (moveDown&&nstackC<2*maxStack) {
2373                      if (!onList) {
2374                        stackC[nstackC]=kcol ;
2375                        saveL[nstackC]=colLower[kcol] ;
2376                        saveU[nstackC]=colUpper[kcol] ;
2377                        assert (saveU[nstackC]>saveL[nstackC]) ;
2378                        assert (nstackC<nCols) ;
2379                        nstackC++ ;
2380                        onList=true ;
2381                      }
2382                      if (newUpper<colsol[kcol]) {
2383                        if (djs[kcol]>0.0) {
2384                          /* should be infeasible */
2385                          assert (colLower[kcol]>newUpper+primalTolerance_) ;
2386                        } else {
2387                          objVal += moveDown*djs[kcol] ;
2388                        }
2389                      }
2390                      if (intVar[kcol])
2391                        newUpper = CoinMin(colUpper[kcol],floor(newUpper+1.0e-4)) ;
2392                      colUpper[kcol]=newUpper ;
2393                      columnGap[kcol] = newUpper-colLower[kcol]-primalTolerance_ ;
2394                      if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
2395                        markC[kcol] = (tightenLower|tightenUpper); // say fixed
2396                      }
2397                      markC[kcol] &= ~(infLower|infUpper) ;
2398                      if (colUpper[kcol]>1.0e10)
2399                        markC[kcol] |= infUpper ;
2400                      if (colLower[kcol]<-1.0e10)
2401                        markC[kcol] |= infLower ;
2402
2403                      if (!updateRowBounds(kcol,moveDown,
2404                               columnStart,columnLength,row,columnElements,
2405                               rowUpper,rowLower,nstackR,stackR,markR,
2406                               minR,maxR,saveMin,saveMax)) {
2407                        colUpper[kcol] = -1.0e50 ;
2408                      }
2409                    }
2410/*
2411  We've propagated the bound changes for this column. Check to see if we
2412  discovered infeasibility. Abort by claiming we've cleared the propagation
2413  stack.
2414*/
2415                    if (colLower[kcol] > colUpper[kcol]+primalTolerance_) {
2416                      notFeasible = true ;
2417                      k = columnStart[jcol]+columnLength[jcol] ;
2418                      istackC = nstackC+1 ;
2419                      break ;
2420                    }
2421                  }   // end if (doUp || doDown) (productive column)
2422                }  // end column not fixed
2423              } // end loop on negative coefficients of this row
2424            } else if (doRowUpN) {
2425/*
2426  The above block propagated change for a variable with a negative coefficient,
2427  where we could work against both the row upper and lower bounds. Now we're
2428  going to do it again, but specialised for the case where we can only work
2429  against the upper bound.
2430*/
2431              // Start neg values loop
2432              for (int kk = rStart ; kk < rEnd ; kk++) {
2433                int kcol =column[kk] ;
2434                int markIt=markC[kcol] ;
2435                if ((markIt&(tightenLower|tightenUpper))!= (tightenLower|tightenUpper)) {
2436                  double value2=rowElements[kk] ;
2437                  double gap = columnGap[kcol]*value2 ;
2438                  if (!(rowUp2 + gap < 0.0))
2439                    continue ;
2440                  double moveUp=0.0 ;
2441                  double newLower=1.0 ;
2442                  if ((markIt&(tightenLower|infUpper))==0) {
2443                    double dbound = colUpper[kcol]+rowUp2/value2 ;
2444                    if (intVar[kcol]) {
2445                      markIt |= tightenLower ;
2446                      newLower = ceil(dbound-primalTolerance_) ;
2447                    } else {
2448                      newLower=dbound ;
2449                      if (newLower+primalTolerance_>colUpper[kcol]&&
2450                          newLower-primalTolerance_<=colUpper[kcol]) {
2451                        newLower=colUpper[kcol] ;
2452                        markIt |= tightenLower ;
2453                        //markIt= (tightenLower|tightenUpper) ;
2454                      } else {
2455                        // avoid problems - fix later ?
2456                        markIt= (tightenLower|tightenUpper) ;
2457                      }
2458                    }
2459                    moveUp = newLower-colLower[kcol] ;
2460                    if (!moveUp)
2461                      continue ;
2462                    bool onList = ((markC[kcol]&(tightenLower|tightenUpper))!=0) ;
2463                    if (nstackC<2*maxStack) {
2464                      markC[kcol] = markIt ;
2465                    }
2466                    if (moveUp&&nstackC<2*maxStack) {
2467                      if (!onList) {
2468                        stackC[nstackC]=kcol ;
2469                        saveL[nstackC]=colLower[kcol] ;
2470                        saveU[nstackC]=colUpper[kcol] ;
2471                        assert (saveU[nstackC]>saveL[nstackC]) ;
2472                        assert (nstackC<nCols) ;
2473                        nstackC++ ;
2474                        onList=true ;
2475                      }
2476                      if (newLower>colsol[kcol]) {
2477                        if (djs[kcol]<0.0) {
2478                          /* should be infeasible */
2479                          assert (newLower>colUpper[kcol]+primalTolerance_) ;
2480                        } else {
2481                          objVal += moveUp*djs[kcol] ;
2482                        }
2483                      }
2484                      if (intVar[kcol]) 
2485                        newLower = CoinMax(colLower[kcol],ceil(newLower-1.0e-4)) ;
2486                      colLower[kcol]=newLower ;
2487                      columnGap[kcol] = colUpper[kcol]-newLower-primalTolerance_ ;
2488                      if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
2489                        markC[kcol]= (tightenLower|tightenUpper); // say fixed
2490                      }
2491                      markC[kcol] &= ~(infLower|infUpper) ;
2492                      if (colUpper[kcol]>1.0e10)
2493                        markC[kcol] |= infUpper ;
2494                      if (colLower[kcol]<-1.0e10)
2495                        markC[kcol] |= infLower ;
2496
2497                      if (!updateRowBounds(kcol,moveUp,
2498                               columnStart,columnLength,row,columnElements,
2499                               rowUpper,rowLower,nstackR,stackR,markR,
2500                               minR,maxR,saveMin,saveMax)) {
2501                        colLower[kcol] = 1.0e50 ;
2502                      }
2503                    }
2504                    if (colLower[kcol] > colUpper[kcol]+primalTolerance_) {
2505                      notFeasible = true ;
2506                      k = columnStart[jcol]+columnLength[jcol] ;
2507                      istackC = nstackC+1 ;
2508                      break ;
2509                    }
2510                  }
2511                }
2512              } // end big loop rStart->rPos
2513            } else if (doRowLoN) {
2514/*
2515  And yet again, for the case where we can only work against the lower bound.
2516*/
2517              // Start neg values loop
2518              for (int kk = rStart ; kk < rEnd ; kk++) {
2519                int kcol =column[kk] ;
2520                if ((markC[kcol]&(tightenLower|tightenUpper))!= (tightenLower|tightenUpper)) {
2521                  double moveDown=0.0 ;
2522                  double newUpper=-1.0 ;
2523                  double value2=rowElements[kk] ;
2524                  int markIt=markC[kcol] ;
2525                  assert (value2<0.0) ;
2526                  double gap = columnGap[kcol]*value2 ;
2527                  bool doDown = (rowLo2 -gap > 0.0) ;
2528                  if (doDown && ((markIt&(tightenUpper|infLower)) == 0)) {
2529                    double dbound = colLower[kcol] + rowLo2/value2 ;
2530                    if (intVar[kcol]) {
2531                      markIt |= tightenUpper ;
2532                      newUpper = floor(dbound+primalTolerance_) ;
2533                    } else {
2534                      newUpper=dbound ;
2535                      if (newUpper-primalTolerance_<colLower[kcol]&&
2536                          newUpper+primalTolerance_>=colLower[kcol]) {
2537                        newUpper=colLower[kcol] ;
2538                        markIt |= tightenUpper ;
2539                        //markIt= (tightenLower|tightenUpper) ;
2540                      } else {
2541                        // avoid problems - fix later ?
2542                        markIt= (tightenLower|tightenUpper) ;
2543                      }
2544                    }
2545                    moveDown = newUpper-colUpper[kcol] ;
2546                    if (!moveDown)
2547                      continue ;
2548                    bool onList = ((markC[kcol]&(tightenLower|tightenUpper))!=0) ;
2549                    if (nstackC<2*maxStack) {
2550                      markC[kcol] = markIt ;
2551                    }
2552                    if (moveDown&&nstackC<2*maxStack) {
2553                      if (!onList) {
2554                        stackC[nstackC]=kcol ;
2555                        saveL[nstackC]=colLower[kcol] ;
2556                        saveU[nstackC]=colUpper[kcol] ;
2557                        assert (saveU[nstackC]>saveL[nstackC]) ;
2558                        assert (nstackC<nCols) ;
2559                        nstackC++ ;
2560                        onList=true ;
2561                      }
2562                      if (newUpper<colsol[kcol]) {
2563                        if (djs[kcol]>0.0) {
2564                          /* should be infeasible */
2565                          assert (colLower[kcol]>newUpper+primalTolerance_) ;
2566                        } else {
2567                          objVal += moveDown*djs[kcol] ;
2568                        }
2569                      }
2570                      if (intVar[kcol])
2571                        newUpper = CoinMin(colUpper[kcol],floor(newUpper+1.0e-4)) ;
2572                      colUpper[kcol]=newUpper ;
2573                      columnGap[kcol] = newUpper-colLower[kcol]-primalTolerance_ ;
2574                      if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
2575                        markC[kcol]= (tightenLower|tightenUpper); // say fixed
2576                      }
2577                      markC[kcol] &= ~(infLower|infUpper) ;
2578                      if (colUpper[kcol]>1.0e10)
2579                        markC[kcol] |= infUpper ;
2580                      if (colLower[kcol]<-1.0e10)
2581                        markC[kcol] |= infLower ;
2582
2583                      if (!updateRowBounds(kcol,moveDown,
2584                               columnStart,columnLength,row,columnElements,
2585                               rowUpper,rowLower,nstackR,stackR,markR,
2586                               minR,maxR,saveMin,saveMax)) {
2587                        colUpper[kcol] = -1.0e50 ;
2588                      }
2589                    }
2590                    if (colLower[kcol] > colUpper[kcol]+primalTolerance_) {
2591                      notFeasible = true ;
2592                      k = columnStart[jcol]+columnLength[jcol] ;
2593                      istackC = nstackC+1 ;
2594                      break ;
2595                    }
2596                  }
2597                }
2598              }  // end big loop rStart->rPos
2599            }
2600/*
2601  We've finished working on the negative coefficients of the row. Advance
2602  rStart and rEnd to cover the positive coefficients and repeat the previous
2603  500 lines.
2604*/
2605            rStart = rEnd ;
2606            rEnd = rowStart[irow+1] ;
2607            if (doRowUpP&&doRowLoP) {
2608              //doRowUpP=doRowLoP=false ;
2609              // Start pos values loop
2610              for (int kk =rStart;kk<rEnd;kk++) {
2611                int kcol=column[kk] ;
2612                int markIt=markC[kcol] ;
2613                if ((markIt&(tightenLower|tightenUpper))!= (tightenLower|tightenUpper)) {
2614                  double value2=rowElements[kk] ;
2615                  assert (value2 > 0.0) ;
2616                  /* positive element */
2617                  double gap = columnGap[kcol]*value2 ;
2618                  bool doDown = (rowLo2 + gap > 0.0) ;
2619                  bool doUp = (rowUp2 - gap < 0.0) ;
2620                  if (doDown||doUp) {
2621                    double moveUp=0.0 ;
2622                    double moveDown=0.0 ;
2623                    double newUpper=-1.0 ;
2624                    double newLower=1.0 ;
2625                    if (doDown && ((markIt&(tightenLower|infUpper)) == 0)) {
2626                      double dbound = colUpper[kcol] + rowLo2/value2 ;
2627                      if (intVar[kcol]) {
2628                        markIt |= tightenLower ;
2629                        newLower = ceil(dbound-primalTolerance_) ;
2630                      } else {
2631                        newLower=dbound ;
2632                        if (newLower+primalTolerance_>colUpper[kcol]&&
2633                            newLower-primalTolerance_<=colUpper[kcol]) {
2634                          newLower=colUpper[kcol] ;
2635                          markIt |= tightenLower ;
2636                          //markIt= (tightenLower|tightenUpper) ;
2637                        } else {
2638                          // avoid problems - fix later ?
2639                          markIt= (tightenLower|tightenUpper) ;
2640                        }
2641                      }
2642                      moveUp = newLower-colLower[kcol] ;
2643                    }
2644                    if (doUp && ((markIt&(tightenUpper|infLower)) == 0)) {
2645                      double dbound = colLower[kcol] + rowUp2/value2 ;
2646                      if (intVar[kcol]) {
2647                        markIt |= tightenUpper ;
2648                        newUpper = floor(dbound+primalTolerance_) ;
2649                      } else {
2650                        newUpper=dbound ;
2651                        if (newUpper-primalTolerance_<colLower[kcol]&&
2652                            newUpper+primalTolerance_>=colLower[kcol]) {
2653                          newUpper=colLower[kcol] ;
2654                          markIt |= tightenUpper ;
2655                          //markIt= (tightenLower|tightenUpper) ;
2656                        } else {
2657                          // avoid problems - fix later ?
2658                          markIt= (tightenLower|tightenUpper) ;
2659                        }
2660                      }
2661                      moveDown = newUpper-colUpper[kcol] ;
2662                    }
2663                    if (!moveUp&&!moveDown)
2664                      continue ;
2665                    bool onList = ((markC[kcol]&(tightenLower|tightenUpper))!=0) ;
2666                    if (nstackC<2*maxStack) {
2667                      markC[kcol] = markIt ;
2668                    }
2669                    if (moveUp&&nstackC<2*maxStack) {
2670                      if (!onList) {
2671                        stackC[nstackC]=kcol ;
2672                        saveL[nstackC]=colLower[kcol] ;
2673                        saveU[nstackC]=colUpper[kcol] ;
2674                        assert (saveU[nstackC]>saveL[nstackC]) ;
2675                        assert (nstackC<nCols) ;
2676                        nstackC++ ;
2677                        onList=true ;
2678                      }
2679                      if (newLower>colsol[kcol]) {
2680                        if (djs[kcol]<0.0) {
2681                          /* should be infeasible */
2682                          assert (newLower>colUpper[kcol]+primalTolerance_) ;
2683                        } else {
2684                          objVal += moveUp*djs[kcol] ;
2685                        }
2686                      }
2687                      if (intVar[kcol]) 
2688                        newLower = CoinMax(colLower[kcol],ceil(newLower-1.0e-4)) ;
2689                      colLower[kcol]=newLower ;
2690                      columnGap[kcol] = colUpper[kcol]-newLower-primalTolerance_ ;
2691                      if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
2692                        markC[kcol]= (tightenLower|tightenUpper); // say fixed
2693                      }
2694                      markC[kcol] &= ~(infLower|infUpper) ;
2695                      if (colUpper[kcol]>1.0e10)
2696                        markC[kcol] |= infUpper ;
2697                      if (colLower[kcol]<-1.0e10)
2698                        markC[kcol] |= infLower ;
2699
2700                      if (!updateRowBounds(kcol,moveUp,
2701                               columnStart,columnLength,row,columnElements,
2702                               rowUpper,rowLower,nstackR,stackR,markR,
2703                               minR,maxR,saveMin,saveMax)) {
2704                        colLower[kcol] = 1.0e50 ;
2705                      }
2706                    }
2707                    if (moveDown&&nstackC<2*maxStack) {
2708                      if (!onList) {
2709                        stackC[nstackC]=kcol ;
2710                        saveL[nstackC]=colLower[kcol] ;
2711                        saveU[nstackC]=colUpper[kcol] ;
2712                        assert (saveU[nstackC]>saveL[nstackC]) ;
2713                        assert (nstackC<nCols) ;
2714                        nstackC++ ;
2715                        onList=true ;
2716                      }
2717                      if (newUpper<colsol[kcol]) {
2718                        if (djs[kcol]>0.0) {
2719                          /* should be infeasible */
2720                          assert (colLower[kcol]>newUpper+primalTolerance_) ;
2721                        } else {
2722                          objVal += moveDown*djs[kcol] ;
2723                        }
2724                      }
2725                      if (intVar[kcol])
2726                        newUpper = CoinMin(colUpper[kcol],floor(newUpper+1.0e-4)) ;
2727                      colUpper[kcol]=newUpper ;
2728                      columnGap[kcol] = newUpper-colLower[kcol]-primalTolerance_ ;
2729                      if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
2730                        markC[kcol]= (tightenLower|tightenUpper); // say fixed
2731                      }
2732                      markC[kcol] &= ~(infLower|infUpper) ;
2733                      if (colUpper[kcol]>1.0e10)
2734                        markC[kcol] |= infUpper ;
2735                      if (colLower[kcol]<-1.0e10)
2736                        markC[kcol] |= infLower ;
2737
2738                      if (!updateRowBounds(kcol,moveDown,
2739                               columnStart,columnLength,row,columnElements,
2740                               rowUpper,rowLower,nstackR,stackR,markR,
2741                               minR,maxR,saveMin,saveMax)) {
2742                        colUpper[kcol] = -1.0e50 ;
2743                      }
2744                    }
2745                    if (colLower[kcol]>colUpper[kcol]+primalTolerance_) {
2746                      notFeasible = true; ;
2747                      k=columnStart[jcol]+columnLength[jcol] ;
2748                      istackC=nstackC+1 ;
2749                      break ;
2750                    }
2751                  }
2752                }
2753              } // end big loop rPos->rEnd
2754            } else if (doRowUpP) {
2755              // Start pos values loop
2756              for (int kk =rStart;kk<rEnd;kk++) {
2757                int kcol =column[kk] ;
2758                int markIt=markC[kcol] ;
2759                if ((markIt&(tightenLower|tightenUpper))!= (tightenLower|tightenUpper)) {
2760                  double value2=rowElements[kk] ;
2761                  assert (value2 > 0.0) ;
2762                  /* positive element */
2763                  double gap = columnGap[kcol]*value2 ;
2764                  bool doUp = (rowUp2 - gap < 0.0) ;
2765                  if (doUp && ((markIt&(tightenUpper|infLower)) == 0)) {
2766                    double newUpper=-1.0 ;
2767                    double dbound = colLower[kcol] + rowUp2/value2 ;
2768                    if (intVar[kcol]) {
2769                      markIt |= tightenUpper ;
2770                      newUpper = floor(dbound+primalTolerance_) ;
2771                    } else {
2772                      newUpper=dbound ;
2773                      if (newUpper-primalTolerance_<colLower[kcol]&&
2774                          newUpper+primalTolerance_>=colLower[kcol]) {
2775                        newUpper=colLower[kcol] ;
2776                        markIt |= tightenUpper ;
2777                        //markIt= (tightenLower|tightenUpper) ;
2778                      } else {
2779                        // avoid problems - fix later ?
2780                        markIt= (tightenLower|tightenUpper) ;
2781                      }
2782                    }
2783                    double moveDown = newUpper-colUpper[kcol] ;
2784                    if (!moveDown)
2785                      continue ;
2786                    bool onList = ((markC[kcol]&(tightenLower|tightenUpper))!=0) ;
2787                    if (nstackC<2*maxStack) {
2788                      markC[kcol] = markIt ;
2789                    }
2790                    if (nstackC<2*maxStack) {
2791                      if (!onList) {
2792                        stackC[nstackC]=kcol ;
2793                        saveL[nstackC]=colLower[kcol] ;
2794                        saveU[nstackC]=colUpper[kcol] ;
2795                        assert (saveU[nstackC]>saveL[nstackC]) ;
2796                        assert (nstackC<nCols) ;
2797                        nstackC++ ;
2798                        onList=true ;
2799                      }
2800                      if (newUpper<colsol[kcol]) {
2801                        if (djs[kcol]>0.0) {
2802                          /* should be infeasible */
2803                          assert (colLower[kcol]>newUpper+primalTolerance_) ;
2804                        } else {
2805                          objVal += moveDown*djs[kcol] ;
2806                        }
2807                      }
2808                      if (intVar[kcol])
2809                        newUpper = CoinMin(colUpper[kcol],floor(newUpper+1.0e-4)) ;
2810                      colUpper[kcol]=newUpper ;
2811                      columnGap[kcol] = newUpper-colLower[kcol]-primalTolerance_ ;
2812                      if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
2813                        markC[kcol]= (tightenLower|tightenUpper); // say fixed
2814                      }
2815                      markC[kcol] &= ~(infLower|infUpper) ;
2816                      if (colUpper[kcol]>1.0e10)
2817                        markC[kcol] |= infUpper ;
2818                      if (colLower[kcol]<-1.0e10)
2819                        markC[kcol] |= infLower ;
2820
2821                      if (!updateRowBounds(kcol,moveDown,
2822                               columnStart,columnLength,row,columnElements,
2823                               rowUpper,rowLower,nstackR,stackR,markR,
2824                               minR,maxR,saveMin,saveMax)) {
2825                        colUpper[kcol] = -1.0e50 ;
2826                      }
2827                    }
2828                    if (colLower[kcol]>colUpper[kcol]+primalTolerance_) {
2829                      notFeasible = true ;
2830                      k=columnStart[jcol]+columnLength[jcol] ;
2831                      istackC=nstackC+1 ;
2832                      break ;
2833                    }
2834                  }
2835                }
2836              } // end big loop rPos->rEnd
2837            } else if (doRowLoP) {
2838              // Start pos values loop
2839              for (int kk =rStart;kk<rEnd;kk++) {
2840                int kcol =column[kk] ;
2841                if ((markC[kcol]&(tightenLower|tightenUpper))!= (tightenLower|tightenUpper)) {
2842                  double value2=rowElements[kk] ;
2843                  int markIt=markC[kcol] ;
2844                  assert (value2 > 0.0) ;
2845                  /* positive element */
2846                  double gap = columnGap[kcol]*value2 ;
2847                  bool doDown = (rowLo2 +gap > 0.0) ;
2848                  if (doDown&&(markIt&(tightenLower|infUpper))==0) {
2849                    double newLower=1.0 ;
2850                    double dbound = colUpper[kcol] + rowLo2/value2 ;
2851                    if (intVar[kcol]) {
2852                      markIt |= tightenLower ;
2853                      newLower = ceil(dbound-primalTolerance_) ;
2854                    } else {
2855                      newLower=dbound ;
2856                      if (newLower+primalTolerance_>colUpper[kcol]&&
2857                          newLower-primalTolerance_<=colUpper[kcol]) {
2858                        newLower=colUpper[kcol] ;
2859                        markIt |= tightenLower ;
2860                        //markIt= (tightenLower|tightenUpper) ;
2861                      } else {
2862                        // avoid problems - fix later ?
2863                        markIt= (tightenLower|tightenUpper) ;
2864                        }
2865                    }
2866                    double moveUp = newLower-colLower[kcol] ;
2867                    if (!moveUp)
2868                      continue ;
2869                    bool onList = ((markC[kcol]&(tightenLower|tightenUpper))!=0) ;
2870                    if (nstackC<2*maxStack) {
2871                      markC[kcol] = markIt ;
2872                    }
2873                    if (nstackC<2*maxStack) {
2874                      if (!onList) {
2875                        stackC[nstackC]=kcol ;
2876                        saveL[nstackC]=colLower[kcol] ;
2877                        saveU[nstackC]=colUpper[kcol] ;
2878                        assert (saveU[nstackC]>saveL[nstackC]) ;
2879                        assert (nstackC<nCols) ;
2880                        nstackC++ ;
2881                        onList=true ;
2882                      }
2883                      if (newLower>colsol[kcol]) {
2884                        if (djs[kcol]<0.0) {
2885                          /* should be infeasible */
2886                          assert (newLower>colUpper[kcol]+primalTolerance_) ;
2887                        } else {
2888                          objVal += moveUp*djs[kcol] ;
2889                        }
2890                      }
2891                      if (intVar[kcol]) 
2892                        newLower = CoinMax(colLower[kcol],ceil(newLower-1.0e-4)) ;
2893                      colLower[kcol]=newLower ;
2894                      columnGap[kcol] = colUpper[kcol]-newLower-primalTolerance_ ;
2895                      if (fabs(colUpper[kcol]-colLower[kcol])<1.0e-6) {
2896                        markC[kcol]= (tightenLower|tightenUpper); // say fixed
2897                      }
2898                      markC[kcol] &= ~(infLower|infUpper) ;
2899                      if (colUpper[kcol]>1.0e10)
2900                        markC[kcol] |= infUpper ;
2901                      if (colLower[kcol]<-1.0e10)
2902                        markC[kcol] |= infLower ;
2903
2904                      if (!updateRowBounds(kcol,moveUp,
2905                               columnStart,columnLength,row,columnElements,
2906                               rowUpper,rowLower,nstackR,stackR,markR,
2907                               minR,maxR,saveMin,saveMax)) {
2908                        colLower[kcol] = 1.0e50 ;
2909                      }
2910                    }
2911                    if (colLower[kcol]>colUpper[kcol]+primalTolerance_) {
2912                      notFeasible = true ;
2913                      k=columnStart[jcol]+columnLength[jcol] ;
2914                      istackC=nstackC+1 ;
2915                      break ;
2916                    }
2917                  }
2918                }
2919              }      // end big loop rPos->rEnd
2920            }    // end processing of positive coefficients of row.
2921          }   // end loop to walk the column of a variable popped off stackC
2922          istackC++ ;
2923        }  // end stackC processing loop
2924/*
2925  PROBE LOOP: END PROPAGATION
2926
2927  End propagation of probe in one direction for a single variable. Hard to
2928  believe 1,000 lines of code can achieve so little.
2929*/
2930/*
2931  PROBE LOOP: INFEASIBILITY
2932
2933  Feasibility check. Primal infeasibility is noted in notFeasible; we still
2934  need to test against the objective bound. If we have shown up and down
2935  infeasibility, we're infeasible, period. Break from the probe loop and
2936  terminate the lookedAt and pass loops by forcing their iteration counts past
2937  the maximum.
2938
2939  If we're not feasible, then we surely can't drive the variable to bound,
2940  so reset goingToTrueBound.
2941
2942  TODO: I'm coming to think forcing goingToTrueBound = 0 is really an overload
2943        to suppress various post-probe activities (cut generation, singleton
2944        motion, etc) that should not be performed when the probe shows
2945        infeasible. It might be that this is not the best approach.
2946        -- lh, 101216 --
2947*/
2948        if (notFeasible || (objVal > cutoff)) {
2949#         if CGL_DEBUG > 1
2950          std::cout << "  Column " << j << " infeasible: " ;
2951          if (notFeasible && (objVal > cutoff))
2952            std::cout << "primal bounds and objective cutoff" ;
2953          else if (notFeasible)
2954            std::cout << "primal bounds" ;
2955          else
2956            std::cout << "objective cutoff" ;
2957          std::cout << "." << std::endl ;
2958#         endif
2959          notFeasible = true ;
2960          if (iway == upProbe && feasRecord == 0) {
2961            ninfeas = 1 ;
2962            j = nCols-1 ;
2963            iLook = numberThisTime_ ;
2964            ipass = maxPass ;
2965            break ;
2966          }
2967          goingToTrueBound = 0 ;
2968        } else {
2969          feasRecord |= feasValue[iway] ; 
2970        }
2971/*
2972  Save the implications? Currently restricted to binary variables. If info
2973  were a CglTreeProbing object, fixes() would save the implication in the
2974  arrays kept there, returning false only when space is exceeded.
2975 
2976  TODO: This would be one of the places that will fail if fixes() and
2977        initializeFixing() are removed from CglTreeInfo.  -- lh, 101127 --
2978*/
2979        if (!notFeasible && saveFixingInfo) {
2980          // save fixing info
2981          assert (j == stackC[0]) ;
2982          int toValue = (way[iway] == probeDown) ? -1 : +1 ;
2983          for (istackC = 1 ; istackC < nstackC ; istackC++) {
2984            int icol = stackC[istackC] ;
2985            // for now back to just 0-1
2986            if (colUpper[icol]-colLower[icol]<1.0e-12 &&
2987                !saveL[istackC] && saveU[istackC]==1.0) {
2988              assert(saveL[istackC] == colLower[icol] ||
2989                     saveU[istackC] == colUpper[icol]) ;
2990              saveFixingInfo =
2991                  info->fixes(j,toValue,icol,
2992                              (colLower[icol] == saveL[istackC])) ;
2993            }
2994          }
2995        }
2996/*
2997  PROBE LOOP: MONOTONE (BREAK AND KEEP)
2998
2999  We can't prove infeasibility. What's left? If we're at the end of the up
3000  probe, we may know enough to prove monotonicity:
3001
3002  * If this is the up probe and we were feasible on this probe but not on
3003    the down pass, we're monotone.
3004  * If this is the second down probe, we were feasible on down, infeasible on
3005    up (monotone) and repeated the down probe to recreate the bounds.
3006
3007  Either way, execute monotoneActions to capture the state and move on to the
3008  next variable. Terminate the probe loop by forcing iway = oneProbeTooMany.
3009
3010  This if-then-else extends to the end of the down/up/down probe loop.
3011
3012  TODO: There's something not right about this `second down pass' business.
3013        Why do we copy stackC to stackC0 on a feasible first pass? Surely that
3014        could be used to restore the state of the first down pass.
3015        -- lh, 101128 --
3016*/
3017        if (iway == downProbe2 ||
3018            (iway == upProbe && feasRecord == upProbeFeas)) {
3019          nfixed++ ;
3020          bool retVal = monotoneActions(primalTolerance_,si,cs,
3021                                        nstackC,stackC,intVar,
3022                                        colLower,colsol,colUpper,
3023                                        index,element) ;
3024          if (retVal) anyColumnCuts = true ;
3025          clearStacks(primalTolerance_,nstackC,stackC,markC,colLower,colUpper,
3026                      nstackR,stackR,markR) ;
3027          break ;
3028        }
3029/*
3030  PROBE LOOP: ITERATE
3031
3032  If we reach here, we may iterate the probe loop.
3033
3034  Cases remaining include
3035
3036  a) End of the first down probe; in this case we will iterate regardless of
3037     feasibility.
3038  b) End of the up probe, up and down feasible; in this case we will not
3039     iterate (we're done, cannot prove monotonicity).
3040  c) End of the up probe, down feasible, up infeasible; in this case we will
3041     iterate to restore the down feasible state.
3042
3043  TODO: Unless I miss my guess, large chunks of this block of code will be
3044        replicated in the code that handles case b), as we generate cuts
3045        related to forcing the variable up.  Further inspection confirms
3046        this. -- lh, 101128 --
3047
3048  TODO: I'm becoming more convinced that the second down probe iteration is
3049        unnecessary. Given that we have stackC0, lo0, and up0, seems like
3050        calling monotoneActions with lo0 = colLower, up0 = colUpper, and
3051        stackC0 = stackC should do the trick.  -- lh, 101216 --
3052
3053  TODO: As I get the control flow sorted out a bit better, I should try to
3054        separate b) (which does not iterate) from a) and c).
3055
3056  The next block deals with the end of the first down probe. If the probe
3057  is feasible, attempt to tighten column bounds with singletonRows, process
3058  stackC to generate disaggregation cuts, and copy the tightened bounds
3059  to stackC0, lo0, and up0 for use after the up probe.  Feasible or not,
3060  we'll need to go on to the up probe, so restore the original bounds.
3061
3062  Note that GoingToTrueBound was set to 0 way back in initial setup unless
3063  the user requested disaggregation cuts (rowCuts&0x01).
3064
3065  I've rearranged the control flow here, moving the call to singletonRows
3066  ahead of the code that captures column bound changes to stackC0, et al.
3067  This makes it possible to use the same code at the end of case b) (down and
3068  up probe both feasible). Also, singletonRows is just another way to tighten
3069  column bounds. In its original position, it was subject to the same
3070  restrictions as coefficient strengthening (coefficient strengthening enabled
3071  and goingToTrueBound for the probing variable). I don't see any reason to
3072  keep them (but it might be good to add a separate control bit precisely for
3073  singletonRows).  -- lh, 110113 --
3074*/
3075        if (iway == downProbe) {
3076          if (!notFeasible) {
3077/*
3078  Attempt to tighten singleton bounds and generate disaggregation cuts.
3079*/
3080            singletonRows(j,primalTolerance_,si,rowCopy,markC,
3081                          nstackC,stackC,saveL,saveU,
3082                          colUpper,colLower,columnGap,
3083                          nstackR,stackR,rowUpper,rowLower,maxR,minR) ;
3084            if (goingToTrueBound == 2 && !justReplace) {
3085              disaggCuts(nstackC,probeDown,primalTolerance_,
3086                         disaggEffectiveness,si,rowCut,stackC,colsol,
3087                         colUpper,colLower,saveU,saveL,index,element) ;
3088            }
3089/*
3090  Copy off the bound changes from the down probe -- we'll need them after the
3091  up probe.
3092*/
3093            nstackC0 = CoinMin(nstackC,maxStack) ;
3094            for (istackC = 0 ; istackC < nstackC0 ; istackC++) {
3095              int icol = stackC[istackC] ;
3096              stackC0[istackC] = icol ;
3097              lo0[istackC] = colLower[icol] ;
3098              up0[istackC] = colUpper[icol] ;
3099            }
3100          } else {
3101            nstackC0 = 0 ;
3102          }
3103/*
3104  Now reset column bounds to their original values in preparation for the
3105  up probe.
3106*/
3107          for (istackC = nstackC-1 ; istackC >= 0 ; istackC--) {
3108            int icol = stackC[istackC] ;
3109            double oldU = saveU[istackC] ;
3110            double oldL = saveL[istackC] ;
3111            colUpper[icol] = oldU ;
3112            colLower[icol] = oldL ;
3113            columnGap[icol] = oldU-oldL-primalTolerance_ ;
3114            markC[icol] = 0 ;
3115            if (oldU > 1.0e10)
3116              markC[icol] |= infUpper ;
3117            if (oldL < -1.0e10)
3118              markC[icol] |= infLower ;
3119          }
3120/*
3121  And now for the rows. Remember, we're still finishing off the first down
3122  probe and readying for the next iteration (up probe). The final three lines
3123  (of about 300) actually restore the row bounds. The rest is cut generation.
3124*/
3125          if ((rowCuts&0x02) != 0 && goingToTrueBound && !anyColumnCuts)
3126            strengthenCoeff(j,primalTolerance_,justReplace,canReplace,
3127                            needEffectiveness,si,rowCut,
3128                            rowCopy,colUpper,colLower,colsol,nstackR,stackR,
3129                            rowUpper,rowLower,maxR,minR,realRows,
3130                            element,index,info) ;
3131/*
3132  Restore row bounds.
3133*/
3134          for (istackR = 0 ; istackR < nstackR ; istackR++) {
3135            int irow = stackR[istackR] ;
3136            minR[irow] = saveMin[istackR] ;
3137            maxR[irow] = saveMax[istackR] ;
3138            markR[irow] = -1 ;
3139          }
3140        }    // end of code to iterate after first down pass
3141/*
3142  PROBE LOOP: DOWN AND UP FEASIBLE
3143
3144  Again, not quite the full story. Cases remaining:
3145  b) End of up probe, up and down feasible; in this case we will not
3146     iterate.
3147  c) End of up probe, down feasible, up infeasible; in this case we will
3148     iterate to restore the down feasible state.
3149
3150  The next block deals with case b). We don't want to iterate the
3151  down/up/down loop, so force iway to oneProbeTooMany. As usual, the move
3152  singleton code consists of two nearly identical blocks, one working off
3153  U(i), the next off L(i).
3154
3155  TODO: The conditions guarding the move singleton code here are different than
3156        those for pass 0: here there's an added guard using strengthenRow.
3157        Just a failure to change all instances? Nor do I see how John's
3158        comment (`just at root') is enforced here, unless it's from outside
3159        via rowCuts or strengthenRow.
3160
3161        The check for any column cuts is not present, so presumably we can
3162        do singleton motion in the presence of column cuts. The check for gap
3163        does not have the high end (gap < 1.0e8) test.
3164
3165        Note also that here the test for cut generation is moved outside the
3166        loop that iterates over rows. Presumably that's because in the
3167        previous case, the loop also restored markR, etc.
3168
3169        Note that the code that handles the bound change is considerably
3170        different than the previous instance in pass 0. We're adding the
3171        changes to stackC rather than stackC0.
3172        -- lh, 101128 --
3173
3174  Moved the call to disaggCuts up from the common code for b) and c), because
3175  c) (up infeasible) sets goingToTrueBound to 0 and we can't execute
3176  disaggCuts. Similarly for strengthenCoeff.
3177*/
3178          else {
3179          if (iway == upProbe &&
3180              feasRecord == (downProbeFeas|upProbeFeas)) {
3181            iway = oneProbeTooMany ;
3182
3183            singletonRows(j,primalTolerance_,si,rowCopy,markC,
3184                          nstackC,stackC,saveL,saveU,
3185                          colUpper,colLower,columnGap,
3186                          nstackR,stackR,rowUpper,rowLower,maxR,minR) ;
3187            if (goingToTrueBound == 2 && !justReplace) {
3188              disaggCuts(nstackC,probeUp,primalTolerance_,
3189                         disaggEffectiveness,si,
3190                         rowCut,stackC,colsol,colUpper,colLower,saveU,saveL,
3191                         index,element) ;
3192            if ((rowCuts&0x02) != 0 && goingToTrueBound && !anyColumnCuts)
3193              strengthenCoeff(j,primalTolerance_,justReplace,canReplace,
3194                              needEffectiveness,si,rowCut,
3195                              rowCopy,colUpper,colLower,colsol,nstackR,stackR,
3196                              rowUpper,rowLower,maxR,minR,realRows,
3197                              element,index,info) ;
3198            }
3199/*
3200  jjf: point back to stack
3201
3202  We're done playing with singletons. Get to the real work here. We don't need
3203  markC any more; coopt it for a cross-reference, var index -> stackC index.
3204  The +1000 offset allows us to distinguish invalid entries (not all variables
3205  are on stackC).
3206*/
3207            for (istackC = nstackC-1 ; istackC >= 0 ; istackC--) {
3208              int icol = stackC[istackC] ;
3209              markC[icol] = istackC+1000 ;
3210            }
3211            OsiColCut cc ;
3212            int nTot = 0, nFix = 0, nInt = 0 ;
3213            bool ifCut = false ;
3214/*
3215  See if we have bounds improvement. Given a lower bound ld<j> from the
3216  down probe, a bound lu<j> from the up probe, and the original bound lo<j>,
3217  the bound we want is max(min(ld<j>,lu<j>),lo<j>). Use a conservative
3218  tolerance.
3219*/
3220            for (istackC = 1 ; istackC < nstackC0 ; istackC++) {
3221              int icol = stackC0[istackC] ;
3222              int istackC1 = markC[icol]-1000 ;
3223              if (istackC1 >= 0) {
3224                if (CoinMin(lo0[istackC],colLower[icol]) >
3225                                              saveL[istackC1]+1.0e-4) {
3226                  saveL[istackC1] = CoinMin(lo0[istackC],colLower[icol]) ;
3227                  if (intVar[icol] /*||!info->inTree*/ ) {
3228                    element[nFix] = saveL[istackC1] ;
3229                    index[nFix++] = icol ;
3230                    nInt++ ;
3231                    if (colsol[icol] < saveL[istackC1]-primalTolerance_)
3232                      ifCut = true ;
3233                  }
3234                  nfixed++ ;
3235                } 
3236              }
3237            }
3238            if (nFix) {
3239              nTot = nFix ;
3240              cc.setLbs(nFix,index,element) ;
3241              nFix = 0 ;
3242            }
3243/*
3244  Repeat for upper bounds. There's an odd bit of asymmetry here; this loop
3245  tries to generate a cut if there's no bound improvement.
3246
3247  TODO: What, pray tell, is a `two cut'?  Binary variables only, it seems.
3248        Judging from the conditions, we're looking for a variable that's
3249        fixed at 0 on a down probe and fixed at 1 on an up probe, or vice
3250        versa.
3251        -- lh, 101128 --
3252*/
3253            for (istackC = 1 ; istackC < nstackC0 ; istackC++) {
3254              int icol = stackC0[istackC] ;
3255              int istackC1 = markC[icol]-1000 ;
3256              if (istackC1 >= 0) {
3257                if (CoinMax(up0[istackC],colUpper[icol]) <
3258                                              saveU[istackC1]-1.0e-4) {
3259                  saveU[istackC1] = CoinMax(up0[istackC],colUpper[icol]) ;
3260                  if (intVar[icol] /*||!info->inTree*/ ) {
3261                    element[nFix] = saveU[istackC1] ;
3262                    index[nFix++] = icol ;
3263                    nInt++ ;
3264                    if (colsol[icol] > saveU[istackC1]+primalTolerance_)
3265                      ifCut = true ;
3266                  }
3267                  nfixed++ ;
3268                } else if (!info->inTree &&
3269                           saveL[0] == 0.0 && saveU[0] == 1.0) {
3270                  // See if can do two cut
3271                  double upperWhenDown = up0[istackC] ;
3272                  double lowerWhenDown = lo0[istackC] ;
3273                  double upperWhenUp = colUpper[icol] ;
3274                  double lowerWhenUp = colLower[icol] ;
3275                  double upperOriginal = saveU[istackC1] ;
3276                  double lowerOriginal = saveL[istackC1] ;
3277                  if (upperWhenDown < lowerOriginal+1.0e-12 &&
3278                      lowerWhenUp > upperOriginal-1.0e-12) {
3279                    OsiRowCut rc ;
3280                    rc.setLb(lowerOriginal) ;
3281                    rc.setUb(lowerOriginal) ;
3282                    rc.setEffectiveness(1.0e-5) ;
3283                    int index[2] ;
3284                    double element[2] ;
3285                    index[0] = j ;
3286                    index[1] = icol ;
3287                    element[0] = -(upperOriginal-lowerOriginal) ;
3288/*
3289  jjf: If zero then - must have been fixed without noticing!
3290
3291  TODO: Uh, fixed without noticing?! That's an algorithm problem.
3292        -- lh, 101128 --
3293*/
3294                    if (fabs(element[0]) > 1.0e-8) {
3295                      element[1] = 1.0 ;
3296                      rc.setRow(2,index,element,false) ;
3297                      cs.insert(rc) ;
3298                    }
3299                  } else if (upperWhenUp < lowerOriginal+1.0e-12 &&
3300                             lowerWhenDown > upperOriginal-1.0e-12) {
3301                    OsiRowCut rc ;
3302                    rc.setLb(upperOriginal) ;
3303                    rc.setUb(upperOriginal) ;
3304                    rc.setEffectiveness(1.0e-5) ;
3305                    int index[2] ;
3306                    double element[2] ;
3307                    index[0] = j ;
3308                    index[1] = icol ;
3309                    element[0] = upperOriginal-lowerOriginal ;
3310                    element[1] = 1.0 ;
3311                    rc.setRow(2,index,element,false) ;
3312                    cs.insert(rc) ;
3313                  } 
3314                }
3315              }
3316            }    // end loop to update upper bounds (w/ two cuts)
3317            if (nFix) {
3318              nTot+=nFix ;
3319              cc.setUbs(nFix,index,element) ;
3320            }
3321            // could tighten continuous as well
3322            if (nInt) {
3323              if (ifCut) {
3324                cc.setEffectiveness(100.0) ;
3325              } else {
3326                cc.setEffectiveness(1.0e-5) ;
3327              }
3328#if CGL_DEBUG > 0
3329              checkBounds(si,cc) ;
3330#endif
3331              cs.insert(cc) ;
3332            }
3333          }    // end of code to handle up and down feasible.
3334/*
3335  PROBE LOOP: DOWN FEASIBLE, UP INFEASIBLE
3336
3337  The only remaining case is down feasible, up infeasible, in which case we
3338  need to reset to the initial state and run the down probe again. Surely
3339  there must be a better way?
3340*/
3341            else {
3342            goingToTrueBound = 0 ;
3343          }
3344/*
3345  PROBE LOOP: RESTORE
3346
3347  And now the code to reset to the initial state.
3348*/
3349          /* restore all */
3350          for (istackC = nstackC-1 ; istackC >= 0 ; istackC--) {
3351            int icol = stackC[istackC] ;
3352            double oldU = saveU[istackC] ;
3353            double oldL = saveL[istackC] ;
3354/*
3355  TODO: This next bit differs in subtle ways from the restore at the end
3356        of the down probe. Here we force markC to show fixed if the bounds are
3357        within 1.0e-4 of one another, a fairly broad tolerance; code for
3358        the pass 0 restore does not restore fixed status.  Also,
3359        we're not distinguishing here between actions for down feasible,
3360        up feasible, vs. down feasible, up infeasible.  -- lh, 101128 --
3361
3362        I still don't see any reason for the restore here to differ from
3363        the restore after the down probe.  -- lh, 110113 --
3364*/
3365            colUpper[icol] = oldU ;
3366            colLower[icol] = oldL ;
3367            columnGap[icol] = oldU-oldL-primalTolerance_ ;
3368            if (oldU > oldL+1.0e-4) {
3369              markC[icol] = 0 ;
3370              if (oldU > 1.0e10)
3371                markC[icol] |= infUpper ;
3372              if (oldL < -1.0e10)
3373                markC[icol] |= infLower ;
3374            } else {
3375              markC[icol] = (tightenUpper|tightenLower) ;
3376            }
3377          }
3378/*
3379  End of column restore. On to rows.
3380*/
3381          for (istackR = 0 ; istackR < nstackR ; istackR++) {
3382            int irow = stackR[istackR] ;
3383            minR[irow] = saveMin[istackR] ;
3384            maxR[irow] = saveMax[istackR] ;
3385            markR[irow] = -1 ;
3386          }
3387        }       // end of PROBE: DOWN AND UP FEASIBLE
3388      }     // PROBE LOOP: END
3389    }    // LOOKEDAT LOOP: END
3390  }   // PASS LOOP: END
3391/*
3392  Free up the space we've been using.
3393*/
3394#ifndef ONE_ARRAY
3395  delete [] stackC0 ;
3396  delete [] lo0 ;
3397  delete [] up0 ;
3398  delete [] columnGap ;
3399  delete [] markC ;
3400  delete [] stackC ;
3401  delete [] stackR ;
3402  delete [] saveL ;
3403  delete [] saveU ;
3404  delete [] saveMin ;
3405  delete [] saveMax ;
3406  delete [] index ;
3407  delete [] element ;
3408  delete [] djs ;
3409  delete [] largestPositiveInRow ;
3410  delete [] largestNegativeInRow ;
3411#endif
3412  delete [] colsol ;
3413/*
3414  jjf:  Add in row cuts
3415
3416  Transfer row cuts from the local container into something that can escape
3417  this method. If we're not in `just replace' mode, simply copy them over to
3418  the container (cs) passed as a parameter using addCuts. If we're in `just
3419  replace' mode, ? need to explore further ?  -- lh, 101125 --
3420
3421  TODO: If there's any possibility of realRows[i] < 0, this code will break
3422        trying to read strengthenRow[]!
3423        -- lh, 101128 --
3424*/
3425  if (!ninfeas) {
3426    if (!justReplace) {
3427      rowCut.addCuts(cs,info->strengthenRow,info->pass) ;
3428    } else {
3429      for (int i = 0 ; i < nRows ; i++) {
3430        int realRow = realRows[i] ;
3431        OsiRowCut *cut = info->strengthenRow[realRow] ;
3432        if (realRow >= 0 && cut) {
3433#ifdef CLP_INVESTIGATE
3434          printf("Row %d, real row %d effectiveness %g\n",
3435                 i,realRow,cut->effectiveness()) ;
3436#endif
3437          cs.insert(cut) ;
3438        }
3439      }
3440    }
3441  }
3442#if 0
3443/*
3444  Debug print.
3445*/
3446  {
3447    int numberRowCutsAfter = cs.sizeRowCuts() ;
3448    int k ;
3449    for (k = 0;k<numberRowCutsAfter;k++) {
3450      OsiRowCut thisCut = cs.rowCut(k) ;
3451      printf("Cut %d is %g <=",k,thisCut.lb()) ;
3452      int n=thisCut.row().getNumElements() ;
3453      const int * column = thisCut.row().getIndices() ;
3454      const double * element = thisCut.row().getElements() ;
3455      assert (n) ;
3456      for (int i=0;i<n;i++) {
3457        printf(" %g*x%d",element[i],column[i]) ;
3458      }
3459      printf(" <= %g\n",thisCut.ub()) ;
3460    }
3461  }
3462#endif
3463
3464# if CGL_DEBUG > 0
3465  std::cout
3466    << "Leaving CglProbing::probe, ninfeas " << ninfeas
3467    << "." << std::endl ;
3468# endif
3469
3470  return (ninfeas) ;
3471}
Note: See TracBrowser for help on using the repository browser.