//
//  Numeric.h
//
//  Qdos - numeric scalar and vector classes
//
//  David Burgess
//  March 2000
//

#ifndef QDOS_NUMERIC_H
#define QDOS_NUMERIC_H

#include <complex>
#include "Exception.h"
#include "Var.h"
#include "Adapter.h"
#include "QObject.h"
#include "Naming.h"
#include "Sequence.h"
#include "Matrix.h"
#include "StringRep.h"

namespace QDOS {

typedef Matrix<string> SMatrix;
 
typedef Scalar<double> RScalar;
typedef Matrix<double> RMatrix;

typedef MatrixSeq<double> RMatrixSeq;
typedef ScalarSeq<double> RScalarSeq;

typedef Scalar<double_complex > CScalar;
typedef Matrix<double_complex > CMatrix;

typedef MatrixSeq<double_complex > CMatrixSeq;
typedef ScalarSeq<double_complex > CScalarSeq;

class QdRScalar;
typedef Var<QdRScalar> QdRScalar_var;

class QdRScalar : public QdNObject, public RScalar {

public:
  QdRScalar() : RScalar() {;}
  QdRScalar(const RMatrix& mat):RScalar(mat){;}
  QdRScalar(const double& val) : RScalar(val) {;}
  QdRScalar(const RScalar& rscalar):RScalar(rscalar){;}
  QdRScalar(const QdRScalar& rscalar):RScalar(rscalar){;}
  
  QdRScalar& operator=(const RMatrix& mat)
  { RScalar& ref=*this; ref=mat.scalar(); return *this;}
    
  template<class U>
  QdRScalar& operator=(const Matrix<U>& mat)
  { RScalar& ref=*this; ref=(double)mat.scalar(); return *this;}
 
  QdRScalar& operator=(const double& tval)
  { RScalar& ref=*this; ref=tval; return *this;}
  
  virtual QdObject_var create() const
  { return new QdRScalar;}
  
  virtual QdObject_var copy() const
  { return new QdRScalar(*this);}
  
  virtual void setValue(const double& d) 
  { *this=d;}
  
  virtual string srep(const StringRep& stringrep,int sequence_index=0) const 
  { return double_to_string( scalar());}
    
  virtual string object_srep( const StringRep& stringrep ) const
  { return "QdRScalar";}
  
  virtual void put(Vostream& v_os);
  virtual void get(Vistream& v_is);

  virtual QdNObject_var negate() const;
  virtual QdNObject_var add(QdNObject_var) const;
  virtual QdNObject_var sub(QdNObject_var) const;
  virtual QdNObject_var mult(QdNObject_var) const;
  virtual QdNObject_var div(QdNObject_var) const;
  virtual QdNObject_var abs() const;
  virtual QdNObject_var acos() const;
  virtual QdNObject_var asin() const;
  virtual QdNObject_var atan() const;
  virtual QdNObject_var atan2(QdNObject_var) const;
  virtual QdNObject_var cos() const;
  virtual QdNObject_var cosh() const;
  virtual QdNObject_var exp() const;
  virtual QdNObject_var invert() const;
  virtual QdNObject_var log() const;
  virtual QdNObject_var log10() const;
  virtual QdNObject_var mod() const;
  virtual QdNObject_var pow(QdNObject_var obj) const;
  virtual QdNObject_var sin() const;
  virtual QdNObject_var sinh() const;
  virtual QdNObject_var sqrt() const;

  virtual bool LT(QdNObject_var ) const;
  virtual bool LE(QdNObject_var ) const;
  virtual bool EQ(QdNObject_var ) const;
  virtual bool GE(QdNObject_var ) const;
  virtual bool GT(QdNObject_var ) const;


  QDOBJECT
};

class QdCScalar : public QdNObject, public CScalar {

public:
  QdCScalar() : CScalar() {;}
  QdCScalar(const CMatrix& mat):CScalar(mat){;}
  QdCScalar(const double_complex& val) : CScalar(val) {;}
  QdCScalar(const CScalar& cscalar):CScalar(cscalar){;}
  QdCScalar(const QdCScalar& cscalar):CScalar(cscalar){;}

  QdCScalar& operator=(const CMatrix& mat)
  { CScalar& ref=*this; ref=mat.scalar(); return *this;}

  template<class U>
  QdCScalar& operator=(const Matrix<U>& mat)
  { CScalar& ref=*this; ref=(double_complex)mat.scalar(); return *this;}

  QdCScalar& operator=(const double_complex& tval)
  { CScalar& ref=*this; ref=tval; return *this;}

  virtual QdObject_var create() const
  { return new QdCScalar;}

  virtual QdObject_var copy() const
  { return new QdCScalar(*this);}

  virtual string srep(const StringRep& stringrep,int sequence_index=0) const
  { return complex_to_string( scalar());}

  virtual string object_srep( const StringRep& stringrep ) const
  { return "QdCScalar";}

