source: trunk/Clp/src/AbcNonLinearCost.hpp

Last change on this file was 2385, checked in by unxusr, 8 months ago

formatting

  • Property svn:keywords set to Id
File size: 8.9 KB
Line 
1/* $Id: AbcNonLinearCost.hpp 2385 2019-01-06 19:43:06Z tkr $ */
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#include "CoinPragma.hpp"
10#include "AbcCommon.hpp"
11
12class AbcSimplex;
13class CoinIndexedVector;
14
15/** Trivial class to deal with non linear costs
16   
17    I don't make any explicit assumptions about convexity but I am
18    sure I do make implicit ones.
19   
20    One interesting idea for normal LP's will be to allow non-basic
21    variables to come into basis as infeasible i.e. if variable at
22    lower bound has very large positive reduced cost (when problem
23    is infeasible) could it reduce overall problem infeasibility more
24    by bringing it into basis below its lower bound.
25   
26    Another feature would be to automatically discover when problems
27    are convex piecewise linear and re-formulate to use non-linear.
28    I did some work on this many years ago on "grade" problems, but
29    while it improved primal interior point algorithms were much better
30    for that particular problem.
31*/
32/* status has original status and current status
33   0 - below lower so stored is upper
34   1 - in range
35   2 - above upper so stored is lower
36   4 - (for current) - same as original
37*/
38#define CLP_BELOW_LOWER 0
39#define CLP_FEASIBLE 1
40#define CLP_ABOVE_UPPER 2
41#define CLP_SAME 4
42#ifndef ClpNonLinearCost_H
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}
70#endif
71class AbcNonLinearCost {
72
73public:
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  /**@name Actual work in primal */
92  //@{
93  /** Changes infeasible costs and computes number and cost of infeas
94      Puts all non-basic (non free) variables to bounds
95      and all free variables to zero if oldTolerance is non-zero
96      - but does not move those <= oldTolerance away*/
97  void checkInfeasibilities(double oldTolerance = 0.0);
98  /** Changes infeasible costs for each variable
99      The indices are row indices and need converting to sequences
100  */
101  void checkInfeasibilities(int numberInArray, const int *index);
102  /** Puts back correct infeasible costs for each variable
103      The input indices are row indices and need converting to sequences
104      for costs.
105      On input array is empty (but indices exist).  On exit just
106      changed costs will be stored as normal CoinIndexedVector
107  */
108  void checkChanged(int numberInArray, CoinIndexedVector *update);
109  /** Goes through one bound for each variable.
110      If multiplier*work[iRow]>0 goes down, otherwise up.
111      The indices are row indices and need converting to sequences
112      Temporary offsets may be set
113      Rhs entries are increased
114  */
115  void goThru(int numberInArray, double multiplier,
116    const int *index, const double *work,
117    double *rhs);
118  /** Takes off last iteration (i.e. offsets closer to 0)
119   */
120  void goBack(int numberInArray, const int *index,
121    double *rhs);
122  /** Puts back correct infeasible costs for each variable
123      The input indices are row indices and need converting to sequences
124      for costs.
125      At the end of this all temporary offsets are zero
126  */
127  void goBackAll(const CoinIndexedVector *update);
128  /// Temporary zeroing of feasible costs
129  void zapCosts();
130  /// Refreshes costs always makes row costs zero
131  void refreshCosts(const double *columnCosts);
132  /// Puts feasible bounds into lower and upper
133  void feasibleBounds();
134  /// Refresh - assuming regions OK
135  void refresh();
136  /// Refresh - from original
137  void refreshFromPerturbed(double tolerance);
138  /** Sets bounds and cost for one variable
139      Returns change in cost
140      May need to be inline for speed */
141  double setOne(int sequence, double solutionValue);
142  /** Sets bounds and cost for one variable
143      Returns change in cost
144      May need to be inline for speed */
145  double setOneBasic(int iRow, double solutionValue);
146  /** Sets bounds and cost for outgoing variable
147      may change value
148      Returns direction */
149  int setOneOutgoing(int sequence, double &solutionValue);
150  /// Returns nearest bound
151  double nearest(int iRow, double solutionValue);
152  /** Returns change in cost - one down if alpha >0.0, up if <0.0
153      Value is current - new
154  */
155  inline double changeInCost(int /*sequence*/, double alpha) const
156  {
157    return (alpha > 0.0) ? infeasibilityWeight_ : -infeasibilityWeight_;
158  }
159  inline double changeUpInCost(int /*sequence*/) const
160  {
161    return -infeasibilityWeight_;
162  }
163  inline double changeDownInCost(int /*sequence*/) const
164  {
165    return infeasibilityWeight_;
166  }
167  /// This also updates next bound
168  inline double changeInCost(int iRow, double alpha, double &rhs)
169  {
170    int sequence = model_->pivotVariable()[iRow];
171    double returnValue = 0.0;
172    unsigned char iStatus = status_[sequence];
173    int iWhere = currentStatus(iStatus);
174    if (iWhere == CLP_SAME)
175      iWhere = originalStatus(iStatus);
176    // rhs always increases
177    if (iWhere == CLP_FEASIBLE) {
178      if (alpha > 0.0) {
179        // going below
180        iWhere = CLP_BELOW_LOWER;
181        rhs = COIN_DBL_MAX;
182      } else {
183        // going above
184        iWhere = CLP_ABOVE_UPPER;
185        rhs = COIN_DBL_MAX;
186      }
187    } else if (iWhere == CLP_BELOW_LOWER) {
188      assert(alpha < 0);
189      // going feasible
190      iWhere = CLP_FEASIBLE;
191      rhs += bound_[sequence] - model_->upperRegion()[sequence];
192    } else {
193      assert(iWhere == CLP_ABOVE_UPPER);
194      // going feasible
195      iWhere = CLP_FEASIBLE;
196      rhs += model_->lowerRegion()[sequence] - bound_[sequence];
197    }
198    setCurrentStatus(status_[sequence], iWhere);
199    returnValue = fabs(alpha) * infeasibilityWeight_;
200    return returnValue;
201  }
202  //@}
203
204  /**@name Gets and sets */
205  //@{
206  /// Number of infeasibilities
207  inline int numberInfeasibilities() const
208  {
209    return numberInfeasibilities_;
210  }
211  /// Change in cost
212  inline double changeInCost() const
213  {
214    return changeCost_;
215  }
216  /// Feasible cost
217  inline double feasibleCost() const
218  {
219    return feasibleCost_;
220  }
221  /// Feasible cost with offset and direction (i.e. for reporting)
222  double feasibleReportCost() const;
223  /// Sum of infeasibilities
224  inline double sumInfeasibilities() const
225  {
226    return sumInfeasibilities_;
227  }
228  /// Largest infeasibility
229  inline double largestInfeasibility() const
230  {
231    return largestInfeasibility_;
232  }
233  /// Average theta
234  inline double averageTheta() const
235  {
236    return averageTheta_;
237  }
238  inline void setAverageTheta(double value)
239  {
240    averageTheta_ = value;
241  }
242  inline void setChangeInCost(double value)
243  {
244    changeCost_ = value;
245  }
246  //@}
247  ///@name Private functions to deal with infeasible regions
248  inline unsigned char *statusArray() const
249  {
250    return status_;
251  }
252  inline int getCurrentStatus(int sequence)
253  {
254    return (status_[sequence] >> 4);
255  }
256  /// For debug
257  void validate();
258  //@}
259
260private:
261  /**@name Data members */
262  //@{
263  /// Change in cost because of infeasibilities
264  double changeCost_;
265  /// Feasible cost
266  double feasibleCost_;
267  /// Current infeasibility weight
268  double infeasibilityWeight_;
269  /// Largest infeasibility
270  double largestInfeasibility_;
271  /// Sum of infeasibilities
272  double sumInfeasibilities_;
273  /// Average theta - kept here as only for primal
274  double averageTheta_;
275  /// Number of rows (mainly for checking and copy)
276  int numberRows_;
277  /// Number of columns (mainly for checking and copy)
278  int numberColumns_;
279  /// Model
280  AbcSimplex *model_;
281  /// Number of infeasibilities found
282  int numberInfeasibilities_;
283  // new stuff
284  /// Contains status at beginning and current
285  unsigned char *status_;
286  /// Bound which has been replaced in lower_ or upper_
287  double *bound_;
288  /// Feasible cost array
289  double *cost_;
290  //@}
291};
292
293#endif
294
295/* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
296*/
Note: See TracBrowser for help on using the repository browser.