source: trunk/Clp/src/Idiot.cpp @ 1525

Last change on this file since 1525 was 1525, checked in by mjs, 10 years ago

Formatted .cpp, .hpp, .c, .h files with "astyle -A4 -p". This matches the formatting used in the grand CBC reorganization.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 79.0 KB
Line 
1/* $Id: Idiot.cpp 1525 2010-02-26 17:27:59Z mjs $ */
2// Copyright (C) 2002, International Business Machines
3// Corporation and others.  All Rights Reserved.
4
5#include "CoinPragma.hpp"
6#include <stdio.h>
7#include <stdarg.h>
8#include <stdlib.h>
9#include <math.h>
10#include "ClpPresolve.hpp"
11#include "Idiot.hpp"
12#include "CoinTime.hpp"
13#include "CoinSort.hpp"
14#include "CoinMessageHandler.hpp"
15#include "CoinHelperFunctions.hpp"
16// Redefine stuff for Clp
17#ifndef OSI_IDIOT
18#include "ClpMessage.hpp"
19#define OsiObjOffset ClpObjOffset
20#endif
21/**** strategy 4 - drop, exitDrop and djTolerance all relative:
22For first two major iterations these are small.  Then:
23
24drop - exit a major iteration if drop over 5*checkFrequency < this is
25used as info->drop*(10.0+fabs(last weighted objective))
26
27exitDrop - exit idiot if feasible and drop < this is
28used as info->exitDrop*(10.0+fabs(last objective))
29
30djExit - exit a major iteration if largest dj (averaged over 5 checks)
31drops below this - used as info->djTolerance*(10.0+fabs(last weighted objective)
32
33djFlag - mostly skip variables with bad dj worse than this => 2*djExit
34
35djTol - only look at variables with dj better than this => 0.01*djExit
36****************/
37
38#define IDIOT_FIX_TOLERANCE 1e-6
39#define SMALL_IDIOT_FIX_TOLERANCE 1e-10
40int
41Idiot::dropping(IdiotResult result,
42                double tolerance,
43                double small,
44                int *nbad)
45{
46     if (result.infeas <= small) {
47          double value = CoinMax(fabs(result.objval), fabs(result.dropThis)) + 1.0;
48          if (result.dropThis > tolerance * value) {
49               *nbad = 0;
50               return 1;
51          } else {
52               (*nbad)++;
53               if (*nbad > 4) {
54                    return 0;
55               } else {
56                    return 1;
57               }
58          }
59     } else {
60          *nbad = 0;
61          return 1;
62     }
63}
64// Deals with whenUsed and slacks
65int
66Idiot::cleanIteration(int iteration, int ordinaryStart, int ordinaryEnd,
67                      double * colsol, const double * lower, const double * upper,
68                      const double * rowLower, const double * rowUpper,
69                      const double * cost, const double * element, double fixTolerance,
70                      double & objValue, double & infValue)
71{
72     int n = 0;
73     if ((strategy_ & 16384) == 0) {
74          for (int i = ordinaryStart; i < ordinaryEnd; i++) {
75               if (colsol[i] > lower[i] + fixTolerance) {
76                    if (colsol[i] < upper[i] - fixTolerance) {
77                         n++;
78                    } else {
79                         colsol[i] = upper[i];
80                    }
81                    whenUsed_[i] = iteration;
82               } else {
83                    colsol[i] = lower[i];
84               }
85          }
86          return n;
87     } else {
88#ifdef COIN_DEVELOP
89          printf("entering inf %g, obj %g\n", infValue, objValue);
90#endif
91          int nrows = model_->getNumRows();
92          int ncols = model_->getNumCols();
93          int * posSlack = whenUsed_ + ncols;
94          int * negSlack = posSlack + nrows;
95          int * nextSlack = negSlack + nrows;
96          double * rowsol = reinterpret_cast<double *> (nextSlack + ncols);
97          memset(rowsol, 0, nrows * sizeof(double));
98#ifdef OSI_IDIOT
99          const CoinPackedMatrix * matrix = model_->getMatrixByCol();
100#else
101          ClpMatrixBase * matrix = model_->clpMatrix();
102#endif
103          const int * row = matrix->getIndices();
104          const CoinBigIndex * columnStart = matrix->getVectorStarts();
105          const int * columnLength = matrix->getVectorLengths();
106          //const double * element = matrix->getElements();
107          int i;
108          objValue = 0.0;
109          infValue = 0.0;
110          for ( i = 0; i < ncols; i++) {
111               if (nextSlack[i] == -1) {
112                    // not a slack
113                    if (colsol[i] > lower[i] + fixTolerance) {
114                         if (colsol[i] < upper[i] - fixTolerance) {
115                              n++;
116                              whenUsed_[i] = iteration;
117                         } else {
118                              colsol[i] = upper[i];
119                         }
120                         whenUsed_[i] = iteration;
121                    } else {
122                         colsol[i] = lower[i];
123                    }
124                    double value = colsol[i];
125                    if (value) {
126                         objValue += cost[i] * value;
127                         CoinBigIndex j;
128                         for (j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
129                              int iRow = row[j];
130                              rowsol[iRow] += value * element[j];
131                         }
132                    }
133               }
134          }
135          // temp fix for infinite lbs - just limit to -1000
136          for (i = 0; i < nrows; i++) {
137               double rowSave = rowsol[i];
138               int iCol;
139               iCol = posSlack[i];
140               if (iCol >= 0) {
141                    // slide all slack down
142                    double rowValue = rowsol[i];
143                    CoinBigIndex j = columnStart[iCol];
144                    double lowerValue = CoinMax(CoinMin(colsol[iCol], 0.0) - 1000.0, lower[iCol]);
145                    rowSave += (colsol[iCol] - lowerValue) * element[j];
146                    colsol[iCol] = lowerValue;
147                    while (nextSlack[iCol] >= 0) {
148                         iCol = nextSlack[iCol];
149                         double lowerValue = CoinMax(CoinMin(colsol[iCol], 0.0) - 1000.0, lower[iCol]);
150                         j = columnStart[iCol];
151                         rowSave += (colsol[iCol] - lowerValue) * element[j];
152                         colsol[iCol] = lowerValue;
153                    }
154                    iCol = posSlack[i];
155                    while (rowValue < rowLower[i] && iCol >= 0) {
156                         // want to increase
157                         double distance = rowLower[i] - rowValue;
158                         double value = element[columnStart[iCol]];
159                         double thisCost = cost[iCol];
160                         if (distance <= value*(upper[iCol] - colsol[iCol])) {
161                              // can get there
162                              double movement = distance / value;
163                              objValue += movement * thisCost;
164                              rowValue = rowLower[i];
165                              colsol[iCol] += movement;
166                         } else {
167                              // can't get there
168                              double movement = upper[iCol] - colsol[iCol];
169                              objValue += movement * thisCost;
170                              rowValue += movement * value;
171                              colsol[iCol] = upper[iCol];
172                              iCol = nextSlack[iCol];
173                         }
174                    }
175                    if (iCol >= 0) {
176                         // may want to carry on - because of cost?
177                         while (iCol >= 0 && cost[iCol] < 0 && rowValue < rowUpper[i]) {
178                              // want to increase
179                              double distance = rowUpper[i] - rowValue;
180                              double value = element[columnStart[iCol]];
181                              double thisCost = cost[iCol];
182                              if (distance <= value*(upper[iCol] - colsol[iCol])) {
183                                   // can get there
184                                   double movement = distance / value;
185                                   objValue += movement * thisCost;
186                                   rowValue = rowUpper[i];
187                                   colsol[iCol] += movement;
188                                   iCol = -1;
189                              } else {
190                                   // can't get there
191                                   double movement = upper[iCol] - colsol[iCol];
192                                   objValue += movement * thisCost;
193                                   rowValue += movement * value;
194                                   colsol[iCol] = upper[iCol];
195                                   iCol = nextSlack[iCol];
196                              }
197                         }
198                         if (iCol >= 0 && colsol[iCol] > lower[iCol] + fixTolerance &&
199                                   colsol[iCol] < upper[iCol] - fixTolerance) {
200                              whenUsed_[i] = iteration;
201                              n++;
202                         }
203                    }
204                    rowsol[i] = rowValue;
205               }
206               iCol = negSlack[i];
207               if (iCol >= 0) {
208                    // slide all slack down
209                    double rowValue = rowsol[i];
210                    CoinBigIndex j = columnStart[iCol];
211                    rowSave += (colsol[iCol] - lower[iCol]) * element[j];
212                    colsol[iCol] = lower[iCol];
213                    assert (lower[iCol] > -1.0e20);
214                    while (nextSlack[iCol] >= 0) {
215                         iCol = nextSlack[iCol];
216                         j = columnStart[iCol];
217                         rowSave += (colsol[iCol] - lower[iCol]) * element[j];
218                         colsol[iCol] = lower[iCol];
219                    }
220                    iCol = negSlack[i];
221                    while (rowValue > rowUpper[i] && iCol >= 0) {
222                         // want to increase
223                         double distance = -(rowUpper[i] - rowValue);
224                         double value = -element[columnStart[iCol]];
225                         double thisCost = cost[iCol];
226                         if (distance <= value*(upper[iCol] - lower[iCol])) {
227                              // can get there
228                              double movement = distance / value;
229                              objValue += movement * thisCost;
230                              rowValue = rowUpper[i];
231                              colsol[iCol] += movement;
232                         } else {
233                              // can't get there
234                              double movement = upper[iCol] - lower[iCol];
235                              objValue += movement * thisCost;
236                              rowValue -= movement * value;
237                              colsol[iCol] = upper[iCol];
238                              iCol = nextSlack[iCol];
239                         }
240                    }
241                    if (iCol >= 0) {
242                         // may want to carry on - because of cost?
243                         while (iCol >= 0 && cost[iCol] < 0 && rowValue > rowLower[i]) {
244                              // want to increase
245                              double distance = -(rowLower[i] - rowValue);
246                              double value = -element[columnStart[iCol]];
247                              double thisCost = cost[iCol];
248                              if (distance <= value*(upper[iCol] - colsol[iCol])) {
249                                   // can get there
250                                   double movement = distance / value;
251                                   objValue += movement * thisCost;
252                                   rowValue = rowLower[i];
253                                   colsol[iCol] += movement;
254                                   iCol = -1;
255                              } else {
256                                   // can't get there
257                                   double movement = upper[iCol] - colsol[iCol];
258                                   objValue += movement * thisCost;
259                                   rowValue -= movement * value;
260                                   colsol[iCol] = upper[iCol];
261                                   iCol = nextSlack[iCol];
262                              }
263                         }
264                         if (iCol >= 0 && colsol[iCol] > lower[iCol] + fixTolerance &&
265                                   colsol[iCol] < upper[iCol] - fixTolerance) {
266                              whenUsed_[i] = iteration;
267                              n++;
268                         }
269                    }
270                    rowsol[i] = rowValue;
271               }
272               infValue += CoinMax(CoinMax(0.0, rowLower[i] - rowsol[i]), rowsol[i] - rowUpper[i]);
273               // just change
274               rowsol[i] -= rowSave;
275          }
276          return n;
277     }
278}
279
280/* returns -1 if none or start of costed slacks or -2 if
281   there are costed slacks but it is messy */
282static int countCostedSlacks(OsiSolverInterface * model)
283{
284#ifdef OSI_IDIOT
285     const CoinPackedMatrix * matrix = model->getMatrixByCol();
286#else
287     ClpMatrixBase * matrix = model->clpMatrix();
288#endif
289     const int * row = matrix->getIndices();
290     const CoinBigIndex * columnStart = matrix->getVectorStarts();
291     const int * columnLength = matrix->getVectorLengths();
292     const double * element = matrix->getElements();
293     const double * rowupper = model->getRowUpper();
294     int nrows = model->getNumRows();
295     int ncols = model->getNumCols();
296     int slackStart = ncols - nrows;
297     int nSlacks = nrows;
298     int i;
299
300     if (ncols <= nrows) return -1;
301     while (1) {
302          for (i = 0; i < nrows; i++) {
303               int j = i + slackStart;
304               CoinBigIndex k = columnStart[j];
305               if (columnLength[j] == 1) {
306                    if (row[k] != i || element[k] != 1.0) {
307                         nSlacks = 0;
308                         break;
309                    }
310               } else {
311                    nSlacks = 0;
312                    break;
313               }
314               if (rowupper[i] <= 0.0) {
315                    nSlacks = 0;
316                    break;
317               }
318          }
319          if (nSlacks || !slackStart) break;
320          slackStart = 0;
321     }
322     if (!nSlacks) slackStart = -1;
323     return slackStart;
324}
325void
326Idiot::crash(int numberPass, CoinMessageHandler * handler,
327             const CoinMessages *messages, bool doCrossover)
328{
329     // lightweight options
330     int numberColumns = model_->getNumCols();
331     const double * objective = model_->getObjCoefficients();
332     int nnzero = 0;
333     double sum = 0.0;
334     int i;
335     for (i = 0; i < numberColumns; i++) {
336          if (objective[i]) {
337               sum += fabs(objective[i]);
338               nnzero++;
339          }
340     }
341     sum /= static_cast<double> (nnzero + 1);
342     if (maxIts_ == 5)
343          maxIts_ = 2;
344     if (numberPass <= 0)
345          majorIterations_ = static_cast<int>(2 + log10(static_cast<double>(numberColumns + 1)));
346     else
347          majorIterations_ = numberPass;
348     // If mu not changed then compute
349     if (mu_ == 1e-4)
350          mu_ = CoinMax(1.0e-3, sum * 1.0e-5);
351     if (maxIts2_ == 100) {
352          if (!lightWeight_) {
353               maxIts2_ = 105;
354          } else if (lightWeight_ == 1) {
355               mu_ *= 1000.0;
356               maxIts2_ = 23;
357          } else if (lightWeight_ == 2) {
358               maxIts2_ = 11;
359          } else {
360               maxIts2_ = 23;
361          }
362     }
363     //printf("setting mu to %g and doing %d passes\n",mu_,majorIterations_);
364     solve2(handler, messages);
365#ifndef OSI_IDIOT
366     if (doCrossover) {
367          double averageInfeas = model_->sumPrimalInfeasibilities() / static_cast<double> (model_->numberRows());
368          if ((averageInfeas < 0.01 && (strategy_ & 512) != 0) || (strategy_ & 8192) != 0)
369               crossOver(16 + 1);
370          else
371               crossOver(majorIterations_ < 1000000 ? 3 : 2);
372     }
373#endif
374}
375
376void
377Idiot::solve()
378{
379     CoinMessages dummy;
380     solve2(NULL, &dummy);
381}
382void
383Idiot::solve2(CoinMessageHandler * handler, const CoinMessages * messages)
384{
385     int strategy = 0;
386     double d2;
387     int i, n;
388     int allOnes = 1;
389     int iteration = 0;
390     int iterationTotal = 0;
391     int nTry = 0; /* number of tries at same weight */
392     double fixTolerance = IDIOT_FIX_TOLERANCE;
393     int maxBigIts = maxBigIts_;
394     int maxIts = maxIts_;
395     int logLevel = logLevel_;
396     int saveMajorIterations = majorIterations_;
397     majorIterations_ = majorIterations_ % 1000000;
398     if (handler) {
399          if (handler->logLevel() > 0 && handler->logLevel() < 3)
400               logLevel = 1;
401          else if (!handler->logLevel())
402               logLevel = 0;
403          else
404               logLevel = 7;
405     }
406     double djExit = djTolerance_;
407     double djFlag = 1.0 + 100.0 * djExit;
408     double djTol = 0.00001;
409     double mu = mu_;
410     double drop = drop_;
411     int maxIts2 = maxIts2_;
412     double factor = muFactor_;
413     double smallInfeas = smallInfeas_;
414     double reasonableInfeas = reasonableInfeas_;
415     double stopMu = stopMu_;
416     double maxmin, offset;
417     double lastWeighted = 1.0e50;
418     double exitDrop = exitDrop_;
419     double fakeSmall = smallInfeas;
420     double firstInfeas;
421     int badIts = 0;
422     int slackStart, slackEnd, ordStart, ordEnd;
423     int checkIteration = 0;
424     int lambdaIteration = 0;
425     int belowReasonable = 0; /* set if ever gone below reasonable infeas */
426     double bestWeighted = 1.0e60;
427     double bestFeasible = 1.0e60; /* best solution while feasible */
428     IdiotResult result, lastResult;
429     int saveStrategy = strategy_;
430     const int strategies[] = {0, 2, 128};
431     double saveLambdaScale = 0.0;
432     if ((saveStrategy & 128) != 0) {
433          fixTolerance = SMALL_IDIOT_FIX_TOLERANCE;
434     }
435#ifdef OSI_IDIOT
436     const CoinPackedMatrix * matrix = model_->getMatrixByCol();
437#else
438     ClpMatrixBase * matrix = model_->clpMatrix();
439#endif
440     const int * row = matrix->getIndices();
441     const CoinBigIndex * columnStart = matrix->getVectorStarts();
442     const int * columnLength = matrix->getVectorLengths();
443     const double * element = matrix->getElements();
444     int nrows = model_->getNumRows();
445     int ncols = model_->getNumCols();
446     double * rowsol, * colsol;
447     double * pi, * dj;
448#ifndef OSI_IDIOT
449     double * cost = model_->objective();
450     double * lower = model_->columnLower();
451     double * upper = model_->columnUpper();
452#else
453     double * cost = new double [ncols];
454     CoinMemcpyN( model_->getObjCoefficients(), ncols, cost);
455     const double * lower = model_->getColLower();
456     const double * upper = model_->getColUpper();
457#endif
458     const double *elemXX;
459     double * saveSol;
460     double * rowupper = new double[nrows]; // not const as modified
461     CoinMemcpyN(model_->getRowUpper(), nrows, rowupper);
462     double * rowlower = new double[nrows]; // not const as modified
463     CoinMemcpyN(model_->getRowLower(), nrows, rowlower);
464     CoinThreadRandom * randomNumberGenerator = model_->randomNumberGenerator();
465     int * whenUsed;
466     double * lambda;
467     saveSol = new double[ncols];
468     lambda = new double [nrows];
469     rowsol = new double[nrows];
470     colsol = new double [ncols];
471     CoinMemcpyN(model_->getColSolution(), ncols, colsol);
472     pi = new double[nrows];
473     dj = new double[ncols];
474     delete [] whenUsed_;
475     bool oddSlacks = false;
476     // See if any costed slacks
477     int numberSlacks = 0;
478     for (i = 0; i < ncols; i++) {
479          if (columnLength[i] == 1)
480               numberSlacks++;
481     }
482     if (!numberSlacks) {
483          whenUsed_ = new int[ncols];
484     } else {
485#ifdef COIN_DEVELOP
486          printf("%d slacks\n", numberSlacks);
487#endif
488          oddSlacks = true;
489          int extra = static_cast<int> (nrows * sizeof(double) / sizeof(int));
490          whenUsed_ = new int[2*ncols+2*nrows+extra];
491          int * posSlack = whenUsed_ + ncols;
492          int * negSlack = posSlack + nrows;
493          int * nextSlack = negSlack + nrows;
494          for (i = 0; i < nrows; i++) {
495               posSlack[i] = -1;
496               negSlack[i] = -1;
497          }
498          for (i = 0; i < ncols; i++)
499               nextSlack[i] = -1;
500          for (i = 0; i < ncols; i++) {
501               if (columnLength[i] == 1) {
502                    CoinBigIndex j = columnStart[i];
503                    int iRow = row[j];
504                    if (element[j] > 0.0) {
505                         if (posSlack[iRow] == -1) {
506                              posSlack[iRow] = i;
507                         } else {
508                              int iCol = posSlack[iRow];
509                              while (nextSlack[iCol] >= 0)
510                                   iCol = nextSlack[iCol];
511                              nextSlack[iCol] = i;
512                         }
513                    } else {
514                         if (negSlack[iRow] == -1) {
515                              negSlack[iRow] = i;
516                         } else {
517                              int iCol = negSlack[iRow];
518                              while (nextSlack[iCol] >= 0)
519                                   iCol = nextSlack[iCol];
520                              nextSlack[iCol] = i;
521                         }
522                    }
523               }
524          }
525          // now sort
526          for (i = 0; i < nrows; i++) {
527               int iCol;
528               iCol = posSlack[i];
529               if (iCol >= 0) {
530                    CoinBigIndex j = columnStart[iCol];
531#ifndef NDEBUG
532                    int iRow = row[j];
533#endif
534                    assert (element[j] > 0.0);
535                    assert (iRow == i);
536                    dj[0] = cost[iCol] / element[j];
537                    whenUsed_[0] = iCol;
538                    int n = 1;
539                    while (nextSlack[iCol] >= 0) {
540                         iCol = nextSlack[iCol];
541                         CoinBigIndex j = columnStart[iCol];
542#ifndef NDEBUG
543                         int iRow = row[j];
544#endif
545                         assert (element[j] > 0.0);
546                         assert (iRow == i);
547                         dj[n] = cost[iCol] / element[j];
548                         whenUsed_[n++] = iCol;
549                    }
550                    for (j = 0; j < n; j++) {
551                         int jCol = whenUsed_[j];
552                         nextSlack[jCol] = -2;
553                    }
554                    CoinSort_2(dj, dj + n, whenUsed_);
555                    // put back
556                    iCol = whenUsed_[0];
557                    posSlack[i] = iCol;
558                    for (j = 1; j < n; j++) {
559                         int jCol = whenUsed_[j];
560                         nextSlack[iCol] = jCol;
561                         iCol = jCol;
562                    }
563               }
564               iCol = negSlack[i];
565               if (iCol >= 0) {
566                    CoinBigIndex j = columnStart[iCol];
567#ifndef NDEBUG
568                    int iRow = row[j];
569#endif
570                    assert (element[j] < 0.0);
571                    assert (iRow == i);
572                    dj[0] = -cost[iCol] / element[j];
573                    whenUsed_[0] = iCol;
574                    int n = 1;
575                    while (nextSlack[iCol] >= 0) {
576                         iCol = nextSlack[iCol];
577                         CoinBigIndex j = columnStart[iCol];
578#ifndef NDEBUG
579                         int iRow = row[j];
580#endif
581                         assert (element[j] < 0.0);
582                         assert (iRow == i);
583                         dj[n] = -cost[iCol] / element[j];
584                         whenUsed_[n++] = iCol;
585                    }
586                    for (j = 0; j < n; j++) {
587                         int jCol = whenUsed_[j];
588                         nextSlack[jCol] = -2;
589                    }
590                    CoinSort_2(dj, dj + n, whenUsed_);
591                    // put back
592                    iCol = whenUsed_[0];
593                    negSlack[i] = iCol;
594                    for (j = 1; j < n; j++) {
595                         int jCol = whenUsed_[j];
596                         nextSlack[iCol] = jCol;
597                         iCol = jCol;
598                    }
599               }
600          }
601     }
602     whenUsed = whenUsed_;
603     if (model_->getObjSense() == -1.0) {
604          maxmin = -1.0;
605     } else {
606          maxmin = 1.0;
607     }
608     model_->getDblParam(OsiObjOffset, offset);
609     if (!maxIts2) maxIts2 = maxIts;
610     strategy = strategy_;
611     strategy &= 3;
612     memset(lambda, 0, nrows * sizeof(double));
613     slackStart = countCostedSlacks(model_);
614     if (slackStart >= 0) {
615          printf("This model has costed slacks\n");
616          slackEnd = slackStart + nrows;
617          if (slackStart) {
618               ordStart = 0;
619               ordEnd = slackStart;
620          } else {
621               ordStart = nrows;
622               ordEnd = ncols;
623          }
624     } else {
625          slackEnd = slackStart;
626          ordStart = 0;
627          ordEnd = ncols;
628     }
629     if (offset && logLevel > 2) {
630          printf("** Objective offset is %g\n", offset);
631     }
632     /* compute reasonable solution cost */
633     for (i = 0; i < nrows; i++) {
634          rowsol[i] = 1.0e31;
635     }
636     for (i = 0; i < ncols; i++) {
637          CoinBigIndex j;
638          for (j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
639               if (element[j] != 1.0) {
640                    allOnes = 0;
641                    break;
642               }
643          }
644     }
645     if (allOnes) {
646          elemXX = NULL;
647     } else {
648          elemXX = element;
649     }
650     // Do scaling if wanted
651     bool scaled = false;
652#ifndef OSI_IDIOT
653     if ((strategy_ & 32) != 0 && !allOnes) {
654          if (model_->scalingFlag() > 0)
655               scaled = model_->clpMatrix()->scale(model_) == 0;
656          if (scaled) {
657               const double * rowScale = model_->rowScale();
658               const double * columnScale = model_->columnScale();
659               double * oldLower = lower;
660               double * oldUpper = upper;
661               double * oldCost = cost;
662               lower = new double[ncols];
663               upper = new double[ncols];
664               cost = new double[ncols];
665               CoinMemcpyN(oldLower, ncols, lower);
666               CoinMemcpyN(oldUpper, ncols, upper);
667               CoinMemcpyN(oldCost, ncols, cost);
668               int icol, irow;
669               for (icol = 0; icol < ncols; icol++) {
670                    double multiplier = 1.0 / columnScale[icol];
671                    if (lower[icol] > -1.0e50)
672                         lower[icol] *= multiplier;
673                    if (upper[icol] < 1.0e50)
674                         upper[icol] *= multiplier;
675                    colsol[icol] *= multiplier;
676                    cost[icol] *= columnScale[icol];
677               }
678               CoinMemcpyN(model_->rowLower(), nrows, rowlower);
679               for (irow = 0; irow < nrows; irow++) {
680                    double multiplier = rowScale[irow];
681                    if (rowlower[irow] > -1.0e50)
682                         rowlower[irow] *= multiplier;
683                    if (rowupper[irow] < 1.0e50)
684                         rowupper[irow] *= multiplier;
685                    rowsol[irow] *= multiplier;
686               }
687               int length = columnStart[ncols-1] + columnLength[ncols-1];
688               double * elemYY = new double[length];
689               for (i = 0; i < ncols; i++) {
690                    CoinBigIndex j;
691                    double scale = columnScale[i];
692                    for (j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
693                         int irow = row[j];
694                         elemYY[j] = element[j] * scale * rowScale[irow];
695                    }
696               }
697               elemXX = elemYY;
698          }
699     }
700#endif
701     for (i = 0; i < ncols; i++) {
702          CoinBigIndex j;
703          double dd = columnLength[i];
704          dd = cost[i] / dd;
705          for (j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
706               int irow = row[j];
707               if (dd < rowsol[irow]) {
708                    rowsol[irow] = dd;
709               }
710          }
711     }
712     d2 = 0.0;
713     for (i = 0; i < nrows; i++) {
714          d2 += rowsol[i];
715     }
716     d2 *= 2.0; /* for luck */
717
718     d2 = d2 / static_cast<double> (4 * nrows + 8000);
719     d2 *= 0.5; /* halve with more flexible method */
720     if (d2 < 5.0) d2 = 5.0;
721     if (djExit == 0.0) {
722          djExit = d2;
723     }
724     if ((saveStrategy & 4) != 0) {
725          /* go to relative tolerances - first small */
726          djExit = 1.0e-10;
727          djFlag = 1.0e-5;
728          drop = 1.0e-10;
729     }
730     memset(whenUsed, 0, ncols * sizeof(int));
731     strategy = strategies[strategy];
732     if ((saveStrategy & 8) != 0) strategy |= 64; /* don't allow large theta */
733     CoinMemcpyN(colsol, ncols, saveSol);
734
735     lastResult = IdiSolve(nrows, ncols, rowsol , colsol, pi,
736                           dj, cost, rowlower, rowupper,
737                           lower, upper, elemXX, row, columnStart, columnLength, lambda,
738                           0, mu, drop,
739                           maxmin, offset, strategy, djTol, djExit, djFlag, randomNumberGenerator);
740     // update whenUsed_
741     n = cleanIteration(iteration, ordStart, ordEnd,
742                        colsol,  lower,  upper,
743                        rowlower, rowupper,
744                        cost, elemXX, fixTolerance, lastResult.objval, lastResult.infeas);
745     if ((strategy_ & 16384) != 0) {
746          int * posSlack = whenUsed_ + ncols;
747          int * negSlack = posSlack + nrows;
748          int * nextSlack = negSlack + nrows;
749          double * rowsol2 = reinterpret_cast<double *> (nextSlack + ncols);
750          for (i = 0; i < nrows; i++)
751               rowsol[i] += rowsol2[i];
752     }
753     if ((logLevel_ & 1) != 0) {
754#ifndef OSI_IDIOT
755          if (!handler) {
756#endif
757               printf("Iteration %d infeasibility %g, objective %g - mu %g, its %d, %d interior\n",
758                      iteration, lastResult.infeas, lastResult.objval, mu, lastResult.iteration, n);
759#ifndef OSI_IDIOT
760          } else {
761               handler->message(CLP_IDIOT_ITERATION, *messages)
762                         << iteration << lastResult.infeas << lastResult.objval << mu << lastResult.iteration << n
763                         << CoinMessageEol;
764          }
765#endif
766     }
767     int numberBaseTrys = 0; // for first time
768     int numberAway = -1;
769     iterationTotal = lastResult.iteration;
770     firstInfeas = lastResult.infeas;
771     if ((strategy_ & 1024) != 0) reasonableInfeas = 0.5 * firstInfeas;
772     if (lastResult.infeas < reasonableInfeas) lastResult.infeas = reasonableInfeas;
773     double keepinfeas = 1.0e31;
774     double lastInfeas = 1.0e31;
775     double bestInfeas = 1.0e31;
776     while ((mu > stopMu && lastResult.infeas > smallInfeas) ||
777               (lastResult.infeas <= smallInfeas &&
778                dropping(lastResult, exitDrop, smallInfeas, &badIts)) ||
779               checkIteration < 2 || lambdaIteration < lambdaIterations_) {
780          if (lastResult.infeas <= exitFeasibility_)
781               break;
782          iteration++;
783          checkIteration++;
784          if (lastResult.infeas <= smallInfeas && lastResult.objval < bestFeasible) {
785               bestFeasible = lastResult.objval;
786          }
787          if (lastResult.infeas + mu * lastResult.objval < bestWeighted) {
788               bestWeighted = lastResult.objval + mu * lastResult.objval;
789          }
790          if ((saveStrategy & 4096)) strategy |= 256;
791          if ((saveStrategy & 4) != 0 && iteration > 2) {
792               /* go to relative tolerances */
793               double weighted = 10.0 + fabs(lastWeighted);
794               djExit = djTolerance_ * weighted;
795               djFlag = 2.0 * djExit;
796               drop = drop_ * weighted;
797               djTol = 0.01 * djExit;
798          }
799          CoinMemcpyN(colsol, ncols, saveSol);
800          result = IdiSolve(nrows, ncols, rowsol , colsol, pi, dj,
801                            cost, rowlower, rowupper,
802                            lower, upper, elemXX, row, columnStart, columnLength, lambda,
803                            maxIts, mu, drop,
804                            maxmin, offset, strategy, djTol, djExit, djFlag, randomNumberGenerator);
805          n = cleanIteration(iteration, ordStart, ordEnd,
806                             colsol,  lower,  upper,
807                             rowlower, rowupper,
808                             cost, elemXX, fixTolerance, result.objval, result.infeas);
809          if ((strategy_ & 16384) != 0) {
810               int * posSlack = whenUsed_ + ncols;
811               int * negSlack = posSlack + nrows;
812               int * nextSlack = negSlack + nrows;
813               double * rowsol2 = reinterpret_cast<double *> (nextSlack + ncols);
814               for (i = 0; i < nrows; i++)
815                    rowsol[i] += rowsol2[i];
816          }
817          if ((logLevel_ & 1) != 0) {
818#ifndef OSI_IDIOT
819               if (!handler) {
820#endif
821                    printf("Iteration %d infeasibility %g, objective %g - mu %g, its %d, %d interior\n",
822                           iteration, result.infeas, result.objval, mu, result.iteration, n);
823#ifndef OSI_IDIOT
824               } else {
825                    handler->message(CLP_IDIOT_ITERATION, *messages)
826                              << iteration << result.infeas << result.objval << mu << result.iteration << n
827                              << CoinMessageEol;
828               }
829#endif
830          }
831          if (iteration > 50 && n == numberAway && result.infeas < 1.0e-4) {
832#ifdef CLP_INVESTIGATE
833               printf("infeas small %g\n", result.infeas);
834#endif
835               break; // not much happening
836          }
837          if (lightWeight_ == 1 && iteration > 10 && result.infeas > 1.0 && maxIts != 7) {
838               if (lastInfeas != bestInfeas && CoinMin(result.infeas, lastInfeas) > 0.95 * bestInfeas)
839                    majorIterations_ = CoinMin(majorIterations_, iteration); // not getting feasible
840          }
841          lastInfeas = result.infeas;
842          numberAway = n;
843          keepinfeas = result.infeas;
844          lastWeighted = result.weighted;
845          iterationTotal += result.iteration;
846          if (iteration == 1) {
847               if ((strategy_ & 1024) != 0 && mu < 1.0e-10)
848                    result.infeas = firstInfeas * 0.8;
849               if (majorIterations_ >= 50 || dropEnoughFeasibility_ <= 0.0)
850                    result.infeas *= 0.8;
851               if (result.infeas > firstInfeas * 0.9
852                         && result.infeas > reasonableInfeas) {
853                    iteration--;
854                    if (majorIterations_ < 50)
855                         mu *= 1.0e-1;
856                    else
857                         mu *= 0.7;
858                    bestFeasible = 1.0e31;
859                    bestWeighted = 1.0e60;
860                    numberBaseTrys++;
861                    if (mu < 1.0e-30 || (numberBaseTrys > 10 && lightWeight_)) {
862                         // back to all slack basis
863                         lightWeight_ = 2;
864                         break;
865                    }
866                    CoinMemcpyN(saveSol, ncols, colsol);
867               } else {
868                    maxIts = maxIts2;
869                    checkIteration = 0;
870                    if ((strategy_ & 1024) != 0) mu *= 1.0e-1;
871               }
872          } else {
873          }
874          bestInfeas = CoinMin(bestInfeas, result.infeas);
875          if (iteration) {
876               /* this code is in to force it to terminate sometime */
877               double changeMu = factor;
878               if ((saveStrategy & 64) != 0) {
879                    keepinfeas = 0.0; /* switch off ranga's increase */
880                    fakeSmall = smallInfeas;
881               } else {
882                    fakeSmall = -1.0;
883               }
884               saveLambdaScale = 0.0;
885               if (result.infeas > reasonableInfeas ||
886                         (nTry + 1 == maxBigIts && result.infeas > fakeSmall)) {
887                    if (result.infeas > lastResult.infeas*(1.0 - dropEnoughFeasibility_) ||
888                              nTry + 1 == maxBigIts ||
889                              (result.infeas > lastResult.infeas * 0.9
890                               && result.weighted > lastResult.weighted
891                               - dropEnoughWeighted_ * CoinMax(fabs(lastResult.weighted), fabs(result.weighted)))) {
892                         mu *= changeMu;
893                         if ((saveStrategy & 32) != 0 && result.infeas < reasonableInfeas && 0) {
894                              reasonableInfeas = CoinMax(smallInfeas, reasonableInfeas * sqrt(changeMu));
895                              printf("reasonable infeas now %g\n", reasonableInfeas);
896                         }
897                         result.weighted = 1.0e60;
898                         nTry = 0;
899                         bestFeasible = 1.0e31;
900                         bestWeighted = 1.0e60;
901                         checkIteration = 0;
902                         lambdaIteration = 0;
903#define LAMBDA
904#ifdef LAMBDA
905                         if ((saveStrategy & 2048) == 0) {
906                              memset(lambda, 0, nrows * sizeof(double));
907                         }
908#else
909                         memset(lambda, 0, nrows * sizeof(double));
910#endif
911                    } else {
912                         nTry++;
913                    }
914               } else if (lambdaIterations_ >= 0) {
915                    /* update lambda  */
916                    double scale = 1.0 / mu;
917                    int i, nnz = 0;
918                    saveLambdaScale = scale;
919                    lambdaIteration++;
920                    if ((saveStrategy & 4) == 0) drop = drop_ / 50.0;
921                    if (lambdaIteration > 4 &&
922                              (((lambdaIteration % 10) == 0 && smallInfeas < keepinfeas) ||
923                               ((lambdaIteration % 5) == 0 && 1.5 * smallInfeas < keepinfeas))) {
924                         //printf(" Increasing smallInfeas from %f to %f\n",smallInfeas,1.5*smallInfeas);
925                         smallInfeas *= 1.5;
926                    }
927                    if ((saveStrategy & 2048) == 0) {
928                         for (i = 0; i < nrows; i++) {
929                              if (lambda[i]) nnz++;
930                              lambda[i] += scale * rowsol[i];
931                         }
932                    } else {
933                         nnz = 1;
934#ifdef LAMBDA
935                         for (i = 0; i < nrows; i++) {
936                              lambda[i] += scale * rowsol[i];
937                         }
938#else
939                         for (i = 0; i < nrows; i++) {
940                              lambda[i] = scale * rowsol[i];
941                         }
942                         for (i = 0; i < ncols; i++) {
943                              CoinBigIndex j;
944                              double value = cost[i] * maxmin;
945                              for (j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
946                                   int irow = row[j];
947                                   value += element[j] * lambda[irow];
948                              }
949                              cost[i] = value * maxmin;
950                         }
951                         for (i = 0; i < nrows; i++) {
952                              offset += lambda[i] * rowupper[i];
953                              lambda[i] = 0.0;
954                         }
955#ifdef DEBUG
956                         printf("offset %g\n", offset);
957#endif
958                         model_->setDblParam(OsiObjOffset, offset);
959#endif
960                    }
961                    nTry++;
962                    if (!nnz) {
963                         bestFeasible = 1.0e32;
964                         bestWeighted = 1.0e60;
965                         checkIteration = 0;
966                         result.weighted = 1.0e31;
967                    }
968#ifdef DEBUG
969                    double trueCost = 0.0;
970                    for (i = 0; i < ncols; i++) {
971                         int j;
972                         trueCost += cost[i] * colsol[i];
973                    }
974                    printf("True objective %g\n", trueCost - offset);
975#endif
976               } else {
977                    nTry++;
978               }
979               lastResult = result;
980               if (result.infeas < reasonableInfeas && !belowReasonable) {
981                    belowReasonable = 1;
982                    bestFeasible = 1.0e32;
983                    bestWeighted = 1.0e60;
984                    checkIteration = 0;
985                    result.weighted = 1.0e31;
986               }
987          }
988          if (iteration >= majorIterations_) {
989               // If not feasible and crash then dive dive dive
990               if (mu > 1.0e-12 && result.infeas > 1.0 && majorIterations_ < 40) {
991                    mu = 1.0e-30;
992                    majorIterations_ = iteration + 1;
993                    stopMu = 0.0;
994               } else {
995                    if (logLevel > 2)
996                         printf("Exiting due to number of major iterations\n");
997                    break;
998               }
999          }
1000     }
1001     majorIterations_ = saveMajorIterations;
1002#ifndef OSI_IDIOT
1003     if (scaled) {
1004          // Scale solution and free arrays
1005          const double * rowScale = model_->rowScale();
1006          const double * columnScale = model_->columnScale();
1007          int icol, irow;
1008          for (icol = 0; icol < ncols; icol++) {
1009               colsol[icol] *= columnScale[icol];
1010               saveSol[icol] *= columnScale[icol];
1011               dj[icol] /= columnScale[icol];
1012          }
1013          for (irow = 0; irow < nrows; irow++) {
1014               rowsol[irow] /= rowScale[irow];
1015               pi[irow] *= rowScale[irow];
1016          }
1017          // Don't know why getting Microsoft problems
1018#if defined (_MSC_VER)
1019          delete [] ( double *) elemXX;
1020#else
1021          delete [] elemXX;
1022#endif
1023          model_->setRowScale(NULL);
1024          model_->setColumnScale(NULL);
1025          delete [] lower;
1026          delete [] upper;
1027          delete [] cost;
1028          lower = model_->columnLower();
1029          upper = model_->columnUpper();
1030          cost = model_->objective();
1031          //rowlower = model_->rowLower();
1032     }
1033#endif
1034#define TRYTHIS
1035#ifdef TRYTHIS
1036     if ((saveStrategy & 2048) != 0) {
1037          double offset;
1038          model_->getDblParam(OsiObjOffset, offset);
1039          for (i = 0; i < ncols; i++) {
1040               CoinBigIndex j;
1041               double djval = cost[i] * maxmin;
1042               for (j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
1043                    int irow = row[j];
1044                    djval -= element[j] * lambda[irow];
1045               }
1046               cost[i] = djval;
1047          }
1048          for (i = 0; i < nrows; i++) {
1049               offset += lambda[i] * rowupper[i];
1050          }
1051          model_->setDblParam(OsiObjOffset, offset);
1052     }
1053#endif
1054     if (saveLambdaScale) {
1055          /* back off last update */
1056          for (i = 0; i < nrows; i++) {
1057               lambda[i] -= saveLambdaScale * rowsol[i];
1058          }
1059     }
1060     muAtExit_ = mu;
1061     // For last iteration make as feasible as possible
1062     if (oddSlacks)
1063          strategy_ |= 16384;
1064     // not scaled
1065     n = cleanIteration(iteration, ordStart, ordEnd,
1066                        colsol,  lower,  upper,
1067                        model_->rowLower(), model_->rowUpper(),
1068                        cost, element, fixTolerance, lastResult.objval, lastResult.infeas);
1069#if 0
1070     if ((logLevel & 1) == 0 || (strategy_ & 16384) != 0) {
1071          printf(
1072               "%d - mu %g, infeasibility %g, objective %g, %d interior\n",
1073               iteration, mu, lastResult.infeas, lastResult.objval, n);
1074     }
1075#endif
1076#ifndef OSI_IDIOT
1077     model_->setSumPrimalInfeasibilities(lastResult.infeas);
1078#endif
1079     // Put back more feasible solution
1080     double saveInfeas[] = {0.0, 0.0};
1081     for (int iSol = 0; iSol < 3; iSol++) {
1082          const double * solution = iSol ? colsol : saveSol;
1083          if (iSol == 2 && saveInfeas[0] < saveInfeas[1]) {
1084               // put back best solution
1085               CoinMemcpyN(saveSol, ncols, colsol);
1086          }
1087          double large = 0.0;
1088          int i;
1089          memset(rowsol, 0, nrows * sizeof(double));
1090          for (i = 0; i < ncols; i++) {
1091               CoinBigIndex j;
1092               double value = solution[i];
1093               for (j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
1094                    int irow = row[j];
1095                    rowsol[irow] += element[j] * value;
1096               }
1097          }
1098          for (i = 0; i < nrows; i++) {
1099               if (rowsol[i] > rowupper[i]) {
1100                    double diff = rowsol[i] - rowupper[i];
1101                    if (diff > large)
1102                         large = diff;
1103               } else if (rowsol[i] < rowlower[i]) {
1104                    double diff = rowlower[i] - rowsol[i];
1105                    if (diff > large)
1106                         large = diff;
1107               }
1108          }
1109          if (iSol < 2)
1110               saveInfeas[iSol] = large;
1111          if (logLevel > 2)
1112               printf("largest infeasibility is %g\n", large);
1113     }
1114     /* subtract out lambda */
1115     for (i = 0; i < nrows; i++) {
1116          pi[i] -= lambda[i];
1117     }
1118     for (i = 0; i < ncols; i++) {
1119          CoinBigIndex j;
1120          double djval = cost[i] * maxmin;
1121          for (j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
1122               int irow = row[j];
1123               djval -= element[j] * pi[irow];
1124          }
1125          dj[i] = djval;
1126     }
1127     if ((strategy_ & 1024) != 0) {
1128          double ratio = static_cast<double> (ncols) / static_cast<double> (nrows);
1129          printf("col/row ratio %g infeas ratio %g\n", ratio, lastResult.infeas / firstInfeas);
1130          if (lastResult.infeas > 0.01 * firstInfeas * ratio) {
1131               strategy_ &= (~1024);
1132               printf(" - layer off\n");
1133          } else {
1134               printf(" - layer on\n");
1135          }
1136     }
1137     delete [] saveSol;
1138     delete [] lambda;
1139     // save solution
1140     // duals not much use - but save anyway
1141#ifndef OSI_IDIOT
1142     CoinMemcpyN(rowsol, nrows, model_->primalRowSolution());
1143     CoinMemcpyN(colsol, ncols, model_->primalColumnSolution());
1144     CoinMemcpyN(pi, nrows, model_->dualRowSolution());
1145     CoinMemcpyN(dj, ncols, model_->dualColumnSolution());
1146#else
1147     model_->setColSolution(colsol);
1148     model_->setRowPrice(pi);
1149     delete [] cost;
1150#endif
1151     delete [] rowsol;
1152     delete [] colsol;
1153     delete [] pi;
1154     delete [] dj;
1155     delete [] rowlower;
1156     delete [] rowupper;
1157     return ;
1158}
1159#ifndef OSI_IDIOT
1160void
1161Idiot::crossOver(int mode)
1162{
1163     if (lightWeight_ == 2) {
1164          // total failure
1165          model_->allSlackBasis();
1166          return;
1167     }
1168     double fixTolerance = IDIOT_FIX_TOLERANCE;
1169#ifdef COIN_DEVELOP
1170     double startTime = CoinCpuTime();
1171#endif
1172     ClpSimplex * saveModel = NULL;
1173     ClpMatrixBase * matrix = model_->clpMatrix();
1174     const int * row = matrix->getIndices();
1175     const CoinBigIndex * columnStart = matrix->getVectorStarts();
1176     const int * columnLength = matrix->getVectorLengths();
1177     const double * element = matrix->getElements();
1178     const double * rowupper = model_->getRowUpper();
1179     int nrows = model_->getNumRows();
1180     int ncols = model_->getNumCols();
1181     double * rowsol, * colsol;
1182     // different for Osi
1183     double * lower = model_->columnLower();
1184     double * upper = model_->columnUpper();
1185     const double * rowlower = model_->getRowLower();
1186     int * whenUsed = whenUsed_;
1187     rowsol = model_->primalRowSolution();
1188     colsol = model_->primalColumnSolution();;
1189     double * cost = model_->objective();
1190
1191     int slackEnd, ordStart, ordEnd;
1192     int slackStart = countCostedSlacks(model_);
1193
1194     int addAll = mode & 7;
1195     int presolve = 0;
1196
1197     double djTolerance = djTolerance_;
1198     if (djTolerance > 0.0 && djTolerance < 1.0)
1199          djTolerance = 1.0;
1200     int iteration;
1201     int i, n = 0;
1202     double ratio = 1.0;
1203     double objValue = 0.0;
1204     if ((strategy_ & 128) != 0) {
1205          fixTolerance = SMALL_IDIOT_FIX_TOLERANCE;
1206     }
1207     if ((mode & 16) != 0 && addAll < 3) presolve = 1;
1208     double * saveUpper = NULL;
1209     double * saveLower = NULL;
1210     double * saveRowUpper = NULL;
1211     double * saveRowLower = NULL;
1212     bool allowInfeasible = ((strategy_ & 8192) != 0) || (majorIterations_ > 1000000);
1213     if (addAll < 3) {
1214          saveUpper = new double [ncols];
1215          saveLower = new double [ncols];
1216          CoinMemcpyN(upper, ncols, saveUpper);
1217          CoinMemcpyN(lower, ncols, saveLower);
1218          if (allowInfeasible) {
1219               saveRowUpper = new double [nrows];
1220               saveRowLower = new double [nrows];
1221               CoinMemcpyN(rowupper, nrows, saveRowUpper);
1222               CoinMemcpyN(rowlower, nrows, saveRowLower);
1223               double averageInfeas = model_->sumPrimalInfeasibilities() / static_cast<double> (model_->numberRows());
1224               fixTolerance = CoinMax(fixTolerance, 1.0e-5 * averageInfeas);
1225          }
1226     }
1227     if (slackStart >= 0) {
1228          slackEnd = slackStart + nrows;
1229          if (slackStart) {
1230               ordStart = 0;
1231               ordEnd = slackStart;
1232          } else {
1233               ordStart = nrows;
1234               ordEnd = ncols;
1235          }
1236     } else {
1237          slackEnd = slackStart;
1238          ordStart = 0;
1239          ordEnd = ncols;
1240     }
1241     /* get correct rowsol (without known slacks) */
1242     memset(rowsol, 0, nrows * sizeof(double));
1243     for (i = ordStart; i < ordEnd; i++) {
1244          CoinBigIndex j;
1245          double value = colsol[i];
1246          if (value < lower[i] + fixTolerance) {
1247               value = lower[i];
1248               colsol[i] = value;
1249          }
1250          for (j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
1251               int irow = row[j];
1252               rowsol[irow] += value * element[j];
1253          }
1254     }
1255     if (slackStart >= 0) {
1256          for (i = 0; i < nrows; i++) {
1257               if (ratio * rowsol[i] > rowlower[i] && rowsol[i] > 1.0e-8) {
1258                    ratio = rowlower[i] / rowsol[i];
1259               }
1260          }
1261          for (i = 0; i < nrows; i++) {
1262               rowsol[i] *= ratio;
1263          }
1264          for (i = ordStart; i < ordEnd; i++) {
1265               double value = colsol[i] * ratio;
1266               colsol[i] = value;
1267               objValue += value * cost[i];
1268          }
1269          for (i = 0; i < nrows; i++) {
1270               double value = rowlower[i] - rowsol[i];
1271               colsol[i+slackStart] = value;
1272               objValue += value * cost[i+slackStart];
1273          }
1274          printf("New objective after scaling %g\n", objValue);
1275     }
1276#if 0
1277     maybe put back - but just get feasible ?
1278     // If not many fixed then just exit
1279     int numberFixed = 0;
1280     for (i = ordStart; i < ordEnd; i++) {
1281          if (colsol[i] < lower[i] + fixTolerance)
1282               numberFixed++;
1283          else if (colsol[i] > upper[i] - fixTolerance)
1284               numberFixed++;
1285     }
1286     if (numberFixed < ncols / 2) {
1287          addAll = 3;
1288          presolve = 0;
1289     }
1290#endif
1291#ifdef FEB_TRY
1292     int savePerturbation = model_->perturbation();
1293     int saveOptions = model_->specialOptions();
1294     model_->setSpecialOptions(saveOptions | 8192);
1295     if (savePerturbation_ == 50)
1296          model_->setPerturbation(56);
1297#endif
1298     model_->createStatus();
1299     /* addAll
1300        0 - chosen,all used, all
1301        1 - chosen, all
1302        2 - all
1303        3 - do not do anything  - maybe basis
1304     */
1305     for (i = ordStart; i < ordEnd; i++) {
1306          if (addAll < 2) {
1307               if (colsol[i] < lower[i] + fixTolerance) {
1308                    upper[i] = lower[i];
1309                    colsol[i] = lower[i];
1310               } else if (colsol[i] > upper[i] - fixTolerance) {
1311                    lower[i] = upper[i];
1312                    colsol[i] = upper[i];
1313               }
1314          }
1315          model_->setColumnStatus(i, ClpSimplex::superBasic);
1316     }
1317     if ((strategy_ & 16384) != 0) {
1318          // put in basis
1319          int * posSlack = whenUsed_ + ncols;
1320          int * negSlack = posSlack + nrows;
1321          int * nextSlack = negSlack + nrows;
1322          /* Laci - try both ways - to see what works -
1323             you can change second part as much as you want */
1324#ifndef LACI_TRY  // was #if 1
1325          // Array for sorting out slack values
1326          double * ratio = new double [ncols];
1327          int * which = new int [ncols];
1328          for (i = 0; i < nrows; i++) {
1329               if (posSlack[i] >= 0 || negSlack[i] >= 0) {
1330                    int iCol;
1331                    int nPlus = 0;
1332                    int nMinus = 0;
1333                    bool possible = true;
1334                    // Get sum
1335                    double sum = 0.0;
1336                    iCol = posSlack[i];
1337                    while (iCol >= 0) {
1338                         double value = element[columnStart[iCol]];
1339                         sum += value * colsol[iCol];
1340                         if (lower[iCol]) {
1341                              possible = false;
1342                              break;
1343                         } else {
1344                              nPlus++;
1345                         }
1346                         iCol = nextSlack[iCol];
1347                    }
1348                    iCol = negSlack[i];
1349                    while (iCol >= 0) {
1350                         double value = -element[columnStart[iCol]];
1351                         sum -= value * colsol[iCol];
1352                         if (lower[iCol]) {
1353                              possible = false;
1354                              break;
1355                         } else {
1356                              nMinus++;
1357                         }
1358                         iCol = nextSlack[iCol];
1359                    }
1360                    //printf("%d plus, %d minus",nPlus,nMinus);
1361                    //printf("\n");
1362                    if ((rowsol[i] - rowlower[i] < 1.0e-7 ||
1363                              rowupper[i] - rowsol[i] < 1.0e-7) &&
1364                              nPlus + nMinus < 2)
1365                         possible = false;
1366                    if (possible) {
1367                         // Amount contributed by other varaibles
1368                         sum = rowsol[i] - sum;
1369                         double lo = rowlower[i];
1370                         if (lo > -1.0e20)
1371                              lo -= sum;
1372                         double up = rowupper[i];
1373                         if (up < 1.0e20)
1374                              up -= sum;
1375                         //printf("row bounds %g %g\n",lo,up);
1376                         if (0) {
1377                              double sum = 0.0;
1378                              double x = 0.0;
1379                              for (int k = 0; k < ncols; k++) {
1380                                   CoinBigIndex j;
1381                                   double value = colsol[k];
1382                                   x += value * cost[k];
1383                                   for (j = columnStart[k]; j < columnStart[k] + columnLength[k]; j++) {
1384                                        int irow = row[j];
1385                                        if (irow == i)
1386                                             sum += element[j] * value;
1387                                   }
1388                              }
1389                              printf("Before sum %g <= %g <= %g cost %.18g\n",
1390                                     rowlower[i], sum, rowupper[i], x);
1391                         }
1392                         // set all to zero
1393                         iCol = posSlack[i];
1394                         while (iCol >= 0) {
1395                              colsol[iCol] = 0.0;
1396                              iCol = nextSlack[iCol];
1397                         }
1398                         iCol = negSlack[i];
1399                         while (iCol >= 0) {
1400                              colsol[iCol] = 0.0;
1401                              iCol = nextSlack[iCol];
1402                         }
1403                         {
1404                              int iCol;
1405                              iCol = posSlack[i];
1406                              while (iCol >= 0) {
1407                                   //printf("col %d el %g sol %g bounds %g %g cost %g\n",
1408                                   //     iCol,element[columnStart[iCol]],
1409                                   //     colsol[iCol],lower[iCol],upper[iCol],cost[iCol]);
1410                                   iCol = nextSlack[iCol];
1411                              }
1412                              iCol = negSlack[i];
1413                              while (iCol >= 0) {
1414                                   //printf("col %d el %g sol %g bounds %g %g cost %g\n",
1415                                   //     iCol,element[columnStart[iCol]],
1416                                   //     colsol[iCol],lower[iCol],upper[iCol],cost[iCol]);
1417                                   iCol = nextSlack[iCol];
1418                              }
1419                         }
1420                         //printf("now what?\n");
1421                         int n = 0;
1422                         bool basic = false;
1423                         if (lo > 0.0) {
1424                              // Add in positive
1425                              iCol = posSlack[i];
1426                              while (iCol >= 0) {
1427                                   double value = element[columnStart[iCol]];
1428                                   ratio[n] = cost[iCol] / value;
1429                                   which[n++] = iCol;
1430                                   iCol = nextSlack[iCol];
1431                              }
1432                              CoinSort_2(ratio, ratio + n, which);
1433                              for (int i = 0; i < n; i++) {
1434                                   iCol = which[i];
1435                                   double value = element[columnStart[iCol]];
1436                                   if (lo >= upper[iCol]*value) {
1437                                        value *= upper[iCol];
1438                                        sum += value;
1439                                        lo -= value;
1440                                        colsol[iCol] = upper[iCol];
1441                                   } else {
1442                                        value = lo / value;
1443                                        sum += lo;
1444                                        lo = 0.0;
1445                                        colsol[iCol] = value;
1446                                        model_->setColumnStatus(iCol, ClpSimplex::basic);
1447                                        basic = true;
1448                                   }
1449                                   if (lo < 1.0e-7)
1450                                        break;
1451                              }
1452                         } else if (up < 0.0) {
1453                              // Use lo so coding is more similar
1454                              lo = -up;
1455                              // Add in negative
1456                              iCol = negSlack[i];
1457                              while (iCol >= 0) {
1458                                   double value = -element[columnStart[iCol]];
1459                                   ratio[n] = cost[iCol] / value;
1460                                   which[n++] = iCol;
1461                                   iCol = nextSlack[iCol];
1462                              }
1463                              CoinSort_2(ratio, ratio + n, which);
1464                              for (int i = 0; i < n; i++) {
1465                                   iCol = which[i];
1466                                   double value = -element[columnStart[iCol]];
1467                                   if (lo >= upper[iCol]*value) {
1468                                        value *= upper[iCol];
1469                                        sum += value;
1470                                        lo -= value;
1471                                        colsol[iCol] = upper[iCol];
1472                                   } else {
1473                                        value = lo / value;
1474                                        sum += lo;
1475                                        lo = 0.0;
1476                                        colsol[iCol] = value;
1477                                        model_->setColumnStatus(iCol, ClpSimplex::basic);
1478                                        basic = true;
1479                                   }
1480                                   if (lo < 1.0e-7)
1481                                        break;
1482                              }
1483                         }
1484                         if (0) {
1485                              double sum2 = 0.0;
1486                              double x = 0.0;
1487                              for (int k = 0; k < ncols; k++) {
1488                                   CoinBigIndex j;
1489                                   double value = colsol[k];
1490                                   x += value * cost[k];
1491                                   for (j = columnStart[k]; j < columnStart[k] + columnLength[k]; j++) {
1492                                        int irow = row[j];
1493                                        if (irow == i)
1494                                             sum2 += element[j] * value;
1495                                   }
1496                              }
1497                              printf("after sum %g <= %g <= %g cost %.18g (sum = %g)\n",
1498                                     rowlower[i], sum2, rowupper[i], x, sum);
1499                         }
1500                         rowsol[i] = sum;
1501                         if (basic) {
1502                              if (fabs(rowsol[i] - rowlower[i]) < fabs(rowsol[i] - rowupper[i]))
1503                                   model_->setRowStatus(i, ClpSimplex::atLowerBound);
1504                              else
1505                                   model_->setRowStatus(i, ClpSimplex::atUpperBound);
1506                         }
1507                    } else {
1508                         int n = 0;
1509                         int iCol;
1510                         iCol = posSlack[i];
1511                         while (iCol >= 0) {
1512                              if (colsol[iCol] > lower[iCol] + 1.0e-8 &&
1513                                        colsol[iCol] < upper[iCol] - 1.0e-8) {
1514                                   model_->setColumnStatus(iCol, ClpSimplex::basic);
1515                                   n++;
1516                              }
1517                              iCol = nextSlack[iCol];
1518                         }
1519                         iCol = negSlack[i];
1520                         while (iCol >= 0) {
1521                              if (colsol[iCol] > lower[iCol] + 1.0e-8 &&
1522                                        colsol[iCol] < upper[iCol] - 1.0e-8) {
1523                                   model_->setColumnStatus(iCol, ClpSimplex::basic);
1524                                   n++;
1525                              }
1526                              iCol = nextSlack[iCol];
1527                         }
1528                         if (n) {
1529                              if (fabs(rowsol[i] - rowlower[i]) < fabs(rowsol[i] - rowupper[i]))
1530                                   model_->setRowStatus(i, ClpSimplex::atLowerBound);
1531                              else
1532                                   model_->setRowStatus(i, ClpSimplex::atUpperBound);
1533#ifdef CLP_INVESTIGATE
1534                              if (n > 1)
1535                                   printf("%d basic on row %d!\n", n, i);
1536#endif
1537                         }
1538                    }
1539               }
1540          }
1541          delete [] ratio;
1542          delete [] which;
1543#else
1544          for (i = 0; i < nrows; i++) {
1545               int n = 0;
1546               int iCol;
1547               iCol = posSlack[i];
1548               while (iCol >= 0) {
1549                    if (colsol[iCol] > lower[iCol] + 1.0e-8 &&
1550                              colsol[iCol] < upper[iCol] - 1.0e-8) {
1551                         model_->setColumnStatus(iCol, ClpSimplex::basic);
1552                         n++;
1553                    }
1554                    iCol = nextSlack[iCol];
1555               }
1556               iCol = negSlack[i];
1557               while (iCol >= 0) {
1558                    if (colsol[iCol] > lower[iCol] + 1.0e-8 &&
1559                              colsol[iCol] < upper[iCol] - 1.0e-8) {
1560                         model_->setColumnStatus(iCol, ClpSimplex::basic);
1561                         n++;
1562                    }
1563                    iCol = nextSlack[iCol];
1564               }
1565               if (n) {
1566                    if (fabs(rowsol[i] - rowlower[i]) < fabs(rowsol[i] - rowupper[i]))
1567                         model_->setRowStatus(i, ClpSimplex::atLowerBound);
1568                    else
1569                         model_->setRowStatus(i, ClpSimplex::atUpperBound);
1570#ifdef CLP_INVESTIGATE
1571                    if (n > 1)
1572                         printf("%d basic on row %d!\n", n, i);
1573#endif
1574               }
1575          }
1576#endif
1577     }
1578     double maxmin;
1579     if (model_->getObjSense() == -1.0) {
1580          maxmin = -1.0;
1581     } else {
1582          maxmin = 1.0;
1583     }
1584     bool justValuesPass = majorIterations_ > 1000000;
1585     if (slackStart >= 0) {
1586          for (i = 0; i < nrows; i++) {
1587               model_->setRowStatus(i, ClpSimplex::superBasic);
1588          }
1589          for (i = slackStart; i < slackEnd; i++) {
1590               model_->setColumnStatus(i, ClpSimplex::basic);
1591          }
1592     } else {
1593          /* still try and put singletons rather than artificials in basis */
1594          int ninbas = 0;
1595          for (i = 0; i < nrows; i++) {
1596               model_->setRowStatus(i, ClpSimplex::basic);
1597          }
1598          for (i = 0; i < ncols; i++) {
1599               if (columnLength[i] == 1 && upper[i] > lower[i] + 1.0e-5) {
1600                    CoinBigIndex j = columnStart[i];
1601                    double value = element[j];
1602                    int irow = row[j];
1603                    double rlo = rowlower[irow];
1604                    double rup = rowupper[irow];
1605                    double clo = lower[i];
1606                    double cup = upper[i];
1607                    double csol = colsol[i];
1608                    /* adjust towards feasibility */
1609                    double move = 0.0;
1610                    if (rowsol[irow] > rup) {
1611                         move = (rup - rowsol[irow]) / value;
1612                         if (value > 0.0) {
1613                              /* reduce */
1614                              if (csol + move < clo) move = CoinMin(0.0, clo - csol);
1615                         } else {
1616                              /* increase */
1617                              if (csol + move > cup) move = CoinMax(0.0, cup - csol);
1618                         }
1619                    } else if (rowsol[irow] < rlo) {
1620                         move = (rlo - rowsol[irow]) / value;
1621                         if (value > 0.0) {
1622                              /* increase */
1623                              if (csol + move > cup) move = CoinMax(0.0, cup - csol);
1624                         } else {
1625                              /* reduce */
1626                              if (csol + move < clo) move = CoinMin(0.0, clo - csol);
1627                         }
1628                    } else {
1629                         /* move to improve objective */
1630                         if (cost[i]*maxmin > 0.0) {
1631                              if (value > 0.0) {
1632                                   move = (rlo - rowsol[irow]) / value;
1633                                   /* reduce */
1634                                   if (csol + move < clo) move = CoinMin(0.0, clo - csol);
1635                              } else {
1636                                   move = (rup - rowsol[irow]) / value;
1637                                   /* increase */
1638                                   if (csol + move > cup) move = CoinMax(0.0, cup - csol);
1639                              }
1640                         } else if (cost[i]*maxmin < 0.0) {
1641                              if (value > 0.0) {
1642                                   move = (rup - rowsol[irow]) / value;
1643                                   /* increase */
1644                                   if (csol + move > cup) move = CoinMax(0.0, cup - csol);
1645                              } else {
1646                                   move = (rlo - rowsol[irow]) / value;
1647                                   /* reduce */
1648                                   if (csol + move < clo) move = CoinMin(0.0, clo - csol);
1649                              }
1650                         }
1651                    }
1652                    rowsol[irow] += move * value;
1653                    colsol[i] += move;
1654                    /* put in basis if row was artificial */
1655                    if (rup - rlo < 1.0e-7 && model_->getRowStatus(irow) == ClpSimplex::basic) {
1656                         model_->setRowStatus(irow, ClpSimplex::superBasic);
1657                         model_->setColumnStatus(i, ClpSimplex::basic);
1658                         ninbas++;
1659                    }
1660               }
1661          }
1662          /*printf("%d in basis\n",ninbas);*/
1663     }
1664     bool wantVector = false;
1665     if (dynamic_cast< ClpPackedMatrix*>(model_->clpMatrix())) {
1666          // See if original wanted vector
1667          ClpPackedMatrix * clpMatrixO = dynamic_cast< ClpPackedMatrix*>(model_->clpMatrix());
1668          wantVector = clpMatrixO->wantsSpecialColumnCopy();
1669     }
1670     if (addAll < 3) {
1671          ClpPresolve pinfo;
1672          if (presolve) {
1673               if (allowInfeasible) {
1674                    // fix up so will be feasible
1675                    double * rhs = new double[nrows];
1676                    memset(rhs, 0, nrows * sizeof(double));
1677                    model_->clpMatrix()->times(1.0, colsol, rhs);
1678                    double * rowupper = model_->rowUpper();
1679                    double * rowlower = model_->rowLower();
1680                    saveRowUpper = CoinCopyOfArray(rowupper, nrows);
1681                    saveRowLower = CoinCopyOfArray(rowlower, nrows);
1682                    double sum = 0.0;
1683                    for (i = 0; i < nrows; i++) {
1684                         if (rhs[i] > rowupper[i]) {
1685                              sum += rhs[i] - rowupper[i];
1686                              rowupper[i] = rhs[i];
1687                         }
1688                         if (rhs[i] < rowlower[i]) {
1689                              sum += rowlower[i] - rhs[i];
1690                              rowlower[i] = rhs[i];
1691                         }
1692                    }
1693                    printf("sum of infeasibilities %g\n", sum);
1694                    delete [] rhs;
1695               }
1696               saveModel = model_;
1697               pinfo.setPresolveActions(pinfo.presolveActions() | 16384);
1698               model_ = pinfo.presolvedModel(*model_, 1.0e-8, false, 5);
1699          }
1700          if (model_) {
1701               if (!wantVector) {
1702                    //#define TWO_GOES
1703#ifndef TWO_GOES
1704                    model_->primal(justValuesPass ? 2 : 1);
1705#else
1706                    model_->primal(1 + 11);
1707#endif
1708               } else {
1709                    ClpMatrixBase * matrix = model_->clpMatrix();
1710                    ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
1711                    assert (clpMatrix);
1712                    clpMatrix->makeSpecialColumnCopy();
1713                    model_->primal(1);
1714                    clpMatrix->releaseSpecialColumnCopy();
1715               }
1716               if (presolve) {
1717                    pinfo.postsolve(true);
1718                    delete model_;
1719                    model_ = saveModel;
1720                    saveModel = NULL;
1721               }
1722          } else {
1723               // not feasible
1724               addAll = 1;
1725               presolve = 0;
1726               model_ = saveModel;
1727               saveModel = NULL;
1728               if (justValuesPass)
1729                    model_->primal(2);
1730          }
1731          if (allowInfeasible) {
1732               CoinMemcpyN(saveRowUpper, nrows, model_->rowUpper());
1733               CoinMemcpyN(saveRowLower, nrows, model_->rowLower());
1734               delete [] saveRowUpper;
1735               delete [] saveRowLower;
1736               saveRowUpper = NULL;
1737               saveRowLower = NULL;
1738          }
1739          if (addAll < 2) {
1740               n = 0;
1741               if (!addAll ) {
1742                    /* could do scans to get a good number */
1743                    iteration = 1;
1744                    for (i = ordStart; i < ordEnd; i++) {
1745                         if (whenUsed[i] >= iteration) {
1746                              if (upper[i] - lower[i] < 1.0e-5 && saveUpper[i] - saveLower[i] > 1.0e-5) {
1747                                   n++;
1748                                   upper[i] = saveUpper[i];
1749                                   lower[i] = saveLower[i];
1750                              }
1751                         }
1752                    }
1753               } else {
1754                    for (i = ordStart; i < ordEnd; i++) {
1755                         if (upper[i] - lower[i] < 1.0e-5 && saveUpper[i] - saveLower[i] > 1.0e-5) {
1756                              n++;
1757                              upper[i] = saveUpper[i];
1758                              lower[i] = saveLower[i];
1759                         }
1760                    }
1761                    delete [] saveUpper;
1762                    delete [] saveLower;
1763                    saveUpper = NULL;
1764                    saveLower = NULL;
1765               }
1766#ifdef COIN_DEVELOP
1767               printf("Time so far %g, %d now added from previous iterations\n",
1768                      CoinCpuTime() - startTime, n);
1769#endif
1770               if (justValuesPass)
1771                    return;
1772               if (addAll)
1773                    presolve = 0;
1774               if (presolve) {
1775                    saveModel = model_;
1776                    model_ = pinfo.presolvedModel(*model_, 1.0e-8, false, 5);
1777               } else {
1778                    presolve = 0;
1779               }
1780               if (!wantVector) {
1781                    model_->primal(1);
1782               } else {
1783                    ClpMatrixBase * matrix = model_->clpMatrix();
1784                    ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
1785                    assert (clpMatrix);
1786                    clpMatrix->makeSpecialColumnCopy();
1787                    model_->primal(1);
1788                    clpMatrix->releaseSpecialColumnCopy();
1789               }
1790               if (presolve) {
1791                    pinfo.postsolve(true);
1792                    delete model_;
1793                    model_ = saveModel;
1794                    saveModel = NULL;
1795               }
1796               if (!addAll) {
1797                    n = 0;
1798                    for (i = ordStart; i < ordEnd; i++) {
1799                         if (upper[i] - lower[i] < 1.0e-5 && saveUpper[i] - saveLower[i] > 1.0e-5) {
1800                              n++;
1801                              upper[i] = saveUpper[i];
1802                              lower[i] = saveLower[i];
1803                         }
1804                    }
1805                    delete [] saveUpper;
1806                    delete [] saveLower;
1807                    saveUpper = NULL;
1808                    saveLower = NULL;
1809#ifdef COIN_DEVELOP
1810                    printf("Time so far %g, %d now added from previous iterations\n",
1811                           CoinCpuTime() - startTime, n);
1812#endif
1813               }
1814               if (presolve) {
1815                    saveModel = model_;
1816                    model_ = pinfo.presolvedModel(*model_, 1.0e-8, false, 5);
1817               } else {
1818                    presolve = 0;
1819               }
1820               if (!wantVector) {
1821                    model_->primal(1);
1822               } else {
1823                    ClpMatrixBase * matrix = model_->clpMatrix();
1824                    ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
1825                    assert (clpMatrix);
1826                    clpMatrix->makeSpecialColumnCopy();
1827                    model_->primal(1);
1828                    clpMatrix->releaseSpecialColumnCopy();
1829               }
1830               if (presolve) {
1831                    pinfo.postsolve(true);
1832                    delete model_;
1833                    model_ = saveModel;
1834                    saveModel = NULL;
1835               }
1836          }
1837#ifdef COIN_DEVELOP
1838          printf("Total time in crossover %g\n", CoinCpuTime() - startTime);
1839#endif
1840          delete [] saveUpper;
1841          delete [] saveLower;
1842     }
1843#ifdef FEB_TRY
1844     model_->setSpecialOptions(saveOptions);
1845     model_->setPerturbation(savePerturbation);
1846#endif
1847     return ;
1848}
1849#endif
1850/*****************************************************************************/
1851
1852// Default contructor
1853Idiot::Idiot()
1854{
1855     model_ = NULL;
1856     maxBigIts_ = 3;
1857     maxIts_ = 5;
1858     logLevel_ = 1;
1859     logFreq_ = 100;
1860     maxIts2_ = 100;
1861     djTolerance_ = 1e-1;
1862     mu_ = 1e-4;
1863     drop_ = 5.0;
1864     exitDrop_ = -1.0e20;
1865     muFactor_ = 0.3333;
1866     stopMu_ = 1e-12;
1867     smallInfeas_ = 1e-1;
1868     reasonableInfeas_ = 1e2;
1869     muAtExit_ = 1.0e31;
1870     strategy_ = 8;
1871     lambdaIterations_ = 0;
1872     checkFrequency_ = 100;
1873     whenUsed_ = NULL;
1874     majorIterations_ = 30;
1875     exitFeasibility_ = -1.0;
1876     dropEnoughFeasibility_ = 0.02;
1877     dropEnoughWeighted_ = 0.01;
1878     // adjust
1879     double nrows = 10000.0;
1880     int baseIts = static_cast<int> (sqrt(static_cast<double>(nrows)));
1881     baseIts = baseIts / 10;
1882     baseIts *= 10;
1883     maxIts2_ = 200 + baseIts + 5;
1884     maxIts2_ = 100;
1885     reasonableInfeas_ = static_cast<double> (nrows) * 0.05;
1886     lightWeight_ = 0;
1887}
1888// Constructor from model
1889Idiot::Idiot(OsiSolverInterface &model)
1890{
1891     model_ = & model;
1892     maxBigIts_ = 3;
1893     maxIts_ = 5;
1894     logLevel_ = 1;
1895     logFreq_ = 100;
1896     maxIts2_ = 100;
1897     djTolerance_ = 1e-1;
1898     mu_ = 1e-4;
1899     drop_ = 5.0;
1900     exitDrop_ = -1.0e20;
1901     muFactor_ = 0.3333;
1902     stopMu_ = 1e-12;
1903     smallInfeas_ = 1e-1;
1904     reasonableInfeas_ = 1e2;
1905     muAtExit_ = 1.0e31;
1906     strategy_ = 8;
1907     lambdaIterations_ = 0;
1908     checkFrequency_ = 100;
1909     whenUsed_ = NULL;
1910     majorIterations_ = 30;
1911     exitFeasibility_ = -1.0;
1912     dropEnoughFeasibility_ = 0.02;
1913     dropEnoughWeighted_ = 0.01;
1914     // adjust
1915     double nrows;
1916     if (model_)
1917          nrows = model_->getNumRows();
1918     else
1919          nrows = 10000.0;
1920     int baseIts = static_cast<int> (sqrt(static_cast<double>(nrows)));
1921     baseIts = baseIts / 10;
1922     baseIts *= 10;
1923     maxIts2_ = 200 + baseIts + 5;
1924     maxIts2_ = 100;
1925     reasonableInfeas_ = static_cast<double> (nrows) * 0.05;
1926     lightWeight_ = 0;
1927}
1928// Copy constructor.
1929Idiot::Idiot(const Idiot &rhs)
1930{
1931     model_ = rhs.model_;
1932     if (model_ && rhs.whenUsed_) {
1933          int numberColumns = model_->getNumCols();
1934          whenUsed_ = new int [numberColumns];
1935          CoinMemcpyN(rhs.whenUsed_, numberColumns, whenUsed_);
1936     } else {
1937          whenUsed_ = NULL;
1938     }
1939     djTolerance_ = rhs.djTolerance_;
1940     mu_ = rhs.mu_;
1941     drop_ = rhs.drop_;
1942     muFactor_ = rhs.muFactor_;
1943     stopMu_ = rhs.stopMu_;
1944     smallInfeas_ = rhs.smallInfeas_;
1945     reasonableInfeas_ = rhs.reasonableInfeas_;
1946     exitDrop_ = rhs.exitDrop_;
1947     muAtExit_ = rhs.muAtExit_;
1948     exitFeasibility_ = rhs.exitFeasibility_;
1949     dropEnoughFeasibility_ = rhs.dropEnoughFeasibility_;
1950     dropEnoughWeighted_ = rhs.dropEnoughWeighted_;
1951     maxBigIts_ = rhs.maxBigIts_;
1952     maxIts_ = rhs.maxIts_;
1953     majorIterations_ = rhs.majorIterations_;
1954     logLevel_ = rhs.logLevel_;
1955     logFreq_ = rhs.logFreq_;
1956     checkFrequency_ = rhs.checkFrequency_;
1957     lambdaIterations_ = rhs.lambdaIterations_;
1958     maxIts2_ = rhs.maxIts2_;
1959     strategy_ = rhs.strategy_;
1960     lightWeight_ = rhs.lightWeight_;
1961}
1962// Assignment operator. This copies the data
1963Idiot &
1964Idiot::operator=(const Idiot & rhs)
1965{
1966     if (this != &rhs) {
1967          delete [] whenUsed_;
1968          model_ = rhs.model_;
1969          if (model_ && rhs.whenUsed_) {
1970               int numberColumns = model_->getNumCols();
1971               whenUsed_ = new int [numberColumns];
1972               CoinMemcpyN(rhs.whenUsed_, numberColumns, whenUsed_);
1973          } else {
1974               whenUsed_ = NULL;
1975          }
1976          djTolerance_ = rhs.djTolerance_;
1977          mu_ = rhs.mu_;
1978          drop_ = rhs.drop_;
1979          muFactor_ = rhs.muFactor_;
1980          stopMu_ = rhs.stopMu_;
1981          smallInfeas_ = rhs.smallInfeas_;
1982          reasonableInfeas_ = rhs.reasonableInfeas_;
1983          exitDrop_ = rhs.exitDrop_;
1984          muAtExit_ = rhs.muAtExit_;
1985          exitFeasibility_ = rhs.exitFeasibility_;
1986          dropEnoughFeasibility_ = rhs.dropEnoughFeasibility_;
1987          dropEnoughWeighted_ = rhs.dropEnoughWeighted_;
1988          maxBigIts_ = rhs.maxBigIts_;
1989          maxIts_ = rhs.maxIts_;
1990          majorIterations_ = rhs.majorIterations_;
1991          logLevel_ = rhs.logLevel_;
1992          logFreq_ = rhs.logFreq_;
1993          checkFrequency_ = rhs.checkFrequency_;
1994          lambdaIterations_ = rhs.lambdaIterations_;
1995          maxIts2_ = rhs.maxIts2_;
1996          strategy_ = rhs.strategy_;
1997          lightWeight_ = rhs.lightWeight_;
1998     }
1999     return *this;
2000}
2001Idiot::~Idiot()
2002{
2003     delete [] whenUsed_;
2004}
Note: See TracBrowser for help on using the repository browser.