  virtual void put(Vostream& v_os);
  virtual void get(Vistream& v_is);

  virtual QdNObject_var negate() const;
  virtual QdNObject_var add(QdNObject_var) const;
  virtual QdNObject_var sub(QdNObject_var) const;
  virtual QdNObject_var mult(QdNObject_var) const;
  virtual QdNObject_var div(QdNObject_var) const;
  virtual QdNObject_var abs() const;
  virtual QdNObject_var cos() const;
  virtual QdNObject_var cosh() const;
  virtual QdNObject_var exp() const;
  virtual QdNObject_var invert() const;
  virtual QdNObject_var log() const;
  virtual QdNObject_var log10() const;
  virtual QdNObject_var mod() const;
  virtual QdNObject_var pow(QdNObject_var obj) const;
  virtual QdNObject_var sin() const;
  virtual QdNObject_var sinh() const;
  virtual QdNObject_var sqrt() const;

  QDOBJECT
};

class QdRScalarSeq;
typedef SeqVar<QdRScalarSeq> QdRScalarSeq_var;

class QdRScalarSeq : public QdNObject, public RScalarSeq {

public:
  QdRScalarSeq() : RScalarSeq() {;}
  QdRScalarSeq( int n, double x=0 ) : RScalarSeq(n,x) {;}
  template<class U>
  QdRScalarSeq( const vector<U>& v ) : RScalarSeq(v) {;}
  QdRScalarSeq( const RScalarSeq& s ) : RScalarSeq(s) {;}
  QdRScalarSeq(const QdRScalarSeq& s) : RScalarSeq(s) {;}

  template<class U>
  QdRScalarSeq& operator=(const vector<U>& v)
  { RScalarSeq& rscalarseq=*this; rscalarseq=v;return *this;}
   
  virtual bool is_sequence(void) const{ return true; }
  int sequence_size(void) const{ return size(); }

  virtual QdObject_var create() const
  { return new QdRScalarSeq;}

  virtual QdObject_var create_elt() const
  { return new QdRScalar;}

  virtual QdObject_var copy() const
  { return new QdRScalarSeq(*this);}

  virtual string srep(const StringRep& stringrep,int sequence_index=0) const
  { string s;
    if(sequence_index<0 || sequence_index>=size())
        return s;
    return double_to_string(elt(sequence_index).scalar());
  }

  virtual string object_srep( const StringRep& stringrep ) const
  { return "QdRScalarSeq";}

  virtual void put(Vostream& v_os);
  virtual void get(Vistream& v_is);

  virtual QdNObject_var negate() const;
  virtual QdNObject_var add(QdNObject_var) const;
  virtual QdNObject_var sub(QdNObject_var) const;
  virtual QdNObject_var mult(QdNObject_var) const;
  virtual QdNObject_var div(QdNObject_var) const;
  virtual QdNObject_var abs() const;
  virtual QdNObject_var acos() const;
  virtual QdNObject_var asin() const;
  virtual QdNObject_var atan() const;
  virtual QdNObject_var atan2(QdNObject_var) const;
  virtual QdNObject_var cos() const;
  virtual QdNObject_var cosh() const;
  virtual QdNObject_var exp() const;
  virtual QdNObject_var invert() const;
  virtual FunctionAdapter* join(QdObject_var x_in, QdObject_var y_in,\
       QdObject_var x_out, QdObject_var y_out,\
       QdObject_var join_method, QdObject_var gap_value, \
       QdObject_var gap_handler) const;
  virtual QdNObject_var log() const;
  virtual QdNObject_var log10() const;
  virtual QdNObject_var mod() const;
  virtual FunctionAdapter* multijoin(QdObjectSeq_var x_in, \
       QdObjectSeq_var y_in,QdObject_var x_out,QdObjectSeq_var y_out,\
       QdObjectSeq_var join_method,\
       QdObjectSeq_var gap_value, QdObjectSeq_var gap_handler,\
       int triage_cnt=-1) const
  { return new MultiJoinAdapter<RScalar>(x_in,y_in,x_out,y_out,\
           join_method,gap_value,gap_handler,triage_cnt);}

  virtual QdNObject_var pow(QdNObject_var obj) const;
  virtual QdNObject_var sin() const;
  virtual QdNObject_var sinh() const;
  virtual QdNObject_var sqrt() const;

  virtual bool LT(QdNObject_var ) const;
  virtual bool LE(QdNObject_var ) const;
  virtual bool EQ(QdNObject_var ) const;
  virtual bool GE(QdNObject_var ) const;
  virtual bool GT(QdNObject_var ) const;

