algebra/triangular_matrix.h

// -*- c++ -*-

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

#ifndef _MATHTL_TRIANGULAR_MATRIX_H
#define _MATHTL_TRIANGULAR_MATRIX_H

#include <iostream>
#include <algebra/vector.h>

// matrix norms, for convenience
#include <algebra/matrix_norms.h>

namespace MathTL
{
  template <class C>
  class LowerTriangularMatrix
  {
  public:
    typedef typename Vector<C>::size_type size_type;

    explicit LowerTriangularMatrix(const size_type n = 0);

    LowerTriangularMatrix(const LowerTriangularMatrix<C>& M);

    LowerTriangularMatrix(const size_type rows, const size_type columns);

    LowerTriangularMatrix(const size_type rows,
                          const size_type columns,
                          const char* str,
                          const bool byrow = true);

    const size_type row_dimension() const;
    
    const size_type column_dimension() const;
    
    const size_type size() const;
    
    void resize(const size_type rows, const size_type columns);

    const size_type memory_consumption() const;
    
    bool empty() const;
    
    const C operator () (const size_type row, const size_type column) const;
    
    const C get_entry(const size_type row, const size_type column) const;

    C& operator () (const size_type row, const size_type column);
    
    void set_entry(const size_type row, const size_type column, const C value);

    template <class C2>
    bool operator == (const LowerTriangularMatrix<C2>& M) const;

    template <class C2>
    bool operator != (const LowerTriangularMatrix<C2>& M) const;

    LowerTriangularMatrix<C>& operator = (const LowerTriangularMatrix<C>& M);

    void scale(const C s);

    void inverse(LowerTriangularMatrix<C>& MInv) const;
    
    template <class VECTOR>
    void apply(const VECTOR& x, VECTOR& Mx) const;

    template <class VECTOR>
    void apply_transposed(const VECTOR& x, VECTOR& Mtx) const;

    void print(std::ostream& os,
               const unsigned int tabwidth = 10,
               const unsigned int precision = 3) const;

  protected:
    static size_type triangle_size(const size_type rows, const size_type columns);

    static size_type triangle_index(const size_type row,
                                    const size_type column,
                                    const size_type rowdim,
                                    const size_type coldim);
   
    Vector<C> entries_;

    size_type rowdim_;

    size_type coldim_;
  };

  template <class C>
  LowerTriangularMatrix<C> operator * (const LowerTriangularMatrix<C>& M,
                                       const LowerTriangularMatrix<C>& N);
  
  template <class C>
  std::ostream& operator << (std::ostream& os, const LowerTriangularMatrix<C>& M);


  template <class C>
  class UpperTriangularMatrix
    : public LowerTriangularMatrix<C>
  {
  public:
    typedef typename LowerTriangularMatrix<C>::size_type size_type;

    explicit UpperTriangularMatrix(const size_type n = 0);

    UpperTriangularMatrix(const UpperTriangularMatrix<C>& M);

    UpperTriangularMatrix(const size_type rows, const size_type columns);

    UpperTriangularMatrix(const size_type rows,
                          const size_type columns,
                          const char* str,
                          const bool byrow = true);

    const size_type row_dimension() const;
    
    const size_type column_dimension() const;
    
    void resize(const size_type rows, const size_type columns);

    const C operator () (const size_type row, const size_type column) const;
    
    const C get_entry(const size_type row, const size_type column) const;

    C& operator () (const size_type row, const size_type column);
    
    void set_entry(const size_type row, const size_type column, const C value);

    template <class C2>
    bool operator == (const UpperTriangularMatrix<C2>& M) const;

    template <class C2>
    bool operator != (const UpperTriangularMatrix<C2>& M) const;

    UpperTriangularMatrix<C>& operator = (const UpperTriangularMatrix<C>& M);

    void inverse(UpperTriangularMatrix<C>& MInv) const;
    
    template <class VECTOR>
    void apply(const VECTOR& x, VECTOR& Mx) const;

    template <class VECTOR>
    void apply_transposed(const VECTOR& x, VECTOR& Mtx) const;

    void print(std::ostream& os,
               const unsigned int tabwidth = 10,
               const unsigned int precision = 3) const;
  };

  template <class C>
  std::ostream& operator << (std::ostream& os, const UpperTriangularMatrix<C>& M);
}

// include implementation of inline functions
#include <algebra/triangular_matrix.cpp>

#endif