source: stable/1.15/Clp/src/AbcNonLinearCost.hpp @ 1989

Last change on this file since 1989 was 1910, checked in by stefan, 7 years ago
  • add configure option --enable-aboca={1,2,3,4,yes,no}
  • compile Aboca source only if --enable-aboca set (instead of compiling empty source files)
  • fix svn properties
  • Property svn:keywords set to Id
File size: 8.8 KB
Line 
1/* $Id: AbcNonLinearCost.hpp 1910 2013-01-27 02:00:13Z forrest $ */
2// Copyright (C) 2002, International Business Machines
3// Corporation and others, Copyright (C) 2012, FasterCoin.  All Rights Reserved.
4// This code is licensed under the terms of the Eclipse Public License (EPL).
5
6#ifndef AbcNonLinearCost_H
7#define AbcNonLinearCost_H
8
9
10#include "CoinPragma.hpp"
11#include "AbcCommon.hpp"
12
13class AbcSimplex;
14class CoinIndexedVector;
15
16/** Trivial class to deal with non linear costs
17   
18    I don't make any explicit assumptions about convexity but I am
19    sure I do make implicit ones.
20   
21    One interesting idea for normal LP's will be to allow non-basic
22    variables to come into basis as infeasible i.e. if variable at
23    lower bound has very large positive reduced cost (when problem
24    is infeasible) could it reduce overall problem infeasibility more
25    by bringing it into basis below its lower bound.
26   
27    Another feature would be to automatically discover when problems
28    are convex piecewise linear and re-formulate to use non-linear.
29    I did some work on this many years ago on "grade" problems, but
30    while it improved primal interior point algorithms were much better
31    for that particular problem.
32*/
33/* status has original status and current status
34   0 - below lower so stored is upper
35   1 - in range
36   2 - above upper so stored is lower
37   4 - (for current) - same as original
38*/
39#define CLP_BELOW_LOWER 0
40#define CLP_FEASIBLE 1
41#define CLP_ABOVE_UPPER 2
42#define CLP_SAME 4
43inline int originalStatus(unsigned char status)
44{
45  return (status & 15);
46}
47inline int currentStatus(unsigned char status)
48{
49  return (status >> 4);
50}
51inline void setOriginalStatus(unsigned char & status, int value)
52{
53  status = static_cast<unsigned char>(status & ~15);
54  status = static_cast<unsigned char>(status | value);
55}
56inline void setCurrentStatus(unsigned char &status, int value)
57{
58  status = static_cast<unsigned char>(status & ~(15 << 4));
59  status = static_cast<unsigned char>(status | (value << 4));
60}
61inline void setInitialStatus(unsigned char &status)
62{
63  status = static_cast<unsigned char>(CLP_FEASIBLE | (CLP_SAME << 4));
64}
65inline void setSameStatus(unsigned char &status)
66{
67  status = static_cast<unsigned char>(status & ~(15 << 4));
68  status = static_cast<unsigned char>(status | (CLP_SAME << 4));
69}
70class AbcNonLinearCost  {
71 
72public:
73 
74  /**@name Constructors, destructor */
75  //@{
76  /// Default constructor.
77  AbcNonLinearCost();
78  /** Constructor from simplex.
79      This will just set up wasteful arrays for linear, but
80      later may do dual analysis and even finding duplicate columns .
81  */
82  AbcNonLinearCost(AbcSimplex * model);
83  /// Destructor
84  ~AbcNonLinearCost();
85  // Copy
86  AbcNonLinearCost(const AbcNonLinearCost&);
87  // Assignment
88  AbcNonLinearCost& operator=(const AbcNonLinearCost&);
89  //@}
90 
91 
92  /**@name Actual work in primal */
93  //@{
94  /** Changes infeasible costs and computes number and cost of infeas
95      Puts all non-basic (non free) variables to bounds
96      and all free variables to zero if oldTolerance is non-zero
97      - but does not move those <= oldTolerance away*/
98  void checkInfeasibilities(double oldTolerance = 0.0);
99  /** Changes infeasible costs for each variable
100      The indices are row indices and need converting to sequences
101  */
102  void checkInfeasibilities(int numberInArray, const int * index);
103  /** Puts back correct infeasible costs for each variable
104      The input indices are row indices and need converting to sequences
105      for costs.
106      On input array is empty (but indices exist).  On exit just
107      changed costs will be stored as normal CoinIndexedVector
108  */
109  void checkChanged(int numberInArray, CoinIndexedVector * update);
110  /** Goes through one bound for each variable.
111      If multiplier*work[iRow]>0 goes down, otherwise up.
112      The indices are row indices and need converting to sequences
113      Temporary offsets may be set
114      Rhs entries are increased
115  */
116  void goThru(int numberInArray, double multiplier,
117              const int * index, const double * work,
118              double * rhs);
119  /** Takes off last iteration (i.e. offsets closer to 0)
120   */
121  void goBack(int numberInArray, const int * index,
122              double * rhs);
123  /** Puts back correct infeasible costs for each variable
124      The input indices are row indices and need converting to sequences
125      for costs.
126      At the end of this all temporary offsets are zero
127  */
128  void goBackAll(const CoinIndexedVector * update);
129  /// Temporary zeroing of feasible costs
130  void zapCosts();
131  /// Refreshes costs always makes row costs zero
132  void refreshCosts(const double * columnCosts);
133  /// Puts feasible bounds into lower and upper
134  void feasibleBounds();
135  /// Refresh - assuming regions OK
136  void refresh();
137  /// Refresh - from original
138  void refreshFromPerturbed(double tolerance);
139  /** Sets bounds and cost for one variable
140      Returns change in cost
141      May need to be inline for speed */
142  double setOne(int sequence, double solutionValue);
143  /** Sets bounds and cost for one variable
144      Returns change in cost
145      May need to be inline for speed */
146  double setOneBasic(int iRow, double solutionValue);
147  /** Sets bounds and cost for outgoing variable
148      may change value
149      Returns direction */
150  int setOneOutgoing(int sequence, double &solutionValue);
151  /// Returns nearest bound
152  double nearest(int iRow, double solutionValue);
153  /** Returns change in cost - one down if alpha >0.0, up if <0.0
154      Value is current - new
155  */
156  inline double changeInCost(int /*sequence*/, double alpha) const {
157    return (alpha > 0.0) ? infeasibilityWeight_ : -infeasibilityWeight_;
158  }
159  inline double changeUpInCost(int /*sequence*/) const {
160    return -infeasibilityWeight_;
161  }
162  inline double changeDownInCost(int /*sequence*/) const {
163    return infeasibilityWeight_;
164  }
165  /// This also updates next bound
166  inline double changeInCost(int iRow, double alpha, double &rhs) {
167    int sequence=model_->pivotVariable()[iRow];
168    double returnValue = 0.0;
169    unsigned char iStatus = status_[sequence];
170    int iWhere = currentStatus(iStatus);
171    if (iWhere == CLP_SAME)
172      iWhere = originalStatus(iStatus);
173    // rhs always increases
174    if (iWhere == CLP_FEASIBLE) {
175      if (alpha > 0.0) {
176        // going below
177        iWhere = CLP_BELOW_LOWER;
178        rhs = COIN_DBL_MAX;
179      } else {
180        // going above
181        iWhere = CLP_ABOVE_UPPER;
182        rhs = COIN_DBL_MAX;
183      }
184    } else if (iWhere == CLP_BELOW_LOWER) {
185      assert (alpha < 0);
186      // going feasible
187      iWhere = CLP_FEASIBLE;
188      rhs += bound_[sequence] - model_->upperRegion()[sequence];
189    } else {
190      assert (iWhere == CLP_ABOVE_UPPER);
191      // going feasible
192      iWhere = CLP_FEASIBLE;
193      rhs += model_->lowerRegion()[sequence] - bound_[sequence];
194    }
195    setCurrentStatus(status_[sequence], iWhere);
196    returnValue = fabs(alpha) * infeasibilityWeight_;
197    return returnValue;
198  }
199  //@}
200 
201 
202  /**@name Gets and sets */
203  //@{
204  /// Number of infeasibilities
205  inline int numberInfeasibilities() const {
206    return numberInfeasibilities_;
207  }
208  /// Change in cost
209  inline double changeInCost() const {
210    return changeCost_;
211  }
212  /// Feasible cost
213  inline double feasibleCost() const {
214    return feasibleCost_;
215  }
216  /// Feasible cost with offset and direction (i.e. for reporting)
217  double feasibleReportCost() const;
218  /// Sum of infeasibilities
219  inline double sumInfeasibilities() const {
220    return sumInfeasibilities_;
221  }
222  /// Largest infeasibility
223  inline double largestInfeasibility() const {
224    return largestInfeasibility_;
225  }
226  /// Average theta
227  inline double averageTheta() const {
228    return averageTheta_;
229  }
230  inline void setAverageTheta(double value) {
231    averageTheta_ = value;
232  }
233  inline void setChangeInCost(double value) {
234    changeCost_ = value;
235  }
236  //@}
237  ///@name Private functions to deal with infeasible regions
238  inline unsigned char * statusArray() const {
239    return status_;
240  }
241  inline int getCurrentStatus(int sequence)
242  {return (status_[sequence] >> 4);}
243  /// For debug
244  void validate();
245  //@}
246 
247private:
248  /**@name Data members */
249  //@{
250  /// Change in cost because of infeasibilities
251  double changeCost_;
252  /// Feasible cost
253  double feasibleCost_;
254  /// Current infeasibility weight
255  double infeasibilityWeight_;
256  /// Largest infeasibility
257  double largestInfeasibility_;
258  /// Sum of infeasibilities
259  double sumInfeasibilities_;
260  /// Average theta - kept here as only for primal
261  double averageTheta_;
262  /// Number of rows (mainly for checking and copy)
263  int numberRows_;
264  /// Number of columns (mainly for checking and copy)
265  int numberColumns_;
266  /// Model
267  AbcSimplex * model_;
268  /// Number of infeasibilities found
269  int numberInfeasibilities_;
270  // new stuff
271  /// Contains status at beginning and current
272  unsigned char * status_;
273  /// Bound which has been replaced in lower_ or upper_
274  double * bound_;
275  /// Feasible cost array
276  double * cost_;
277  //@}
278};
279
280#endif
Note: See TracBrowser for help on using the repository browser.