  QDOBJECT
};

class QdCScalarSeq : public QdNObject, public CScalarSeq {

public:
  QdCScalarSeq() : CScalarSeq() {;}
  QdCScalarSeq( int n, double_complex x=double_complex(0.,0.) ) : \
      CScalarSeq(n,x) {;}
  template<class U>
  QdCScalarSeq( const vector<U>& v ) : CScalarSeq(v) {;}
  QdCScalarSeq( const CScalarSeq& s ) : CScalarSeq(s) {;}
  QdCScalarSeq(const QdCScalarSeq& s) : CScalarSeq(s) {;}

  template<class U>
  QdCScalarSeq& operator=(const vector<U>& v)
  { CScalarSeq& cscalarseq=*this; cscalarseq=v;return *this;}

  virtual bool is_sequence(void) const{ return true; }
  int sequence_size(void) const{ return size(); }

  virtual QdObject_var create() const
  { return new QdCScalarSeq;}

  virtual QdObject_var create_elt() const
  { return new QdCScalar;}

  virtual QdObject_var copy() const
  { return new QdCScalarSeq(*this);}

  virtual string srep(const StringRep& stringrep,int sequence_index=0) const
  { string s;
    if(sequence_index<0 || sequence_index>=size())
        return s;
    return complex_to_string(elt(sequence_index).scalar());
  }

  virtual string object_srep( const StringRep& stringrep ) const
  { return "QdCScalarSeq";}

  virtual void put(Vostream& v_os);
  virtual void get(Vistream& v_is);

  virtual QdNObject_var negate() const;
  virtual QdNObject_var add(QdNObject_var) const;
  virtual QdNObject_var sub(QdNObject_var) const;
  virtual QdNObject_var mult(QdNObject_var) const;
  virtual QdNObject_var div(QdNObject_var) const;
  virtual QdNObject_var abs() const;
  virtual QdNObject_var cos() const;
  virtual QdNObject_var cosh() const;
  virtual QdNObject_var exp() const;
  virtual QdNObject_var invert() const;
  virtual FunctionAdapter* join(QdObject_var x_in, QdObject_var y_in,\
      QdObject_var x_out, QdObject_var y_out, \
      QdObject_var join_method, QdObject_var gap_value, \
      QdObject_var gap_handler) const;
  virtual QdNObject_var log() const;
  virtual QdNObject_var log10() const;
  virtual QdNObject_var mod() const;
  virtual QdNObject_var pow(QdNObject_var obj) const;
  virtual QdNObject_var sin() const;
  virtual QdNObject_var sinh() const;
  virtual QdNObject_var sqrt() const;

  QDOBJECT
};

typedef Var<QdCScalar> QdCScalar_var;
typedef SeqVar<QdCScalarSeq> QdCScalarSeq_var;

class QdSMatrix;
typedef Var<QdSMatrix> QdSMatrix_var;

class QdSMatrix: public QdNObject, public SMatrix
{
public:
    QdSMatrix():SMatrix(){;}
    //QdSMatrix(const char* p):SMatrix(p){;}
    QdSMatrix(const string& s):SMatrix(s){;}
    QdSMatrix(const Dim& d):SMatrix(d){;}
    QdSMatrix(const Dim& d,const string& s):SMatrix(d,s){;}
    QdSMatrix(const Dim& d,const char* p):SMatrix(d,string(p)){;}
    //QdSMatrix(const Dim& d,size_t n, const char** p):SMatrix(d,n,p){;}
    QdSMatrix(const Dim& d,const vector<string,allocator<string> >& v):\
       SMatrix(d,v){;}
    //QdSMatrix(const Dim& d,const vector<char*,allocator<char*> >& v):\
    //    SMatrix(d,v){;}
    QdSMatrix(const Matrix<string>& m):SMatrix(m){;}
    template<class U> QdSMatrix(const Matrix<U>& m):SMatrix(m){;}
    QdSMatrix(const QdSMatrix& smatrix):SMatrix(smatrix){;}

    ListInitialization operator=(const string& s)
    { return SMatrix::operator=(s);}
    ListInitialization operator=(const string& s) const
    { return SMatrix::operator=(s);}
    QdSMatrix& operator=(const Matrix<string>& m)
    { assign_matrix(m);return *this;}
    const QdSMatrix& operator=(const Matrix<string>& m) const
    { assign_matrix(m);return *this;}
    template<class U>
    QdSMatrix& operator=(const vector<U>& v)
    { assign_vector(v); return *this;}
    template<class U>
    const QdSMatrix& operator=(const vector<U>& v) const
    { assign_vector(v); return *this;}
    template<class U>
    QdSMatrix& operator=(const Matrix<U>& m)
    { assign_matix(m); return *this;}
    template<class U>
    const QdSMatrix& operator=(const Matrix<U>& m) const
    { assign_matrix(m); return *this;}

    virtual QdObject_var create() const
    { return new QdSMatrix(dimdata());}

    virtual QdObject_var copy() const
    { return new QdSMatrix(*this);}

