supermatrix.hpp

#ifndef SUPER_MATRIX_CPP_H
#define SUPER_MATRIX_CPP_H

#include "parallel.hpp"
#include <cstdlib>
#include <exception>
#include <iostream>
#include <fstream>
#include <list>
#include <map>
#include <vector>
#include <memory>

#include "matrix.hpp"

class super_Matrix;

class super_Matrix_iterator_factory{
public:
    virtual size_t nnz() const = 0;
};

class super_Matrix_strategy{
private:
    super_Matrix_strategy(){};
    super_Matrix_strategy(const super_Matrix_strategy&){};
protected:

    super_Matrix* my_sM;

    void fail(const char* msg) const{
    throw ( super_Matrix_Invalid_request(msg));
    }

    super_Matrix_strategy(super_Matrix* m): my_sM(m){};

public:
    typedef STD_TR1::shared_ptr<super_Matrix_iterator_factory> p_iterator_factory;
    typedef STD_TR1::shared_ptr<const super_Matrix_iterator_factory> p_const_iterator_factory;
    
    enum status {SetUp, Baked, Multiply, Factorized, Factorized_dist};

    virtual void addmyselfto(Matrix& m) const = 0;

    virtual ~super_Matrix_strategy(){};

    virtual void obtain_Diagonal(complex *) = 0;

    virtual void Insert(row, column, complex) = 0;
    virtual void Add(row, column, complex) = 0;
    virtual complex Get(row, column) const = 0;

    virtual void Apply(const complex* , 
                     complex* ) const = 0;
    virtual void Apply(complex*) const = 0;

    virtual void Scale(complex) = 0;

    virtual void Output(std::ostream& os) const = 0;

    virtual p_iterator_factory provide_iterator() = 0;
    virtual p_const_iterator_factory provide_const_iterator() const = 0;

    virtual std::unique_ptr<super_Matrix_strategy> deep_copy(
    STD_TR1::shared_ptr<super_Matrix> new_mat) const = 0;

    virtual status my_status() const = 0;
};


class super_Matrix: public calculation_Matrix{
public:
    typedef STD_TR1::shared_ptr<super_Matrix_iterator_factory> p_iterator_factory;
    typedef STD_TR1::shared_ptr<const super_Matrix_iterator_factory> p_const_iterator_factory;

    typedef std::pair<int,int> indexpair;
    
    typedef super_Matrix_strategy::status status;

    typedef STD_TR1::shared_ptr<super_Matrix> ref;
    typedef STD_TR1::shared_ptr<const super_Matrix> const_ref;

    explicit super_Matrix(MPI_Comm comm, 
              const int N);
    virtual ~super_Matrix() throw();

    status get_status(){return strategy->my_status();};

    p_iterator_factory provide_iterator(){ return strategy->provide_iterator();};    
    p_const_iterator_factory provide_const_iterator() const { return strategy->provide_const_iterator();};    
    
    void set_param(const std::string&, double);
    void set_param(const std::string&, int);

    bool get_double_param(const std::string&, double&) const;
    bool get_int_param(const std::string&, int&) const;
    
    void obtain_Diagonal(complex *);

    void Add_To(Matrix& m) const;

    void Output(std::ostream& os) const {strategy->Output(os);};

    void Bake_to_Multiply(bool use_mkl = true);
    void Bake_to_Distributed_Multiply();
    void Bake_to_Factorize(bool distributed = false);

    void Factorize();
    void Reset();

    ref deep_copy() const;

    MPI_Comm get_comm() const { return comm;};


private:
    int NN;

    MPI_Comm comm;

    std::map<std::string, int> iparams;
    std::map<std::string, double> dparams;
    
    std::unique_ptr<super_Matrix_strategy> strategy;

//    super_Matrix(const super_Matrix&);

    void scale(complex v) { strategy->Scale(v);};
    void insert(row r, column c, complex v){ strategy->Insert(r,c,v);};
    void add(row r, column c, complex v){ strategy->Add(r,c,v);};
    complex get(row r, column c) const { return strategy->Get(r,c);};

    void apply_two(const complex* phi, 
           complex* result) const {strategy->Apply(phi, result);};
    void apply_one(complex* phi) const {strategy->Apply(phi);};

    size_t num_rows() const {return NN;};
    size_t num_cols() const {return NN;};
};

/*-----------------*
 * SETUP STRATEGY  *
 *-----------------*/

class super_Matrix_SetUp_iterator: public super_Matrix_iterator_factory{    
public:
    typedef std::pair<int,int> indexpair;
private:
    const std::map<indexpair, complex>& data;
public:
    super_Matrix_SetUp_iterator(const std::map<indexpair, complex>& d):
    data(d){}

    size_t nnz() const{ return data.size();};
    
    std::map<std::pair<int,int>,complex>::const_iterator begin() const {return data.begin();};
    std::map<std::pair<int,int>,complex>::const_iterator end() const {return data.end();};
};

class super_Matrix_SetUp: public super_Matrix_strategy{


    typedef std::pair<int,int> indexpair;
    std::map<indexpair, complex> data;
    int NN;

    indexpair index(row r, column c) const { return std::pair<int,int>(r.get(),c.get());};

    int get_row(const indexpair& i) const {return i.first;};
    int get_col(const indexpair& i) const {return i.second;};

    friend class matrix_factory;

    void Scale(complex v);


public:

    p_iterator_factory provide_iterator();
    p_const_iterator_factory provide_const_iterator() const;
    
    super_Matrix_SetUp(super_Matrix* m): 
    super_Matrix_strategy(m), 
    NN(m->num_of_cols()){};

    ~super_Matrix_SetUp(){};

    void addmyselfto(Matrix&) const;
  
    void Insert(row, column, complex);
    void Add(row, column, complex);
    complex Get(row, column) const;

    void obtain_Diagonal(complex* diag);

    void Apply(const complex*, 
                     complex*) const{
    fail("Apply not valid in SetUp Phase!\n");
    };
    void Apply(      complex*) const{
    fail("Apply not valid in SetUp Phase!\n");
    };

    std::unique_ptr<super_Matrix_strategy> deep_copy(
    STD_TR1::shared_ptr<super_Matrix> new_mat
    ) const;
    
    void Output(std::ostream& os) const;

    status my_status() const {return SetUp;};
};
#endif