matrice.hpp

/*
    Copyright 2017 Philippe Grandclement
 
    This file is part of Kadath.
 
    Kadath 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 3 of the License, or
    (at your option) any later version.
 
    Kadath 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 Kadath.  If not, see <http://www.gnu.org/licenses/>.
*/
 
#ifndef __MATRICE_HPP_
#define __MATRICE_HPP_
//fichiers includes
#include "headcpp.hpp"
#include "array.hpp"
 
namespace Kadath {
    class Matrice : public Memory_mapped {
        //Elements
    private:
 
        Dim_array sizes ;
        Array<double>* std ; 
 
 
        mutable int ku ;    
        mutable int kl ;    
 
        mutable Array<double>* band ;
 
 
        mutable Array<double>* lu ;   
        mutable Array<int>* permute ;   
 
        // Constructeurs destructeurs
    public:
        Matrice (int size1, int size2 ) ;
 
        Matrice (const Matrice& ) ; 
 
        Matrice (const Array<double>& tab) ;
 
        ~Matrice() ; 
 
            //Gestion memoire
    private:
        void del_deriv() ; 
 
    public:
         void annule() ;
 
    public:
        int get_dim(int i) const {return sizes(i) ;};
        Array<double>* get_lu() {return lu ;} ; 
        Array<int>* get_permute() {return permute ;} ;
        const Array<double>& get_array() const {return *std; } ; 
        Array<double>& get_array() {return *std;}
 
        // affectation
    public:
        void operator=(double x) ;
 
        void operator=(const Matrice& ) ; 
        void operator=(const Array<double>& ) ; 
 
        friend ostream& operator<<(ostream& , const Matrice& ) ;
 
        // extraction d'un element :
    public:
        double& set(int i, int j) ;
 
        void copy_inside (int i, int j, const Matrice& so) ;
 
        double operator()(int i, int j) const ;
 
        // Passage matrice a bande
        void set_band (int up, int low) const ;
 
        // Decomposition LU
        void set_lu () const ;
 
        // Inversion de la matrice
        Array<double> solve (const Array<double>& sec_membre) const ;
 
        // Les valeurs propres :
        Array<double> val_propre() const ;
        Matrice vect_propre() const ;
 
        Matrice transpose() const ;
 
 
        // Member arithmetics
        // ------------------
        public:
        void operator+=(const Matrice &) ; 
        void operator+=(double) ; 
        void operator-=(const Matrice &) ; 
        void operator-=(double) ; 
        void operator*=(double) ; 
        void operator/=(double) ; 
 
        // Operateurs amis
        friend Matrice operator+ (const Matrice&, const Matrice& ) ; 
        friend Matrice operator- (const Matrice&, const Matrice& ) ; 
        friend Matrice operator- (const Matrice& ) ; 
        friend Matrice operator* (const Matrice&, double ) ; 
        friend Matrice operator* (double, const Matrice& ) ; 
        friend Matrice operator* (const Matrice&, const Matrice& ) ; 
        friend Matrice operator/ (const Matrice&,  double ) ; 
    } ;
    ostream& operator<<(ostream& , const Matrice& ) ;
 
 
    Matrice operator+ (const Matrice&, const Matrice& ) ;
    Matrice operator- (const Matrice&, const Matrice& ) ;
    Matrice operator- (const Matrice& ) ;
    Matrice operator* (const Matrice&, double ) ;
    Matrice operator* (double, const Matrice& ) ;
    Matrice operator* (const Matrice&, const Matrice& ) ;
    Matrice operator/ (const Matrice&,  double ) ;
 
    // Lapack routines :
    #define F77_dswap dswap_
    #define F77_dgbtrf dgbtrf_
    #define F77_dgbtrs dgbtrs_
    #define F77_dgetrf dgetrf_
    #define F77_dgetrs dgetrs_
    #define F77_dgeev dgeev_
 
    extern "C" {
    void F77_dgbtrf(int*, int*, int*, int*, double[], int*, int[], int *) ;
    void F77_dgbtrs(char*, int*, int*, int*, int*,
                double[], int*, int[], double [], int*, int *) ;
    void F77_dgetrf(int*, int*, double[], int*, int[], int *) ;
    void F77_dgetrs(char*, int*, int*, double[], int*, int[], double [], int*, int* ) ;
    void F77_dgeev(char*, char*, int*, double[], int*, double[], double[],
                double[], int*, double[], int*, double[], int*, int*) ;
    }
 
}
#endif