summaryrefslogtreecommitdiff
path: root/thirdparty/linux/include/coin/coin/ClpCholeskyBase.hpp
blob: 815af014ce328f08f0a66d7fde18a09ee0482dbc (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
/* $Id: ClpCholeskyBase.hpp 1722 2011-04-17 09:58:37Z stefan $ */
// Copyright (C) 2003, International Business Machines
// Corporation and others.  All Rights Reserved.
// This code is licensed under the terms of the Eclipse Public License (EPL).

#ifndef ClpCholeskyBase_H
#define ClpCholeskyBase_H

#include "CoinPragma.hpp"
#include "CoinTypes.hpp"
//#define CLP_LONG_CHOLESKY 0
#ifndef CLP_LONG_CHOLESKY
#define CLP_LONG_CHOLESKY 0
#endif
/* valid combinations are
   CLP_LONG_CHOLESKY 0 and COIN_LONG_WORK 0
   CLP_LONG_CHOLESKY 1 and COIN_LONG_WORK 1
   CLP_LONG_CHOLESKY 2 and COIN_LONG_WORK 1
*/
#if COIN_LONG_WORK==0
#if CLP_LONG_CHOLESKY>0
#define CHOLESKY_BAD_COMBINATION
#endif
#else
#if CLP_LONG_CHOLESKY==0
#define CHOLESKY_BAD_COMBINATION
#endif
#endif
#ifdef CHOLESKY_BAD_COMBINATION
#  warning("Bad combination of CLP_LONG_CHOLESKY and COIN_BIG_DOUBLE/COIN_LONG_WORK");
"Bad combination of CLP_LONG_CHOLESKY and COIN_LONG_WORK"
#endif
#if CLP_LONG_CHOLESKY>1
typedef long double longDouble;
#define CHOL_SMALL_VALUE 1.0e-15
#elif CLP_LONG_CHOLESKY==1
typedef double longDouble;
#define CHOL_SMALL_VALUE 1.0e-11
#else
typedef double longDouble;
#define CHOL_SMALL_VALUE 1.0e-11
#endif
class ClpInterior;
class ClpCholeskyDense;
class ClpMatrixBase;

/** Base class for Clp Cholesky factorization
    Will do better factorization.  very crude ordering

    Derived classes may be using more sophisticated methods
*/

class ClpCholeskyBase  {

public:
     /**@name Virtual methods that the derived classes may provide  */
     //@{
     /** Orders rows and saves pointer to matrix.and model.
      returns non-zero if not enough memory.
      You can use preOrder to set up ADAT
      If using default symbolic etc then must set sizeFactor_ to
      size of input matrix to order (and to symbolic).
      Also just permute_ and permuteInverse_ should be created */
     virtual int order(ClpInterior * model);
     /** Does Symbolic factorization given permutation.
         This is called immediately after order.  If user provides this then
         user must provide factorize and solve.  Otherwise the default factorization is used
         returns non-zero if not enough memory */
     virtual int symbolic();
     /** Factorize - filling in rowsDropped and returning number dropped.
         If return code negative then out of memory */
     virtual int factorize(const CoinWorkDouble * diagonal, int * rowsDropped) ;
     /** Uses factorization to solve. */
     virtual void solve (CoinWorkDouble * region) ;
     /** Uses factorization to solve. - given as if KKT.
      region1 is rows+columns, region2 is rows */
     virtual void solveKKT (CoinWorkDouble * region1, CoinWorkDouble * region2, const CoinWorkDouble * diagonal,
                            CoinWorkDouble diagonalScaleFactor);
private:
     /// AMD ordering
     int orderAMD();
public:
     //@}

     /**@name Gets */
     //@{
     /// status.  Returns status
     inline int status() const {
          return status_;
     }
     /// numberRowsDropped.  Number of rows gone
     inline int numberRowsDropped() const {
          return numberRowsDropped_;
     }
     /// reset numberRowsDropped and rowsDropped.
     void resetRowsDropped();
     /// rowsDropped - which rows are gone
     inline char * rowsDropped() const {
          return rowsDropped_;
     }
     /// choleskyCondition.
     inline double choleskyCondition() const {
          return choleskyCondition_;
     }
     /// goDense i.e. use dense factoriaztion if > this (default 0.7).
     inline double goDense() const {
          return goDense_;
     }
     /// goDense i.e. use dense factoriaztion if > this (default 0.7).
     inline void setGoDense(double value) {
          goDense_ = value;
     }
     /// rank.  Returns rank
     inline int rank() const {
          return numberRows_ - numberRowsDropped_;
     }
     /// Return number of rows
     inline int numberRows() const {
          return numberRows_;
     }
     /// Return size
     inline CoinBigIndex size() const {
          return sizeFactor_;
     }
     /// Return sparseFactor
     inline longDouble * sparseFactor() const {
          return sparseFactor_;
     }
     /// Return diagonal
     inline longDouble * diagonal() const {
          return diagonal_;
     }
     /// Return workDouble
     inline longDouble * workDouble() const {
          return workDouble_;
     }
     /// If KKT on
     inline bool kkt() const {
          return doKKT_;
     }
     /// Set KKT
     inline void setKKT(bool yesNo) {
          doKKT_ = yesNo;
     }
     /// Set integer parameter
     inline void setIntegerParameter(int i, int value) {
          integerParameters_[i] = value;
     }
     /// get integer parameter
     inline int getIntegerParameter(int i) {
          return integerParameters_[i];
     }
     /// Set double parameter
     inline void setDoubleParameter(int i, double value) {
          doubleParameters_[i] = value;
     }
     /// get double parameter
     inline double getDoubleParameter(int i) {
          return doubleParameters_[i];
     }
     //@}


public:

     /**@name Constructors, destructor
      */
     //@{
     /** Constructor which has dense columns activated.
         Default is off. */
     ClpCholeskyBase(int denseThreshold = -1);
     /** Destructor (has to be public) */
     virtual ~ClpCholeskyBase();
     /// Copy
     ClpCholeskyBase(const ClpCholeskyBase&);
     /// Assignment
     ClpCholeskyBase& operator=(const ClpCholeskyBase&);
     //@}
     //@{
     ///@name Other
     /// Clone
     virtual ClpCholeskyBase * clone() const;

     /// Returns type
     inline int type() const {
          if (doKKT_) return 100;
          else return type_;
     }
protected:
     /// Sets type
     inline void setType(int type) {
          type_ = type;
     }
     /// model.
     inline void setModel(ClpInterior * model) {
          model_ = model;
     }
     //@}

     /**@name Symbolic, factor and solve */
     //@{
     /** Symbolic1  - works out size without clever stuff.
         Uses upper triangular as much easier.
         Returns size
      */
     int symbolic1(const CoinBigIndex * Astart, const int * Arow);
     /** Symbolic2  - Fills in indices
         Uses lower triangular so can do cliques etc
      */
     void symbolic2(const CoinBigIndex * Astart, const int * Arow);
     /** Factorize - filling in rowsDropped and returning number dropped
         in integerParam.
      */
     void factorizePart2(int * rowsDropped) ;
     /** solve - 1 just first half, 2 just second half - 3 both.
     If 1 and 2 then diagonal has sqrt of inverse otherwise inverse
     */
     void solve(CoinWorkDouble * region, int type);
     /// Forms ADAT - returns nonzero if not enough memory
     int preOrder(bool lowerTriangular, bool includeDiagonal, bool doKKT);
     /// Updates dense part (broken out for profiling)
     void updateDense(longDouble * d, /*longDouble * work,*/ int * first);
     //@}

protected:
     /**@name Data members
        The data members are protected to allow access for derived classes. */
     //@{
     /// type (may be useful) if > 20 do KKT
     int type_;
     /// Doing full KKT (only used if default symbolic and factorization)
     bool doKKT_;
     /// Go dense at this fraction
     double goDense_;
     /// choleskyCondition.
     double choleskyCondition_;
     /// model.
     ClpInterior * model_;
     /// numberTrials.  Number of trials before rejection
     int numberTrials_;
     /// numberRows.  Number of Rows in factorization
     int numberRows_;
     /// status.  Status of factorization
     int status_;
     /// rowsDropped
     char * rowsDropped_;
     /// permute inverse.
     int * permuteInverse_;
     /// main permute.
     int * permute_;
     /// numberRowsDropped.  Number of rows gone
     int numberRowsDropped_;
     /// sparseFactor.
     longDouble * sparseFactor_;
     /// choleskyStart - element starts
     CoinBigIndex * choleskyStart_;
     /// choleskyRow (can be shorter than sparsefactor)
     int * choleskyRow_;
     /// Index starts
     CoinBigIndex * indexStart_;
     /// Diagonal
     longDouble * diagonal_;
     /// double work array
     longDouble * workDouble_;
     /// link array
     int * link_;
     // Integer work array
     CoinBigIndex * workInteger_;
     // Clique information
     int * clique_;
     /// sizeFactor.
     CoinBigIndex sizeFactor_;
     /// Size of index array
     CoinBigIndex sizeIndex_;
     /// First dense row
     int firstDense_;
     /// integerParameters
     int integerParameters_[64];
     /// doubleParameters;
     double doubleParameters_[64];
     /// Row copy of matrix
     ClpMatrixBase * rowCopy_;
     /// Dense indicators
     char * whichDense_;
     /// Dense columns (updated)
     longDouble * denseColumn_;
     /// Dense cholesky
     ClpCholeskyDense * dense_;
     /// Dense threshold (for taking out of Cholesky)
     int denseThreshold_;
     //@}
};

#endif