numerics/matrix_decomp.h

// -*- c++ -*-

// +--------------------------------------------------------------------+
// | This file is part of MathTL - the Mathematical Template Library    |
// |                                                                    |
// | Copyright (c) 2002-2009                                            |
// | Thorsten Raasch, Manuel Werner                                     |
// +--------------------------------------------------------------------+

#ifndef _MATHTL_MATRIX_DECOMP_H
#define _MATHTL_MATRIX_DECOMP_H

#include <algebra/vector.h>
#include <algebra/matrix.h>
#include <algebra/triangular_matrix.h>
#include <algebra/symmetric_matrix.h>

namespace MathTL
{

  template <class C>
  bool CholeskyDecomposition(const SymmetricMatrix<C>& A,
                             LowerTriangularMatrix<C>& L);

  template <class C>
  class QUDecomposition
  {
  public:
    QUDecomposition(const Matrix<C>& A);

    bool hasFullRank() const;

    void getU(UpperTriangularMatrix<C>& U) const;

    void getQ(Matrix<C>& Q) const;

    void solve(const Vector<C>& b, Vector<C>& x) const;

    void inverse(Matrix<C>& AInv) const;

  protected:
    typename Matrix<C>::size_type rowdim_, coldim_;
    Matrix<C> QU_;
    Vector<C> Udiag_;
  };


  template <class C>
  class SVD
  {
  public:
    SVD(const Matrix<C>& A);

    void getU(Matrix<C>& U) const;

    void getUS(Matrix<C>& US) const;

    void getS(Vector<C>& S) const { S = Sdiag_; }

    void getV(Matrix<C>& V) const;

  protected:
    typename Matrix<C>::size_type rowdim_, coldim_;
    Matrix<C> U_, V_;
    Vector<C> Sdiag_;
  };

  template <class C>
  class Hessenberg
    : public Matrix<C>
  {
  public:
    Hessenberg(const Matrix<C>& A);
  };
}

#include <numerics/matrix_decomp.cpp>

#endif