/**
 * \file Mmn.h
 *
 * \brief Matrizenmanipulation beliebiger Matrizen
 *
 * Stellt eine Klasse zur verfügung, mit der sich alle denkbaren
 * algebraischen Manipulationen an Matrizen beliebiger größe
 * ausführen lassen, incl. lösen eines Gleichungssystems dessen
 * erweiterte Koeffizientenmatrix diese Matrix ist.
 *
 * \author Georg Steffers <georg@steffers.org> [GS]
 *
 * \date 06.12.2003
 *
 * \version ..2002 erste funktionierende Version [GS]
 * \version 06.12.2003 beginn der Dokumentation via doxygen [GS]
 */

/* 
 * Copyright (C)2003 Georg Steffers
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 * 
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 */

#ifndef __MMN_H__
#define __MMN_H__

using namespace std;

#include <cmath>
#include <iostream>

#include "sin_cos.h"

typedef enum gauss_typ {TREPPENNORMALE, INVERSE}; //!< \brief Was soll mit dem
                                                  //! Gauss-Algorithmus
                                                  //! ermittelt werden?

template <class K>
class Mmn;

//! Skalarmultiplikation Skalar * Matrize
template <class K>
Mmn<K> operator*(const K&, const Mmn<K>&);

/**
 * \brief Darstellung einer Matrize und der auf ihr m&ouml;glichen
 *        algebraischen Operationen
 *
 * Diese Klasse stellt meinen Lernerfolg in Lineare Algebra I dar.
 */
template <class K>
class Mmn {
    protected:
        K** _a;         //!< \brief Die Matrize als doppelter Pointer auf
                        //! den skalaren Typ K
        unsigned m, n;  //!< Zeilen, Spalten der Matrize
        unsigned Rg;    //!< \brief Der Rang der Matrize
                        //!
                        //! wird nur gefüllt wenn die Matrize &uuml;ber 
                        //! Gauss in Treppennormalform gebracht wurde.
        bool Det_F;     //!< enth&auml;lt Det einen g&uuml;ltigen Wert?
        K Det;          //!< \brief Die Determinante der Matrize
                        //!
                        //! Sobald sie errechnet wurde wird sie hier
                        //! gespeichert, bis die Matrize ver&auml;ndert
                        //! wird, damit sie nicht immer wieder neu berechnet
                        //! werden mu&szlig;

        void freeall(void);             //!< die matrize wieder freigeben
        K sarrus(void) const;  //!< \brief Sarrusregel
                               //!
                               //! zum ermitteln
                               //! der Determinante einer 3x3 Matrize
        K laplace(void) const; //!< \brief Laplace'scher Entwicklungssatz
                               //!
                               //! zum ermitteln der Determinante bei
                               //! gr&ouml;&szlig;eren Matrizen als 3x3

    public:
        // Konstruktoren
        //! Defaultkontruktor
        Mmn() : _a(NULL), m(0), n(0), Rg(0), Det_F(false) {}
        Mmn(unsigned);                   //!< EinheitsMatrix
        Mmn(unsigned, unsigned);         //!< NullMatrix
        Mmn(K*, unsigned, unsigned);     //!< \brief Konstruktor für eine 
                                         //!         beliebige Matrix
        Mmn(const Mmn<K>&);              //!< Copykonstruktor
        virtual ~Mmn();                  //!< Destruktor

        //! Zuweisungsoperator
        const Mmn<K>& operator=(const Mmn<K>&);

        //! Skalarmultiplikation (Matrix * Skalar)
        Mmn<K> operator*(const K&) const;
        //! Als friend damit auch Parameter 1 ein Skalar sein kann
        friend Mmn<K> operator*<>(const K&, const Mmn<K>&);

        // Hier können eh nur Matrizen verwendet werden. Daher sind
        // keine weiteren friends nötig
        Mmn<K> operator%(const Mmn<K>&) const;   //!< Matrizenprodukt
        Mmn<K> operator+(const Mmn<K>&) const;   //!< Matrizensumme
        Mmn<K> operator-(void) const;            //!< Negation
        Mmn<K> operator-(const Mmn<K>&) const;   //!< Subtraktion
       
        //! eine NullMatrix derselben Groesse
        Mmn<K> null(void) const { return Mmn<K>(m,n); }

        bool is_null(void) const {
           for(unsigned i=0; i<m; i++)
              for(unsigned j=0; j<n; j++)
                 if(_a[i][j]!=(K)0)
                    return false;

           return true;
        }

        // Einheitsmatrizen, entweder zu Anzahl Spalten, oder Zeilen
        //! Einheitsmatrizen zu Anzahl Spalten
        Mmn<K> In(void) const { return Mmn<K>(n); }
        //! Einheitsmatrizen zu Anzahl Zeilen
        Mmn<K> Im(void) const { return Mmn<K>(m); }

        //! Ein Element aij der Matrize
        K& a(unsigned i, unsigned j) { Det_F=false; return _a[i][j]; } 
        const K& a(unsigned i, unsigned j) const { return _a[i][j]; } 

        //! Determinante
        K det(void);
        K det(void) const;
        
        //! Transponierte
        Mmn<K> T(void) const;

        // Zeile und/oder Spalte aus matrize enftfernen
        //! Zeile aus matrize enftfernen
        Mmn<K> Ai(unsigned) const;
        //! Spalte aus matrize enftfernen
        Mmn<K> Aj(unsigned) const;
        //! Zeile und Spalte aus matrize enftfernen
        Mmn<K> Aij(unsigned, unsigned) const;

        // Adjunkten ermitteln (2 verschiedene Verfahren)
        Mmn<K> Ad_det(void) const;  //!< Adjunkte ueber Determinanten entwickelt
        Mmn<K> Ad_inv(void) const;  //!< Inverse=>Adjunkte / gauss=>Adj.-Satz

        // Inverse ermitteln (2 verschiedene Verfahren)
        Mmn<K> gauss(gauss_typ) const; //!< \brief inverse Matrix ueber 
                                       //!         Gauss-Algorithmus
        Mmn<K> inv_Ad(void) const;     //!< \brief Adjunkte=>Inverse / 
                                       //!         Ad_det=>Adj.-Satz

        // lineare Gleichungssysteme
        Mmn<K> solve(void) const;         //!< Gleichungssystem loesen
        Mmn<K> special_solve(void) const; //!< \brief eine Spezielle Loesung des
                                          //!         mit an=1, damit alle xi!=0

        // formatierte Ausgabe
        void print(ostream& = cout) const;    //!< Daten der Matrize ausgeben
};

template <class K>
Mmn<K> operator*(const K& a, const Mmn<K>& A) {
    return A * a;
}

template <class K>
inline ostream& operator<<(ostream& strm, const Mmn<K>& A) {
    A.print(strm);

    return strm;
}

template <class K>
inline void Mmn<K>::print(ostream& strm) const {
    strm << "lines: " << m << "\n";
    strm << "cols: " << n << "\n";

    for(unsigned i=0; i<m; i++) {
        for(unsigned j=0; j<n; j++)
            strm << _a[i][j] << "\t";
        strm << "\n";
    }
}

#include "Mmn_tpl.h"
#endif // __MMN_H__