source: trunk/ClpQuadraticObjective.cpp @ 225

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

This should break everything

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 12.3 KB
Line 
1// Copyright (C) 2003, International Business Machines
2// Corporation and others.  All Rights Reserved.
3
4#include "CoinPragma.hpp"
5#include "ClpModel.hpp"
6#include "ClpQuadraticObjective.hpp"
7
8//#############################################################################
9// Constructors / Destructor / Assignment
10//#############################################################################
11
12//-------------------------------------------------------------------
13// Default Constructor
14//-------------------------------------------------------------------
15ClpQuadraticObjective::ClpQuadraticObjective () 
16: ClpObjective()
17{
18  type_=2;
19  objective_=NULL;
20  quadraticObjective_=NULL;
21  gradient_ = NULL;
22  numberColumns_=0;
23  numberExtendedColumns_=0;
24}
25
26//-------------------------------------------------------------------
27// Useful Constructor
28//-------------------------------------------------------------------
29ClpQuadraticObjective::ClpQuadraticObjective (const double * objective , 
30                                              int numberColumns,
31                                              const CoinBigIndex * start,
32                                              const int * column, const double * element,
33                                              int numberExtendedColumns)
34  : ClpObjective()
35{
36  type_=2;
37  numberColumns_ = numberColumns;
38  if (numberExtendedColumns>=0)
39    numberExtendedColumns_= max(numberColumns_,numberExtendedColumns);
40  else
41    numberExtendedColumns_= numberColumns_;
42  if (objective) {
43    objective_ = new double [numberExtendedColumns_];
44    memcpy(objective_,objective,numberColumns_*sizeof(double));
45    memset(objective_+numberColumns_,0,(numberExtendedColumns_-numberColumns_)*sizeof(double));
46  } else {
47    objective_ = new double [numberExtendedColumns_];
48    memset(objective_,0,numberExtendedColumns_*sizeof(double));
49  }
50  if (start) 
51    quadraticObjective_ = new CoinPackedMatrix(true,numberColumns,numberColumns,
52                                             start[numberColumns],element,column,start,NULL);
53  else
54  quadraticObjective_=NULL;
55  gradient_ = NULL;
56}
57
58//-------------------------------------------------------------------
59// Copy constructor
60//-------------------------------------------------------------------
61ClpQuadraticObjective::ClpQuadraticObjective (const ClpQuadraticObjective & rhs) 
62: ClpObjective(rhs)
63{ 
64  numberColumns_=rhs.numberColumns_;
65  numberExtendedColumns_=rhs.numberExtendedColumns_;
66  if (rhs.objective_) {
67    objective_ = new double [numberExtendedColumns_];
68    memcpy(objective_,rhs.objective_,numberExtendedColumns_*sizeof(double));
69  } else {
70    objective_=NULL;
71  }
72  if (rhs.gradient_) {
73    gradient_ = new double [numberExtendedColumns_];
74    memcpy(gradient_,rhs.gradient_,numberExtendedColumns_*sizeof(double));
75  } else {
76    gradient_=NULL;
77  }
78  if (rhs.quadraticObjective_) 
79    quadraticObjective_ = new CoinPackedMatrix(*rhs.quadraticObjective_);
80  else 
81    quadraticObjective_=NULL;
82}
83/* Subset constructor.  Duplicates are allowed
84   and order is as given.
85*/
86ClpQuadraticObjective::ClpQuadraticObjective (const ClpQuadraticObjective &rhs,
87                                        int numberColumns, 
88                                        const int * whichColumn) 
89: ClpObjective(rhs)
90{
91  objective_=NULL;
92  int extra = rhs.numberExtendedColumns_-rhs.numberColumns_;
93  numberColumns_=0;
94  if (numberColumns>0) {
95    // check valid lists
96    int numberBad=0;
97    int i;
98    for (i=0;i<numberColumns;i++)
99      if (whichColumn[i]<0||whichColumn[i]>=rhs.numberColumns_)
100        numberBad++;
101    if (numberBad)
102      throw CoinError("bad column list", "subset constructor", 
103                      "ClpQuadraticObjective");
104    numberColumns_ = numberColumns;
105    numberExtendedColumns_ = numberColumns+extra;
106    objective_ = new double[numberExtendedColumns_];
107    for (i=0;i<numberColumns_;i++) 
108      objective_[i]=rhs.objective_[whichColumn[i]];
109    memcpy(objective_+numberColumns_,rhs.objective_+rhs.numberColumns_,
110           (numberExtendedColumns_-numberColumns_)*sizeof(double));
111    if (rhs.gradient_) {
112      gradient_ = new double[numberExtendedColumns_];
113      for (i=0;i<numberColumns_;i++) 
114        gradient_[i]=rhs.gradient_[whichColumn[i]];
115      memcpy(gradient_+numberColumns_,rhs.gradient_+rhs.numberColumns_,
116             (numberExtendedColumns_-numberColumns_)*sizeof(double));
117    }
118  }
119  if (rhs.quadraticObjective_) {
120    quadraticObjective_ = new CoinPackedMatrix(*rhs.quadraticObjective_,
121                                               numberColumns,whichColumn,
122                                               numberColumns,whichColumn);
123  } else {
124    quadraticObjective_=NULL;
125  }
126}
127 
128
129//-------------------------------------------------------------------
130// Destructor
131//-------------------------------------------------------------------
132ClpQuadraticObjective::~ClpQuadraticObjective ()
133{
134  delete [] objective_;
135  delete [] gradient_;
136  delete quadraticObjective_;
137}
138
139//----------------------------------------------------------------
140// Assignment operator
141//-------------------------------------------------------------------
142ClpQuadraticObjective &
143ClpQuadraticObjective::operator=(const ClpQuadraticObjective& rhs)
144{
145  if (this != &rhs) {
146    delete quadraticObjective_;
147    quadraticObjective_ = NULL;
148    ClpObjective::operator=(rhs);
149    numberColumns_=rhs.numberColumns_;
150    numberExtendedColumns_=rhs.numberExtendedColumns_;
151    if (rhs.objective_) {
152      objective_ = new double [numberExtendedColumns_];
153      memcpy(objective_,rhs.objective_,numberExtendedColumns_*sizeof(double));
154    } else {
155      objective_=NULL;
156    }
157    if (rhs.gradient_) {
158      gradient_ = new double [numberExtendedColumns_];
159      memcpy(gradient_,rhs.gradient_,numberExtendedColumns_*sizeof(double));
160    } else {
161      gradient_=NULL;
162    }
163    if (rhs.quadraticObjective_) {
164      quadraticObjective_ = new CoinPackedMatrix(*rhs.quadraticObjective_);
165    } else {
166      quadraticObjective_=NULL;
167    }
168  }
169  return *this;
170}
171
172// Returns gradient
173double * 
174ClpQuadraticObjective::gradient(const double * solution, double & offset,bool refresh)
175{
176  offset=0.0;
177  if (!quadraticObjective_||!solution) {
178    return objective_;
179  } else {
180    if (refresh||!gradient_) {
181      if (!gradient_) 
182        gradient_ = new double[numberExtendedColumns_];
183      const int * columnQuadratic = quadraticObjective_->getIndices();
184      const int * columnQuadraticStart = quadraticObjective_->getVectorStarts();
185      const int * columnQuadraticLength = quadraticObjective_->getVectorLengths();
186      const double * quadraticElement = quadraticObjective_->getElements();
187      double offset=0.0;
188      memcpy(gradient_,objective_,numberExtendedColumns_*sizeof(double));
189      int iColumn;
190      for (iColumn=0;iColumn<numberColumns_;iColumn++) {
191        double valueI = solution[iColumn];
192        int j;
193        for (j=columnQuadraticStart[iColumn];
194             j<columnQuadraticStart[iColumn]+columnQuadraticLength[iColumn];j++) {
195          int jColumn = columnQuadratic[j];
196          double valueJ = solution[jColumn];
197          double elementValue = quadraticElement[j];
198          elementValue *= 0.5;
199          offset += valueI*valueJ*elementValue;
200          double gradientI = valueJ*elementValue;
201          double gradientJ = valueI*elementValue;
202          offset -= gradientI*valueI;
203          gradient_[iColumn] += gradientI;
204          offset -= gradientJ*valueJ;
205          gradient_[jColumn] += gradientJ;
206        }
207      }
208    }
209    return gradient_;
210  }
211}
212 
213//-------------------------------------------------------------------
214// Clone
215//-------------------------------------------------------------------
216ClpObjective * ClpQuadraticObjective::clone() const
217{
218  return new ClpQuadraticObjective(*this);
219}
220/* Subset clone.  Duplicates are allowed
221   and order is as given.
222*/
223ClpObjective * 
224ClpQuadraticObjective::subsetClone (int numberColumns, 
225                           const int * whichColumns) const
226{
227  return new ClpQuadraticObjective(*this, numberColumns, whichColumns);
228}
229// Resize objective
230void 
231ClpQuadraticObjective::resize(int newNumberColumns)
232{
233  if (numberColumns_!=newNumberColumns) {
234    int newExtended = newNumberColumns + (numberExtendedColumns_-numberColumns_);
235    int i;
236    double * newArray = new double[newExtended];
237    if (objective_)
238      memcpy(newArray,objective_,
239             min(newExtended,numberExtendedColumns_)*sizeof(double));
240    delete [] objective_;
241    objective_ = newArray;
242    for (i=numberColumns_;i<newNumberColumns;i++) 
243      objective_[i]=0.0;
244    if (gradient_) {
245      newArray = new double[newExtended];
246      if (gradient_)
247        memcpy(newArray,gradient_,
248               min(newExtended,numberExtendedColumns_)*sizeof(double));
249      delete [] gradient_;
250      gradient_ = newArray;
251      for (i=numberColumns_;i<newNumberColumns;i++) 
252        gradient_[i]=0.0;
253    }
254    if (quadraticObjective_) {
255      if (newNumberColumns<numberColumns_) {
256        int * which = new int[numberColumns_-newNumberColumns];
257        int i;
258        for (i=newNumberColumns;i<numberColumns_;i++) 
259          which[i-newNumberColumns]=i;
260        quadraticObjective_->deleteCols(numberColumns_-newNumberColumns,which);
261        quadraticObjective_->deleteRows(numberColumns_-newNumberColumns,which);
262        delete [] which;
263      } else {
264        quadraticObjective_->setDimensions(newNumberColumns,newNumberColumns);
265      }
266    }
267    numberColumns_ = newNumberColumns;
268    numberExtendedColumns_ = newExtended;
269  } 
270 
271}
272// Delete columns in  objective
273void 
274ClpQuadraticObjective::deleteSome(int numberToDelete, const int * which) 
275{
276  int newNumberColumns = numberColumns_-numberToDelete;
277  int newExtended = numberExtendedColumns_ - numberToDelete;
278  if (objective_) {
279    int i ;
280    char * deleted = new char[numberColumns_];
281    int numberDeleted=0;
282    memset(deleted,0,numberColumns_*sizeof(char));
283    for (i=0;i<numberToDelete;i++) {
284      int j = which[i];
285      if (j>=0&&j<numberColumns_&&!deleted[j]) {
286        numberDeleted++;
287        deleted[j]=1;
288      }
289    }
290    newNumberColumns = numberColumns_-numberDeleted;
291    newExtended = numberExtendedColumns_ - numberDeleted;
292    double * newArray = new double[newExtended];
293    int put=0;
294    for (i=0;i<numberColumns_;i++) {
295      if (!deleted[i]) {
296        newArray[put++]=objective_[i];
297      }
298    }
299    delete [] objective_;
300    objective_ = newArray;
301    delete [] deleted;
302    memcpy(objective_+newNumberColumns,objective_+numberColumns_,
303           (numberExtendedColumns_-numberColumns_)*sizeof(double));
304  }
305  if (gradient_) {
306    int i ;
307    char * deleted = new char[numberColumns_];
308    int numberDeleted=0;
309    memset(deleted,0,numberColumns_*sizeof(char));
310    for (i=0;i<numberToDelete;i++) {
311      int j = which[i];
312      if (j>=0&&j<numberColumns_&&!deleted[j]) {
313        numberDeleted++;
314        deleted[j]=1;
315      }
316    }
317    newNumberColumns = numberColumns_-numberDeleted;
318    newExtended = numberExtendedColumns_ - numberDeleted;
319    double * newArray = new double[newExtended];
320    int put=0;
321    for (i=0;i<numberColumns_;i++) {
322      if (!deleted[i]) {
323        newArray[put++]=gradient_[i];
324      }
325    }
326    delete [] gradient_;
327    gradient_ = newArray;
328    delete [] deleted;
329    memcpy(gradient_+newNumberColumns,gradient_+numberColumns_,
330           (numberExtendedColumns_-numberColumns_)*sizeof(double));
331  }
332  numberColumns_ = newNumberColumns;
333  numberExtendedColumns_ = newExtended;
334  if (quadraticObjective_) {
335    quadraticObjective_->deleteCols(numberToDelete,which);
336    quadraticObjective_->deleteRows(numberToDelete,which);
337  }
338}
339
340// Load up quadratic objective
341void 
342ClpQuadraticObjective::loadQuadraticObjective(const int numberColumns, const CoinBigIndex * start,
343                              const int * column, const double * element,int numberExtended)
344{
345  delete quadraticObjective_;
346  quadraticObjective_ = new CoinPackedMatrix(true,numberColumns,numberColumns,
347                                             start[numberColumns],element,column,start,NULL);
348  numberColumns_=numberColumns;
349  if (numberExtended>numberExtendedColumns_) {
350    if (objective_) {
351      // make correct size
352      double * newArray = new double[numberExtended];
353      memcpy(newArray,objective_,numberColumns_*sizeof(double));
354      delete [] objective_;
355      objective_ = newArray;
356      memset(objective_+numberColumns_,0,(numberExtended-numberColumns_)*sizeof(double));
357    }
358    if (gradient_) {
359      // make correct size
360      double * newArray = new double[numberExtended];
361      memcpy(newArray,gradient_,numberColumns_*sizeof(double));
362      delete [] gradient_;
363      gradient_ = newArray;
364      memset(gradient_+numberColumns_,0,(numberExtended-numberColumns_)*sizeof(double));
365    }
366    numberExtendedColumns_ = numberExtended;
367  } else {
368    numberExtendedColumns_ = numberColumns_;
369  }
370}
371void 
372ClpQuadraticObjective::loadQuadraticObjective (  const CoinPackedMatrix& matrix)
373{
374  delete quadraticObjective_;
375  quadraticObjective_ = new CoinPackedMatrix(matrix);
376}
377// Get rid of quadratic objective
378void 
379ClpQuadraticObjective::deleteQuadraticObjective()
380{
381  delete quadraticObjective_;
382  quadraticObjective_ = NULL;
383}
Note: See TracBrowser for help on using the repository browser.