    virtual string srep(const StringRep& stringrep,int sequence_index=0) const
    { string s(stringrep.get_token(StringRep::SR_TOKEN_BEGIN_MATRIX));
      string separator=stringrep.get_token(StringRep::SR_TOKEN_SEPARATOR);
      int sz = size();
      iterator itr = begin();
      for( int i=0; i<sz; i++,itr++ )
      {
          if( i!=0 ) s += separator;
          s+= *itr;
      }
      s+=stringrep.get_token(StringRep::SR_TOKEN_END_MATRIX);
      return s;
     }

    virtual string object_srep( const StringRep& stringrep) const
    { return "QdSMatrix";}

    virtual string size_srep(const StringRep& stringrep, int sz=-1) const
    { return stringrep.matrix_dims_to_string(dimdata());}

    virtual void put(Vostream& v_os);
    virtual void get(Vistream& v_is);

    virtual QdNObject_var add(QdNObject_var) const;

   QDOBJECT
};

class QdRMatrix;
typedef Var<QdRMatrix> QdRMatrix_var;

class QdRMatrix: public QdNObject, public RMatrix
{
public:
    QdRMatrix():RMatrix(){;}
    QdRMatrix(const double& val):RMatrix(val){;}
    QdRMatrix(const vector<double>& v):\
        RMatrix(v){;}
    QdRMatrix(const Dim& d):RMatrix(d){;}
    QdRMatrix(const Dim& d,const double& val):RMatrix(d,val){;}
    QdRMatrix(const Dim& d,size_t n, const double* p):RMatrix(d,n,p){;}
    QdRMatrix(const Dim& d,const vector<double>& v):\
        RMatrix(d,v){;}
    //QdRMatrix(const Dim& d,const string& str):RMatrix(d,str){;}
    QdRMatrix(const Matrix<double>& m):RMatrix(m){;}
    template<class U> QdRMatrix(const Matrix<U>& m):RMatrix(m){;}
    QdRMatrix(const QdRMatrix& rmatrix):RMatrix(rmatrix){;}

    ListInitialization operator=(const double& d)
    { return RMatrix::operator=(d);}
    ListInitialization operator=(const double& d) const
    { return RMatrix::operator=(d);}
    QdRMatrix& operator=(const Matrix<double>& m)
    { assign_matrix(m);return *this;}
    const QdRMatrix& operator=(const Matrix<double>& m) const
    { assign_matrix(m);return *this;}
    template<class U>
    QdRMatrix& operator=(const vector<U>& v)
    { assign_vector(v); return *this;}
    template<class U>
    const QdRMatrix& operator=(const vector<U>& v) const
    { assign_vector(v); return *this;}
    template<class U>
    QdRMatrix& operator=(const Matrix<U>& m)
    { assign_matrix(m); return *this;}
    template<class U>
    const QdRMatrix& operator=(const Matrix<U>& m) const
    { assign_matrix(m); return *this;}

    virtual QdObject_var create() const
    { return new QdRMatrix(dimdata());}

    virtual QdObject_var copy() const
    { return new QdRMatrix(*this);}

    virtual void setValue(const double& d) 
    { *this=d;}

    virtual string srep(const StringRep& stringrep,int sequence_index=0) const
    { string s(stringrep.get_token(StringRep::SR_TOKEN_BEGIN_MATRIX));
      string separator=stringrep.get_token(StringRep::SR_TOKEN_SEPARATOR);
      int sz = size();
      iterator itr = begin();
      for( int i=0; i<sz; i++,itr++ )
      {
          if( i!=0 ) s += separator;
          s+= double_to_string( (*itr) );
      }
      s+=stringrep.get_token(StringRep::SR_TOKEN_END_MATRIX);
      return s;
     }

    virtual string object_srep( const StringRep& stringrep) const
    { return "QdRMatrix";}

    virtual string size_srep(const StringRep& stringrep, int sz=-1) const
    { return stringrep.matrix_dims_to_string(dimdata());}

    virtual void put(Vostream& v_os);
    virtual void get(Vistream& v_is);

    virtual QdNObject_var negate() const;
    virtual QdNObject_var add(QdNObject_var) const;
    virtual QdNObject_var sub(QdNObject_var) const;
    virtual QdNObject_var mult(QdNObject_var) const;
    virtual QdNObject_var div(QdNObject_var) const;
    virtual QdNObject_var abs() const;
    virtual QdNObject_var acos() const;
    virtual QdNObject_var asin() const;
    virtual QdNObject_var atan() const;
    virtual QdNObject_var atan2(QdNObject_var) const;
    virtual QdNObject_var cos() const;
    virtual QdNObject_var cosh() const;
    virtual QdNObject_var exp() const;
    virtual QdNObject_var inner(QdNObject_var) const;
    virtual QdNObject_var invert() const;
    virtual QdNObject_var log() const;
    virtual QdNObject_var log10() const;
    virtual QdNObject_var mod() const;
    virtual QdNObject_var outer(QdNObject_var) const;
    virtual QdNObject_var pow(QdNObject_var obj) const;
    virtual QdNObject_var sin() const;
    virtual QdNObject_var sinh() const;
    virtual QdNObject_var sqrt() const;

