/* $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