source: trunk/Clp/src/ClpCholeskyBase.hpp @ 1367

Last change on this file since 1367 was 1367, checked in by forrest, 11 years ago

try and improve stability - also long double interior point

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 7.8 KB
Line 
1// Copyright (C) 2003, International Business Machines
2// Corporation and others.  All Rights Reserved.
3#ifndef ClpCholeskyBase_H
4#define ClpCholeskyBase_H
5
6#include "CoinPragma.hpp"
7#include "CoinFinite.hpp"
8//#define CLP_LONG_CHOLESKY 0
9#ifndef CLP_LONG_CHOLESKY
10#define CLP_LONG_CHOLESKY 0
11#endif
12/* valid combinations are
13   CLP_LONG_CHOLESKY 0 and COIN_LONG_WORK 0
14   CLP_LONG_CHOLESKY 1 and COIN_LONG_WORK 1
15   CLP_LONG_CHOLESKY 2 and COIN_LONG_WORK 1
16*/
17#if COIN_LONG_WORK==0
18#if CLP_LONG_CHOLESKY>0
19#define CHOLESKY_BAD_COMBINATION
20#endif
21#else
22#if CLP_LONG_CHOLESKY==0
23#define CHOLESKY_BAD_COMBINATION
24#endif
25#endif
26#ifdef CHOLESKY_BAD_COMBINATION
27#  warning("Bad combination of CLP_LONG_CHOLESKY and COIN_BIG_DOUBLE/COIN_LONG_WORK");
28"Bad combination of CLP_LONG_CHOLESKY and COIN_LONG_WORK"
29#endif
30#if CLP_LONG_CHOLESKY>1
31typedef long double longDouble;
32#define CHOL_SMALL_VALUE 1.0e-15
33#elif CLP_LONG_CHOLESKY==1
34typedef double longDouble;
35#define CHOL_SMALL_VALUE 1.0e-11
36#else
37typedef double longDouble;
38#define CHOL_SMALL_VALUE 1.0e-11
39#endif
40class ClpInterior;
41class ClpCholeskyDense;
42class ClpMatrixBase;
43
44/** Base class for Clp Cholesky factorization
45    Will do better factorization.  very crude ordering
46
47    Derived classes may be using more sophisticated methods
48*/
49
50class ClpCholeskyBase  {
51 
52public:
53   /**@name Virtual methods that the derived classes may provide  */
54   //@{
55  /** Orders rows and saves pointer to matrix.and model.
56   returns non-zero if not enough memory.
57   You can use preOrder to set up ADAT
58   If using default symbolic etc then must set sizeFactor_ to
59   size of input matrix to order (and to symbolic).
60   Also just permute_ and permuteInverse_ should be created */
61  virtual int order(ClpInterior * model);
62  /** Does Symbolic factorization given permutation.
63      This is called immediately after order.  If user provides this then
64      user must provide factorize and solve.  Otherwise the default factorization is used
65      returns non-zero if not enough memory */
66  virtual int symbolic();
67  /** Factorize - filling in rowsDropped and returning number dropped.
68      If return code negative then out of memory */
69  virtual int factorize(const CoinWorkDouble * diagonal, int * rowsDropped) ;
70  /** Uses factorization to solve. */
71  virtual void solve (double * region) ;
72  /** Uses factorization to solve. - given as if KKT.
73   region1 is rows+columns, region2 is rows */
74  virtual void solveKKT (CoinWorkDouble * region1, CoinWorkDouble * region2, const CoinWorkDouble * diagonal,
75                         CoinWorkDouble diagonalScaleFactor);
76#if CLP_LONG_CHOLESKY
77  /** Uses factorization to solve. */
78  virtual void solve (CoinWorkDouble * region) ;
79#endif
80  //@}
81
82  /**@name Gets */
83  //@{
84  /// status.  Returns status
85  inline int status() const 
86  {return status_;}
87  /// numberRowsDropped.  Number of rows gone
88  inline int numberRowsDropped() const 
89  {return numberRowsDropped_;}
90  /// reset numberRowsDropped and rowsDropped.
91  void resetRowsDropped();
92  /// rowsDropped - which rows are gone
93  inline char * rowsDropped() const 
94  {return rowsDropped_;}
95  /// choleskyCondition.
96  inline double choleskyCondition() const 
97  {return choleskyCondition_;}
98  /// goDense i.e. use dense factoriaztion if > this (default 0.7).
99  inline double goDense() const 
100  {return goDense_;}
101  /// goDense i.e. use dense factoriaztion if > this (default 0.7).
102  inline void setGoDense(double value)
103  {goDense_=value;}
104  /// rank.  Returns rank
105  inline int rank() const 
106  {return numberRows_-numberRowsDropped_;}
107  /// Return number of rows
108  inline int numberRows() const 
109  {return numberRows_;}
110  /// Return size
111  inline CoinBigIndex size() const
112  { return sizeFactor_;}
113  /// Return sparseFactor
114  inline longDouble * sparseFactor() const
115  { return sparseFactor_;}
116  /// Return diagonal
117  inline longDouble * diagonal() const
118  { return diagonal_;}
119  /// Return workDouble
120  inline longDouble * workDouble() const
121  { return workDouble_;}
122  /// If KKT on
123  inline bool kkt() const
124  { return doKKT_;}
125  /// Set KKT
126  inline void setKKT(bool yesNo)
127  { doKKT_ = yesNo;}
128  /// Set integer parameter
129  inline void setIntegerParameter(int i,int value)
130  { integerParameters_[i]=value;}
131  /// get integer parameter
132  inline int getIntegerParameter(int i)
133  { return integerParameters_[i];}
134  /// Set double parameter
135  inline void setDoubleParameter(int i,double value)
136  { doubleParameters_[i]=value;}
137  /// get double parameter
138  inline double getDoubleParameter(int i)
139  { return doubleParameters_[i];}
140   //@}
141 
142 
143public:
144
145   /**@name Constructors, destructor
146    */
147   //@{
148  /** Constructor which has dense columns activated.
149      Default is off. */
150  ClpCholeskyBase(int denseThreshold=-1);
151   /** Destructor (has to be public) */
152   virtual ~ClpCholeskyBase();
153  /// Copy
154   ClpCholeskyBase(const ClpCholeskyBase&);
155  /// Assignment
156   ClpCholeskyBase& operator=(const ClpCholeskyBase&);
157   //@}
158  //@{
159  ///@name Other
160  /// Clone
161  virtual ClpCholeskyBase * clone() const;
162 
163  /// Returns type
164  inline int type() const
165  { if (doKKT_) return 100; else return type_;}
166protected:
167  /// Sets type
168  void setType(int type) {type_=type;}
169   //@}
170   
171  /**@name Symbolic, factor and solve */
172  //@{
173  /** Symbolic1  - works out size without clever stuff.
174      Uses upper triangular as much easier.
175      Returns size
176   */
177  int symbolic1(const CoinBigIndex * Astart, const int * Arow);
178  /** Symbolic2  - Fills in indices
179      Uses lower triangular so can do cliques etc
180   */
181  void symbolic2(const CoinBigIndex * Astart, const int * Arow);
182  /** Factorize - filling in rowsDropped and returning number dropped
183      in integerParam.
184   */
185  void factorizePart2(int * rowsDropped) ;
186  /** solve - 1 just first half, 2 just second half - 3 both.
187  If 1 and 2 then diagonal has sqrt of inverse otherwise inverse
188  */
189  void solve(double * region, int type);
190  /// Forms ADAT - returns nonzero if not enough memory
191  int preOrder(bool lowerTriangular, bool includeDiagonal, bool doKKT);
192  /// Updates dense part (broken out for profiling)
193  void updateDense(longDouble * d, longDouble * work, int * first);
194  //@}
195   
196protected:
197  /**@name Data members
198     The data members are protected to allow access for derived classes. */
199  //@{
200  /// type (may be useful) if > 20 do KKT
201  int type_;
202  /// Doing full KKT (only used if default symbolic and factorization)
203  bool doKKT_;
204  /// Go dense at this fraction
205  double goDense_;
206  /// choleskyCondition.
207  double choleskyCondition_;
208  /// model.
209  ClpInterior * model_;
210  /// numberTrials.  Number of trials before rejection
211  int numberTrials_;
212  /// numberRows.  Number of Rows in factorization
213  int numberRows_;
214  /// status.  Status of factorization
215  int status_;
216  /// rowsDropped
217  char * rowsDropped_;
218  /// permute inverse.
219  int * permuteInverse_;
220  /// main permute.
221  int * permute_;
222  /// numberRowsDropped.  Number of rows gone
223  int numberRowsDropped_;
224  /// sparseFactor.
225  longDouble * sparseFactor_;
226  /// choleskyStart - element starts
227  CoinBigIndex * choleskyStart_;
228  /// choleskyRow (can be shorter than sparsefactor)
229  int * choleskyRow_;
230  /// Index starts
231  CoinBigIndex * indexStart_;
232  /// Diagonal
233  longDouble * diagonal_;
234  /// double work array
235  longDouble * workDouble_;
236  /// link array
237  int * link_;
238  // Integer work array
239  CoinBigIndex * workInteger_;
240  // Clique information
241  int * clique_;
242  /// sizeFactor.
243  CoinBigIndex sizeFactor_;
244  /// Size of index array
245  CoinBigIndex sizeIndex_;
246  /// First dense row
247  int firstDense_;
248  /// integerParameters
249  int integerParameters_[64];
250  /// doubleParameters;
251  double doubleParameters_[64];
252  /// Row copy of matrix
253  ClpMatrixBase * rowCopy_;
254  /// Dense indicators
255  char * whichDense_;
256  /// Dense columns (updated)
257  longDouble * denseColumn_;
258  /// Dense cholesky
259  ClpCholeskyDense * dense_;
260  /// Dense threshold (for taking out of Cholesky)
261  int denseThreshold_;
262  //@}
263};
264
265#endif
Note: See TracBrowser for help on using the repository browser.