    virtual bool LT(QdNObject_var ) const;
    virtual bool LE(QdNObject_var ) const;
    virtual bool EQ(QdNObject_var ) const;
    virtual bool GE(QdNObject_var ) const;
    virtual bool GT(QdNObject_var ) const;

    QDOBJECT
};

class QdCMatrix: public QdNObject, public CMatrix
{
public:
    QdCMatrix():CMatrix(){;}
    QdCMatrix(const double_complex& val):CMatrix(val){;}
    QdCMatrix(const vector<double_complex >& v):CMatrix(v){;}
    QdCMatrix(const Dim& d):CMatrix(d){;}
    QdCMatrix(const Dim& d,const double_complex& val):CMatrix(d,val){;}
    QdCMatrix(const Dim& d,size_t n, const double_complex* p):CMatrix(d,n,p){;}
    QdCMatrix(const Dim& d,const vector<double_complex >& v):\
        CMatrix(d,v){;}
    //QdCMatrix(const Dim& d,const string& str):CMatrix(d,str){;}
    QdCMatrix(const Matrix<double_complex >& m):CMatrix(m){;}
    template<class U> QdCMatrix(const Matrix<U>& m):CMatrix(m){;}
    QdCMatrix(const QdCMatrix& cmatrix):CMatrix(cmatrix){;}

    ListInitialization operator=(const double_complex& c)
    { return CMatrix::operator=(c);}
    ListInitialization operator=(const double_complex& c) const
    { return CMatrix::operator=(c);}
    QdCMatrix& operator=(const Matrix<double_complex>& m)
    { assign_matrix(m);return *this;}
    const QdCMatrix& operator=(const Matrix<double_complex>& m) const
    { assign_matrix(m);return *this;}
    template<class U>
    QdCMatrix& operator=(const vector<U>& v)
    { assign_vector(v); return *this;}
    template<class U>
    const QdCMatrix& operator=(const vector<U>& v) const
    { assign_vector(v); return *this;}
    template<class U>
    QdCMatrix& operator=(const Matrix<U>& m)
    { assign_matrix(m); return *this;}
    template<class U>
    const QdCMatrix& operator=(const Matrix<U>& m) const
    { assign_matrix(m); return *this;}

    virtual QdObject_var create() const
    { return new QdCMatrix(dimdata());}

    virtual QdObject_var copy() const
    { return new QdCMatrix(*this);}

    virtual string srep(const StringRep& stringrep,int sequence_index=0) const
    { string s(stringrep.get_token(StringRep::SR_TOKEN_BEGIN_MATRIX));
      string separator=stringrep.get_token(StringRep::SR_TOKEN_SEPARATOR);
      int sz = size();
      iterator itr = begin();
      for( int i=0; i<sz; i++,itr++ )
      {
          if( i!=0 ) s += separator;
          s+= complex_to_string( (*itr) );
      }
      s+=stringrep.get_token(StringRep::SR_TOKEN_END_MATRIX);
      return s;
     }

    virtual string object_srep( const StringRep& stringrep) const
    { return "QdCMatrix";}

    virtual string size_srep(const StringRep& stringrep, int sz=-1) const
    { return stringrep.matrix_dims_to_string(dimdata());}

    virtual void put(Vostream& v_os);
    virtual void get(Vistream& v_is);

    virtual QdNObject_var negate() const;
    virtual QdNObject_var add(QdNObject_var) const;
    virtual QdNObject_var sub(QdNObject_var) const;
    virtual QdNObject_var mult(QdNObject_var) const;
    virtual QdNObject_var div(QdNObject_var) const;
    virtual QdNObject_var abs() const;
    virtual QdNObject_var cos() const;
    virtual QdNObject_var cosh() const;
    virtual QdNObject_var exp() const;
    virtual QdNObject_var inner(QdNObject_var) const;
    virtual QdNObject_var invert() const;
    virtual QdNObject_var log() const;
    virtual QdNObject_var log10() const;
    virtual QdNObject_var mod() const;
    virtual QdNObject_var outer(QdNObject_var) const;
    virtual QdNObject_var pow(QdNObject_var obj) const;
    virtual QdNObject_var sin() const;
    virtual QdNObject_var sinh() const;
    virtual QdNObject_var sqrt() const;

    QDOBJECT
};

typedef Sequence<SMatrix> SMatrixSeq;
class QdSMatrixSeq;
typedef SeqVar<QdSMatrixSeq> QdSMatrixSeq_var;

class QdSMatrixSeq: public QdNObject, public SMatrixSeq
{
public:
    QdSMatrixSeq(int num=0):SMatrixSeq(num){;}
    QdSMatrixSeq(int num, const Dim& d)
    {
		Dim_var dv=new Dim(d);
        SMatrix m(dv);
        m.set_copyoption(DEEP);
		this->assign(num, m);
    }

