source: branches/pre/ClpInterior.cpp @ 1926

Last change on this file since 1926 was 217, checked in by forrest, 17 years ago

Nearly there

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 11.3 KB
Line 
1// Copyright (C) 2002, International Business Machines
2// Corporation and others.  All Rights Reserved.
3
4
5
6
7#include "CoinPragma.hpp"
8
9#include <math.h>
10
11#include "CoinHelperFunctions.hpp"
12#include "ClpInterior.hpp"
13#include "ClpMatrixBase.hpp"
14#include "ClpLsqr.hpp"
15#include "ClpPdcoBase.hpp"
16#include "CoinDenseVector.hpp"
17#include "ClpMessage.hpp"
18#include "ClpLinearObjective.hpp"
19#include <cfloat>
20
21#include <string>
22#include <stdio.h>
23#include <iostream>
24//#############################################################################
25
26ClpInterior::ClpInterior () :
27
28  ClpModel(),
29  largestPrimalError_(0.0),
30  largestDualError_(0.0),
31  sumDualInfeasibilities_(0.0),
32  sumPrimalInfeasibilities_(0.0),
33  xsize_(0.0),
34  zsize_(0.0),
35  lower_(NULL),
36  rowLowerWork_(NULL),
37  columnLowerWork_(NULL),
38  upper_(NULL),
39  rowUpperWork_(NULL),
40  columnUpperWork_(NULL),
41  cost_(NULL),
42  rhs_(NULL),
43   x_(NULL),
44   y_(NULL),
45  dj_(NULL),
46  lsqrObject_(NULL),
47  pdcoStuff_(NULL)
48{
49  solveType_=2; // say interior based life form
50}
51
52// Subproblem constructor
53ClpInterior::ClpInterior ( const ClpModel * rhs,
54                           int numberRows, const int * whichRow,
55                           int numberColumns, const int * whichColumn,
56                           bool dropNames, bool dropIntegers)
57  : ClpModel(rhs, numberRows, whichRow,
58             numberColumns,whichColumn,dropNames,dropIntegers),
59    largestPrimalError_(0.0),
60    largestDualError_(0.0),
61    sumDualInfeasibilities_(0.0),
62    sumPrimalInfeasibilities_(0.0),
63    xsize_(0.0),
64    zsize_(0.0),
65    lower_(NULL),
66    rowLowerWork_(NULL),
67    columnLowerWork_(NULL),
68    upper_(NULL),
69    rowUpperWork_(NULL),
70    columnUpperWork_(NULL),
71    cost_(NULL),
72    rhs_(NULL),
73    x_(NULL),
74    y_(NULL),
75    dj_(NULL),
76    lsqrObject_(NULL),
77    pdcoStuff_(NULL)
78{
79  solveType_=2; // say interior based life form
80 
81}
82
83//-----------------------------------------------------------------------------
84
85ClpInterior::~ClpInterior ()
86{
87  gutsOfDelete();
88}
89//#############################################################################
90/*
91   This does housekeeping
92*/
93int 
94ClpInterior::housekeeping()
95{
96  numberIterations_++;
97  return 0;
98}
99// Copy constructor.
100ClpInterior::ClpInterior(const ClpInterior &rhs) :
101  ClpModel(rhs),
102  largestPrimalError_(0.0),
103  largestDualError_(0.0),
104  sumDualInfeasibilities_(0.0),
105  sumPrimalInfeasibilities_(0.0),
106  xsize_(0.0),
107  zsize_(0.0),
108  lower_(NULL),
109  rowLowerWork_(NULL),
110  columnLowerWork_(NULL),
111  upper_(NULL),
112  rowUpperWork_(NULL),
113  columnUpperWork_(NULL),
114  cost_(NULL),
115  rhs_(NULL),
116   x_(NULL),
117   y_(NULL),
118  dj_(NULL),
119  lsqrObject_(NULL),
120  pdcoStuff_(NULL)
121{
122  gutsOfDelete();
123  gutsOfCopy(rhs);
124  solveType_=2; // say interior based life form
125}
126// Copy constructor from model
127ClpInterior::ClpInterior(const ClpModel &rhs) :
128  ClpModel(rhs),
129  largestPrimalError_(0.0),
130  largestDualError_(0.0),
131  sumDualInfeasibilities_(0.0),
132  sumPrimalInfeasibilities_(0.0),
133  xsize_(0.0),
134  zsize_(0.0),
135  lower_(NULL),
136  rowLowerWork_(NULL),
137  columnLowerWork_(NULL),
138  upper_(NULL),
139  rowUpperWork_(NULL),
140  columnUpperWork_(NULL),
141  cost_(NULL),
142  rhs_(NULL),
143   x_(NULL),
144   y_(NULL),
145  dj_(NULL),
146  lsqrObject_(NULL),
147  pdcoStuff_(NULL)
148{
149  solveType_=2; // say interior based life form
150}
151// Assignment operator. This copies the data
152ClpInterior & 
153ClpInterior::operator=(const ClpInterior & rhs)
154{
155  if (this != &rhs) {
156    gutsOfDelete();
157    ClpModel::operator=(rhs);
158    gutsOfCopy(rhs);
159  }
160  return *this;
161}
162void 
163ClpInterior::gutsOfCopy(const ClpInterior & rhs)
164{
165  lower_ = ClpCopyOfArray(rhs.lower_,numberColumns_+numberRows_);
166  rowLowerWork_ = lower_+numberColumns_;
167  columnLowerWork_ = lower_;
168  upper_ = ClpCopyOfArray(rhs.upper_,numberColumns_+numberRows_);
169  rowUpperWork_ = upper_+numberColumns_;
170  columnUpperWork_ = upper_;
171  //cost_ = ClpCopyOfArray(rhs.cost_,2*(numberColumns_+numberRows_));
172  cost_ = ClpCopyOfArray(rhs.cost_,numberColumns_);
173  rhs_ = ClpCopyOfArray(rhs.rhs_,numberRows_);
174   x_ = ClpCopyOfArray(rhs.x_,numberColumns_);
175   y_ = ClpCopyOfArray(rhs.y_,numberRows_);
176  dj_ = ClpCopyOfArray(rhs.dj_,numberColumns_);
177  lsqrObject_= new ClpLsqr(*rhs.lsqrObject_);
178  pdcoStuff_ = rhs.pdcoStuff_->clone();
179  largestPrimalError_ = rhs.largestPrimalError_;
180  largestDualError_ = rhs.largestDualError_;
181  sumDualInfeasibilities_ = rhs.sumDualInfeasibilities_;
182  sumPrimalInfeasibilities_ = rhs.sumPrimalInfeasibilities_;
183  xsize_ = rhs.xsize_;
184  zsize_ = rhs.zsize_;
185  solveType_=rhs.solveType_;
186}
187// type == 0 do everything, most + pivot data, 2 factorization data as well
188void 
189ClpInterior::gutsOfDelete()
190{
191  delete [] lower_;
192  lower_=NULL;
193  rowLowerWork_=NULL;
194  columnLowerWork_=NULL;
195  delete [] upper_;
196  upper_=NULL;
197  rowUpperWork_=NULL;
198  columnUpperWork_=NULL;
199  delete [] cost_;
200  cost_=NULL;
201  delete [] rhs_;
202  rhs_ = NULL;
203  delete [] x_;
204  x_ = NULL;
205  delete [] y_;
206  y_ = NULL;
207  delete [] dj_;
208  dj_ = NULL;
209  delete lsqrObject_;
210  lsqrObject_ = NULL;
211  //delete pdcoStuff_;
212  pdcoStuff_=NULL;
213}
214bool
215ClpInterior::createWorkingData()
216{
217  bool goodMatrix=true;
218  //check matrix
219  if (!matrix_->allElementsInRange(this,1.0e-12,1.0e20)) {
220    problemStatus_=4;
221    goodMatrix= false;
222  }
223  // check rim of problem okay
224  if (!sanityCheck())
225    goodMatrix=false;
226  return goodMatrix;
227}
228void
229ClpInterior::deleteWorkingData()
230{
231}
232// Sanity check on input data - returns true if okay
233bool 
234ClpInterior::sanityCheck()
235{
236  // bad if empty
237  if (!numberRows_||!numberColumns_||!matrix_->getNumElements()) {
238    handler_->message(CLP_EMPTY_PROBLEM,messages_)
239      <<numberRows_
240      <<numberColumns_
241      <<matrix_->getNumElements()
242      <<CoinMessageEol;
243    problemStatus_=4;
244    return false;
245  }
246  int numberBad ;
247  double largestBound, smallestBound, minimumGap;
248  double smallestObj, largestObj;
249  int firstBad;
250  int modifiedBounds=0;
251  int i;
252  numberBad=0;
253  firstBad=-1;
254  minimumGap=1.0e100;
255  smallestBound=1.0e100;
256  largestBound=0.0;
257  smallestObj=1.0e100;
258  largestObj=0.0;
259  // If bounds are too close - fix
260  double fixTolerance = 1.1*primalTolerance();
261  for (i=numberColumns_;i<numberColumns_+numberRows_;i++) {
262    double value;
263    value = fabs(cost_[i]);
264    if (value>1.0e50) {
265      numberBad++;
266      if (firstBad<0)
267        firstBad=i;
268    } else if (value) {
269      if (value>largestObj)
270        largestObj=value;
271      if (value<smallestObj)
272        smallestObj=value;
273    }
274    value=upper_[i]-lower_[i];
275    if (value<-primalTolerance()) {
276      numberBad++;
277      if (firstBad<0)
278        firstBad=i;
279    } else if (value<=fixTolerance) {
280      if (value) {
281        // modify
282        upper_[i] = lower_[i];
283        modifiedBounds++;
284      }
285    } else {
286      if (value<minimumGap)
287        minimumGap=value;
288    }
289    if (lower_[i]>-1.0e100&&lower_[i]) {
290      value = fabs(lower_[i]);
291      if (value>largestBound)
292        largestBound=value;
293      if (value<smallestBound)
294        smallestBound=value;
295    }
296    if (upper_[i]<1.0e100&&upper_[i]) {
297      value = fabs(upper_[i]);
298      if (value>largestBound)
299        largestBound=value;
300      if (value<smallestBound)
301        smallestBound=value;
302    }
303  }
304  if (largestBound)
305    handler_->message(CLP_RIMSTATISTICS3,messages_)
306      <<smallestBound
307      <<largestBound
308      <<minimumGap
309      <<CoinMessageEol;
310  minimumGap=1.0e100;
311  smallestBound=1.0e100;
312  largestBound=0.0;
313  for (i=0;i<numberColumns_;i++) {
314    double value;
315    value = fabs(cost_[i]);
316    if (value>1.0e50) {
317      numberBad++;
318      if (firstBad<0)
319        firstBad=i;
320    } else if (value) {
321      if (value>largestObj)
322        largestObj=value;
323      if (value<smallestObj)
324        smallestObj=value;
325    }
326    value=upper_[i]-lower_[i];
327    if (value<-primalTolerance()) {
328      numberBad++;
329      if (firstBad<0)
330        firstBad=i;
331    } else if (value<=fixTolerance) {
332      if (value) {
333        // modify
334        upper_[i] = lower_[i];
335        modifiedBounds++;
336      }
337    } else {
338      if (value<minimumGap)
339        minimumGap=value;
340    }
341    if (lower_[i]>-1.0e100&&lower_[i]) {
342      value = fabs(lower_[i]);
343      if (value>largestBound)
344        largestBound=value;
345      if (value<smallestBound)
346        smallestBound=value;
347    }
348    if (upper_[i]<1.0e100&&upper_[i]) {
349      value = fabs(upper_[i]);
350      if (value>largestBound)
351        largestBound=value;
352      if (value<smallestBound)
353        smallestBound=value;
354    }
355  }
356  char rowcol[]={'R','C'};
357  if (numberBad) {
358    handler_->message(CLP_BAD_BOUNDS,messages_)
359      <<numberBad
360      <<rowcol[isColumn(firstBad)]<<sequenceWithin(firstBad)
361      <<CoinMessageEol;
362    problemStatus_=4;
363    return false;
364  }
365  if (modifiedBounds)
366    handler_->message(CLP_MODIFIEDBOUNDS,messages_)
367      <<modifiedBounds
368      <<CoinMessageEol;
369  handler_->message(CLP_RIMSTATISTICS1,messages_)
370    <<smallestObj
371    <<largestObj
372    <<CoinMessageEol;  if (largestBound)
373    handler_->message(CLP_RIMSTATISTICS2,messages_)
374      <<smallestBound
375      <<largestBound
376      <<minimumGap
377      <<CoinMessageEol;
378  return true;
379}
380/* Loads a problem (the constraints on the
381   rows are given by lower and upper bounds). If a pointer is 0 then the
382   following values are the default:
383   <ul>
384   <li> <code>colub</code>: all columns have upper bound infinity
385   <li> <code>collb</code>: all columns have lower bound 0
386   <li> <code>rowub</code>: all rows have upper bound infinity
387   <li> <code>rowlb</code>: all rows have lower bound -infinity
388   <li> <code>obj</code>: all variables have 0 objective coefficient
389   </ul>
390*/
391void 
392ClpInterior::loadProblem (  const ClpMatrixBase& matrix,
393                    const double* collb, const double* colub,   
394                    const double* obj,
395                    const double* rowlb, const double* rowub,
396                    const double * rowObjective)
397{
398  ClpModel::loadProblem(matrix, collb, colub, obj, rowlb, rowub,
399                        rowObjective);
400}
401void 
402ClpInterior::loadProblem (  const CoinPackedMatrix& matrix,
403                    const double* collb, const double* colub,   
404                    const double* obj,
405                    const double* rowlb, const double* rowub,
406                    const double * rowObjective)
407{
408  ClpModel::loadProblem(matrix, collb, colub, obj, rowlb, rowub,
409                        rowObjective);
410}
411
412/* Just like the other loadProblem() method except that the matrix is
413   given in a standard column major ordered format (without gaps). */
414void 
415ClpInterior::loadProblem (  const int numcols, const int numrows,
416                    const CoinBigIndex* start, const int* index,
417                    const double* value,
418                    const double* collb, const double* colub,   
419                    const double* obj,
420                    const double* rowlb, const double* rowub,
421                    const double * rowObjective)
422{
423  ClpModel::loadProblem(numcols, numrows, start, index, value,
424                          collb, colub, obj, rowlb, rowub,
425                          rowObjective);
426}
427void 
428ClpInterior::loadProblem (  const int numcols, const int numrows,
429                           const CoinBigIndex* start, const int* index,
430                           const double* value,const int * length,
431                           const double* collb, const double* colub,   
432                           const double* obj,
433                           const double* rowlb, const double* rowub,
434                           const double * rowObjective)
435{
436  ClpModel::loadProblem(numcols, numrows, start, index, value, length,
437                          collb, colub, obj, rowlb, rowub,
438                          rowObjective);
439}
440// Read an mps file from the given filename
441int 
442ClpInterior::readMps(const char *filename,
443            bool keepNames,
444            bool ignoreErrors)
445{
446  int status = ClpModel::readMps(filename,keepNames,ignoreErrors);
447  return status;
448}
449#include "ClpPdco.hpp"
450/* Pdco algorithm - see ClpPdco.hpp for method */
451int 
452ClpInterior::pdco()
453{
454  return ((ClpPdco *) this)->pdco();
455}
456// ** Temporary version
457int 
458ClpInterior::pdco( ClpPdcoBase * stuff, Options &options, Info &info, Outfo &outfo)
459{
460  return ((ClpPdco *) this)->pdco(stuff,options,info,outfo);
461}
Note: See TracBrowser for help on using the repository browser.