source: branches/pre/ClpInterior.cpp @ 214

Last change on this file since 214 was 214, checked in by forrest, 16 years ago

More pdco

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