    QdSMatrixSeq(int num,const SMatrix& m):SMatrixSeq(num,m){;}
    template<class U>
    QdSMatrixSeq(const vector<U>& v):SMatrixSeq(v){;}
    QdSMatrixSeq(const SMatrixSeq& seq):SMatrixSeq(seq){;}
    QdSMatrixSeq(const QdSMatrixSeq& seq):SMatrixSeq(seq){;}

    virtual bool is_sequence(void) const { return true; }
    int sequence_size(void) const { return size(); }

    virtual QdObject_var create() const
    { return new QdSMatrixSeq;}

    virtual QdObject_var create_elt() const
    { if(empty())
          return new QdSMatrix;
      return new QdSMatrix(elt(0).dimdata());}

    virtual QdObject_var copy() const
    { return new QdSMatrixSeq(*this);}

    virtual string srep(const StringRep& stringrep,int sequence_index=0) const
    {
      string s;
      if(sequence_index<0 || sequence_index>=size())
          return s;

      s+=stringrep.matrix_dims_to_string(elt(sequence_index).dimdata());

      s+=stringrep.get_token(StringRep::SR_TOKEN_BEGIN_MATRIX);
      string separator=stringrep.get_token(StringRep::SR_TOKEN_SEPARATOR);
      int sz = elt(sequence_index).size();
      SMatrix::iterator itr = elt(sequence_index).begin();
      for( int i=0; i<sz; i++,itr++ )
      {
          if( i!=0 ) s += separator;
          s+= *itr;
      }
      s+=stringrep.get_token(StringRep::SR_TOKEN_END_MATRIX);
      return s;
     }

    virtual string object_srep( const StringRep& stringrep) const
    { return "QdSMatrixSeq";}

    virtual string size_srep(const StringRep& stringrep, int sz=-1) const
    { if(sz==-1) sz=sequence_size();
      if(sz==0)
          return stringrep.get_token(StringRep::SR_TOKEN_BEGIN_DIMS)+\
             stringrep.get_token(StringRep::SR_TOKEN_UNKNOWN_DIMS)+\
             stringrep.get_token(StringRep::SR_TOKEN_END_DIMS);
      string s;
      if(sz>1)
      {
        s+=stringrep.get_token(StringRep::SR_TOKEN_BEGIN_DIMS)+\
              int_to_string(sz)+\
              stringrep.get_token(StringRep::SR_TOKEN_END_DIMS);
      }
      s+=stringrep.get_token(StringRep::SR_TOKEN_BEGIN_DIMS)+\
        stringrep.get_token(StringRep::SR_TOKEN_UNKNOWN_DIMS)+\
        stringrep.get_token(StringRep::SR_TOKEN_END_DIMS);
      return s;
    }

    virtual void put(Vostream& v_os);
    virtual void get(Vistream& v_is);

    virtual QdNObject_var add(QdNObject_var) const;

    QDOBJECT
};

class QdRMatrixSeq;
typedef SeqVar<QdRMatrixSeq> QdRMatrixSeq_var;

class QdRMatrixSeq: public QdNObject, public RMatrixSeq
{
public:
    QdRMatrixSeq(int num=0):RMatrixSeq(num){;}
    QdRMatrixSeq(int num, const Dim& d):RMatrixSeq(num,d){;}
    QdRMatrixSeq(int num,const RMatrix& m):RMatrixSeq(num,m){;}
    template<class U>
    QdRMatrixSeq(const vector<U>& v):RMatrixSeq(v){;}
    QdRMatrixSeq(const RMatrixSeq& seq):RMatrixSeq(seq){;}
    QdRMatrixSeq(const QdRMatrixSeq& seq):RMatrixSeq(seq){;}

    template<class U>
    QdRMatrixSeq& operator=(const vector<U>& v)
    { RMatrixSeq& rmatrixseq=*this; rmatrixseq=v;return *this;}

    virtual bool is_sequence(void) const { return true; }
    int sequence_size(void) const { return size(); }

    virtual QdObject_var create() const
    { return new QdRMatrixSeq;}

    virtual QdObject_var create_elt() const
    { if(empty())
          return new QdRMatrix;
      return new QdRMatrix(elt(0).dimdata());}

    virtual QdObject_var copy() const
    { return new QdRMatrixSeq(*this);}

    virtual string srep(const StringRep& stringrep,int sequence_index=0) const
    {
      string s;
      if(sequence_index<0 || sequence_index>=size())
          return s;

      if(!is_homogeneous_seq())
          s+=stringrep.matrix_dims_to_string(elt(sequence_index).dimdata());

      s+=stringrep.get_token(StringRep::SR_TOKEN_BEGIN_MATRIX);
      string separator=stringrep.get_token(StringRep::SR_TOKEN_SEPARATOR);
      int sz = elt(sequence_index).size();
      RMatrix::iterator itr = elt(sequence_index).begin();
      for( int i=0; i<sz; i++,itr++ )
      {
          if( i!=0 ) s += separator;
          s+= double_to_string( (*itr) );
      }
      s+=stringrep.get_token(StringRep::SR_TOKEN_END_MATRIX);
      return s;
     }

