#ifndef FILE_CHOLESKY
#define FILE_CHOLESKY

/****************************************************************************/
/* File:   cholesky.hpp                                                     */
/* Author: Joachim Schoeberl                                                */
/* Date:   25. Mar. 2000, 16. June 2002                                     */
/****************************************************************************/


/**
   The Cholesky-factorization of a symmetric dense matrix.
   A = L D L^T
 */
template <class T>
class CholeskyFactors
{
protected:
  /// matrix size
  int n;
  /// left factor
  T * lfact;
  /// inverse diagonal
  T * diag;
public:
  typedef typename mat_traits<T>::TV_COL TV;
  /// Factor the matrix A
  CholeskyFactors (const FlatMatrix<T> & a);
  /// Delete memory
  ~CholeskyFactors ();
  /// Multiply with the inverse of A 
  void Mult (const FlatVector<TV> & x, FlatVector<TV> & y) const;
  /// Print factorization
  ostream & Print (ostream & ost) const;
private:
  /// first element in row
  T * PRow (int i) { return lfact + (i*(i-1)) / 2; }
  /// first element in row
  const T * PRow (int i) const { return lfact + (i*(i-1)) / 2; }
};


///
template<typename T>
inline std::ostream & operator<< (std::ostream & s, const CholeskyFactors<T> & m)
{
  m.Print (s);
  return s;
}
#endif