#ifndef SMTX_H_SEEN 
#define SMTX_H_SEEN 

#include "machdefs.h"
#include <iostream>

#include "ndatablock.h"
#include "sopnamsp.h"
#include "ctimer.h"
#include "timing.h"

#include "jet.h"


//========================================================
//   Test d'expression template sur tableaux 
//       Reza - Avril 2003
//========================================================


// ---------------------------------------------------
//------- Expression/iterateur sur matrice -----------
// ---------------------------------------------------
template <class T>
class JETExpMtx {
public:
  JETExpMtx() : db_(NULL), msz_(NULL), pd_(NULL) { } 
  //cout << " DBG-AA" << hex << db_ << " - " << msz_ << dec << endl;}
  JETExpMtx(NDataBlock<T>   * db, MtxSize const * msz) : 
    db_(db), msz_(msz), pd_(NULL) { } 
  // cout << " DBG-BBB" << hex << db_ << " - " << msz_ << dec << endl; }
  JETExpMtx(JETExpMtx const & a) : db_(a.db_), msz_(a.msz_), pd_(a.pd_) { } 
  // cout << " DBG-CCC" << hex << db_ << " - " << msz_ << dec << endl;}
  inline T EvaluateAndStep() const { return (*pd_++); }
  inline void Reset() const { pd_ = db_->Begin(); } 
  //{ if (db_ == NULL) cout << " AAA-BUG db_" << endl;
  //pd_ = db_->Begin(); if (pd_ == NULL) cout << " BBB-BUG db_" << endl;} 
  inline MtxSize const & RefSize()  const { return *msz_; } 
  // if (msz_ == NULL) cout << " BUG msz_ " << endl; return *msz_; }
  void JETExpMtxSet(NDataBlock<T> * db, MtxSize const * msz)
  { db_ = db; msz_ = msz; pd_ = NULL; }
protected:
  NDataBlock<T> * db_;
  MtxSize const * msz_;
  mutable T * pd_;
};


// ---------------------------------------------------
// ------------- Classe de matrice -------------------
// ---------------------------------------------------

template <class T>
class SimpleMatrix : public AnyDataObj, public MtxSize, 
		     public JETExpression<T, JETExpMtx<T> > {
public:
		SimpleMatrix();
		SimpleMatrix(sa_size_t nr, sa_size_t nc);
		SimpleMatrix(const SimpleMatrix<T> & m);
		SimpleMatrix(const SimpleMatrix<T> & m, bool share);
  virtual       ~SimpleMatrix();

  virtual void ReSize(const MtxSize & a);
  
  virtual SimpleMatrix<T>& Set(const SimpleMatrix<T> & a);

  inline  SimpleMatrix<T>& Set(JETEvaluator const &  exp)  
  {  exp.EvaluateAll(*this);  return(*this); }

  virtual SimpleMatrix<T>& Init(T c);
  virtual SimpleMatrix<T>& Init(T c, T step);

  inline SimpleMatrix<T>& operator = (const SimpleMatrix<T> & a)
  { return Set(a); }
  inline SimpleMatrix<T>& operator = (JETEvaluator const & a)
  { return Set(a); }
  inline SimpleMatrix<T>& operator = (T c)
  { return Init(c); }
 
  inline T operator()(int r, int c) const 
                { return data_p_[r*ncol_+c]; }
  inline T& operator()(int r, int c) 
                { return data_p_[r*ncol_+c]; }

  inline sa_size_t NRows() const {return nrow_; }
  inline sa_size_t NCols() const {return ncol_; }
  inline sa_size_t Size() const {return nrow_*ncol_; }

  virtual void Print(ostream& os) const;

  inline NDataBlock<T> const & DataBlock() const { return data_; }
  inline NDataBlock<T>& DataBlock() { return data_; }

  SimpleMatrix<T>& AddCst(T c);
  SimpleMatrix<T>& MulCst(T c);

  SimpleMatrix<T>& AddElt(const SimpleMatrix<T> & b);
  SimpleMatrix<T>& MulElt(const SimpleMatrix<T> & b);

  inline SimpleMatrix<T>& operator += (T c) 
  { return AddCst(c); }
  inline SimpleMatrix<T>& operator *= (T c) 
  { return MulCst(c); }
  inline SimpleMatrix<T>& operator += (const SimpleMatrix<T> & b)
  { return AddElt(b); }
  inline SimpleMatrix<T>& operator *= (const SimpleMatrix<T> & b)
  { return MulElt(b); }

protected:
  NDataBlock<T> data_;
  T* data_p_; 
  JETExpMtx<T> dbv_;
};


template <class T>
inline ostream& operator << (ostream& os, const SimpleMatrix<T> & m)
{ m.Print(os); return(os); }

// -------------------------------------------------------------
template <class T, class E>
void JETExpression<T,E>::EvaluateAll(AnyDataObj& dest) const
  {
    //  cout << " JETExpression<T,E>::EvaluateAll(AnyDataObj& dest) " << endl;
    SimpleMatrix<T>* dmx = dynamic_cast<SimpleMatrix<T> *>(&dest);
    if (dmx == NULL) throw ParmError("SMExprMult<T>::EvaluateAll() - dest not a SimpleMatrix");
    //  cout << " --- Doing Reset() " << endl;
    e_.Reset();
    //  cout << " --- Reset() OK " << endl;
    dmx->ReSize(e_.RefSize());
    //  cout << " --- Resize OK " << endl;
    for(sa_size_t k=0; k<dmx->Size(); k++) dmx->DataBlock()(k) = e_.EvaluateAndStep();
  }


#endif