    virtual string object_srep( const StringRep& stringrep) const
    { return "QdRMatrixSeq";}

    virtual string size_srep(const StringRep& stringrep, int sz=-1) const
    { if(sz==-1) sz=sequence_size();
      if(sz==0)
          return stringrep.get_token(StringRep::SR_TOKEN_BEGIN_DIMS)+\
             stringrep.get_token(StringRep::SR_TOKEN_UNKNOWN_DIMS)+\
             stringrep.get_token(StringRep::SR_TOKEN_END_DIMS);
      string s;
      if(sz>1)
      {
        s+=stringrep.get_token(StringRep::SR_TOKEN_BEGIN_DIMS)+\
              int_to_string(sz)+\
              stringrep.get_token(StringRep::SR_TOKEN_END_DIMS);
      }
      if(is_homogeneous_seq())
         s+=stringrep.matrix_dims_to_string(elt(0).dimdata());
      else
         s+=stringrep.get_token(StringRep::SR_TOKEN_BEGIN_DIMS)+\
             stringrep.get_token(StringRep::SR_TOKEN_UNKNOWN_DIMS)+\
             stringrep.get_token(StringRep::SR_TOKEN_END_DIMS);
      return s;
    }

    virtual void put(Vostream& v_os);
    virtual void get(Vistream& v_is);

    virtual QdNObject_var negate() const;
    virtual QdNObject_var add(QdNObject_var) const;
    virtual QdNObject_var sub(QdNObject_var) const;
    virtual QdNObject_var mult(QdNObject_var) const;
    virtual QdNObject_var div(QdNObject_var) const;
    virtual QdNObject_var abs() const;
    virtual QdNObject_var acos() const;
    virtual QdNObject_var asin() const;
    virtual QdNObject_var atan() const;
    virtual QdNObject_var atan2(QdNObject_var) const;
    virtual QdNObject_var cos() const;
    virtual QdNObject_var cosh() const;
    virtual QdNObject_var exp() const;
    virtual QdNObject_var inner(QdNObject_var) const;
    virtual QdNObject_var invert() const;
    virtual FunctionAdapter* join(QdObject_var x_in, QdObject_var y_in,\
        QdObject_var x_out, QdObject_var y_out, \
        QdObject_var join_method, QdObject_var gap_value, \
        QdObject_var gap_handler) const;
    virtual QdNObject_var log() const;
    virtual QdNObject_var log10() const;
    virtual QdNObject_var mod() const;
    virtual QdNObject_var outer(QdNObject_var) const;
    virtual QdNObject_var pow(QdNObject_var obj) const;
    virtual QdNObject_var sin() const;
    virtual QdNObject_var sinh() const;
    virtual QdNObject_var sqrt() const;

    QDOBJECT
};

class QdCMatrixSeq: public QdNObject, public CMatrixSeq
{
public:
    QdCMatrixSeq(int num=0):CMatrixSeq(num){;}
    QdCMatrixSeq(int num,const Dim& d):CMatrixSeq(num,d){;}
    QdCMatrixSeq(int num,const CMatrix& m):CMatrixSeq(num,m){;}
    template<class U>
    QdCMatrixSeq(const vector<U>& v):CMatrixSeq(v){;}
    QdCMatrixSeq(const CMatrixSeq& seq):CMatrixSeq(seq){;}
    QdCMatrixSeq(const QdCMatrixSeq& seq):CMatrixSeq(seq){;}

    template<class U>
    QdCMatrixSeq& operator=(const vector<U>& v)
    { CMatrixSeq& cmatrixseq=*this; cmatrixseq=v;return *this;}

    virtual bool is_sequence(void) const { return true; }
    virtual int sequence_size(void) const { return size(); }

    virtual QdObject_var create() const
    { return new QdCMatrixSeq();}

    virtual QdObject_var create_elt() const
    { if(empty())
          return new QdCMatrix;
      return new QdCMatrix(elt(0).dimdata());}

    virtual QdObject_var copy() const
    { return new QdCMatrixSeq(*this);}

    virtual string srep(const StringRep& stringrep,int sequence_index=0) const
    {
      string s;
      if(sequence_index<0 || sequence_index>=size())
          return s;
      if(!is_homogeneous_seq())
          s+=stringrep.matrix_dims_to_string(elt(sequence_index).dimdata());

      s+=stringrep.get_token(StringRep::SR_TOKEN_BEGIN_MATRIX);
      string separator=stringrep.get_token(StringRep::SR_TOKEN_SEPARATOR);
      int sz = elt(sequence_index).size();
      CMatrix::iterator itr = elt(sequence_index).begin();
      for( int i=0; i<sz; i++,itr++ )
      {
          if( i!=0 ) s += separator;
          s+= complex_to_string( (*itr) );
      }
      s+=stringrep.get_token(StringRep::SR_TOKEN_END_MATRIX);
      return s;
     }

    virtual string object_srep( const StringRep& stringrep) const
    { return "QdCMatrixSeq";}

    virtual string size_srep(const StringRep& stringrep, int sz=-1) const
    { if(sz==-1) sz=sequence_size();
      if(sz==0)
          return stringrep.get_token(StringRep::SR_TOKEN_BEGIN_DIMS)+\
             stringrep.get_token(StringRep::SR_TOKEN_UNKNOWN_DIMS)+\
             stringrep.get_token(StringRep::SR_TOKEN_END_DIMS);
      string s;
      if(sz>1)
      {
        s+=stringrep.get_token(StringRep::SR_TOKEN_BEGIN_DIMS)+\
              int_to_string(sz)+\
              stringrep.get_token(StringRep::SR_TOKEN_END_DIMS);
      }
      if(is_homogeneous_seq())
         s+=stringrep.matrix_dims_to_string(elt(0).dimdata());
      else
         s+=stringrep.get_token(StringRep::SR_TOKEN_BEGIN_DIMS)+\
             stringrep.get_token(StringRep::SR_TOKEN_UNKNOWN_DIMS)+\
             stringrep.get_token(StringRep::SR_TOKEN_END_DIMS);
      return s;
    }

    virtual void put(Vostream& v_os);
    virtual void get(Vistream& v_is);

    virtual QdNObject_var negate() const;
    virtual QdNObject_var add(QdNObject_var) const;
    virtual QdNObject_var sub(QdNObject_var) const;
    virtual QdNObject_var mult(QdNObject_var) const;
    virtual QdNObject_var div(QdNObject_var) const;
    virtual QdNObject_var abs() const;
    virtual QdNObject_var cos() const;
    virtual QdNObject_var cosh() const;
    virtual QdNObject_var exp() const;
    virtual QdNObject_var inner(QdNObject_var) const;
    virtual QdNObject_var invert() const;
    virtual FunctionAdapter* join(QdObject_var x_in, QdObject_var y_in,\
        QdObject_var x_out, QdObject_var y_out, \
        QdObject_var join_method, QdObject_var gap_value, \
        QdObject_var gap_handler) const;
    virtual QdNObject_var outer(QdNObject_var) const;
    virtual QdNObject_var log() const;
    virtual QdNObject_var log10() const;
    virtual QdNObject_var mod() const;
    virtual QdNObject_var pow(QdNObject_var obj) const;
    virtual QdNObject_var sin() const;
    virtual QdNObject_var sinh() const;
    virtual QdNObject_var sqrt() const;

    QDOBJECT
};

typedef Var<QdCMatrix> QdCMatrix_var;
typedef SeqVar<QdCMatrixSeq> QdCMatrixSeq_var;

template<> class Var<QdRMatrix>
{
  VAROBJECT(QdRMatrix)

  // static because does not need access to members of a particular
  // instance of this class

  template<class T2>
  static Var<QdRMatrix> narrow(Var<T2> obj)
  {
     QdRMatrixSeq_var rmats=QdRMatrixSeq_var::narrow(obj);
     QdRMatrix_var res;
     if(!rmats.is_nil())
     {
       if(rmats->sequence_size()==1)
       {
			RMatrix m(rmats[0]);
			res=Var<QdRMatrix>(new QdRMatrix(m));
			return res;
       }
     }
     res=Var<QdRMatrix>((dynamic_cast<QdRMatrix*>(obj.ptr_())));
     return res;
   }
};

template<> class Var<QdRScalar>
{
  VAROBJECT(QdRScalar)

  // static because does not need access to members of a particular
  // instance of this class

  template<class T2>
  static Var<QdRScalar> narrow(Var<T2> obj)
  {
     QdRMatrix_var rmat=QdRMatrix_var::narrow(obj);
     QdRScalar_var res;
     if(!rmat.is_nil())
     {
       if(rmat->size()==1)
       {
			RMatrix m(*rmat);
			while(m.rank()>1)
				m=m[0];
			res=Var<QdRScalar>(new QdRScalar(m[0]));
			return res;
       }
     }
     QdRScalarSeq_var ss=QdRScalarSeq_var::narrow(obj);
     if(!ss.is_nil())
     {
       if(ss->sequence_size()==1)
       {			
			res=Var<QdRScalar>(new QdRScalar(ss[0]));
			return res;
       }
     }
     res=Var<QdRScalar>((dynamic_cast<QdRScalar*>(obj.ptr_())));
     return res;
   }
};

} // namespace QDOS


#endif //  #ifndef QDOS_NUMERIC_H

