source: Sophya/trunk/SophyaLib/TArray/basarr.h@ 890

// This may look like C code, but it is really -*- C++ -*-
//      Base array class - Memory organisation management
//                     R. Ansari, C.Magneville   03/2000

#ifndef BaseArray_SEEN
#define BaseArray_SEEN

#include "machdefs.h"
#include <math.h>
#include <iostream.h>
#include "anydataobj.h"
#include "dvlist.h"


//! Maximum number of dimensions for an array
#define BASEARRAY_MAXNDIMS  5

namespace SOPHYA {  

//   ------------ classe template Array ----------- 
//! Base class for template arrays
/*!
  Define base methods, enum and defaults for TArray , TMatrix and TVector
*/
class BaseArray : public AnyDataObj {
public:
  //! To define Array or Matrix memory mapping
  enum MemoryMapping {
    AutoMemoryMapping = -1,  //!< define Auto Memory Mapping
    SameMemoryMapping = 0,   //!< define Same Memory Mapping
    CMemoryMapping = 1,      //!< define C Memory Mapping
    FortranMemoryMapping = 2 //!< define Fortran Memory Mapping
    };
  //! To define Vector type
  enum VectorType {
    AutoVectorType = -1,  //!< define Auto Vector Type
    SameVectorType = 0,   //!< define Same Vector Type
    ColumnVector = 1,     //!< define Column Vector Type
    RowVector = 2         //!< define Row Vector Type
    };

  // threshold for parallel routine call
  static void    SetOpenMPSizeThreshold(uint_8 thr=200000);
  //! Get Size threshold for parallel routine call
  static inline uint_8 GetOpenMPSizeThreshold() { return openmp_size_threshold; }

  static void    SetMaxPrint(uint_4 nprt=50, uint_4 lev=0);   
  //! Get maximum number of printed elements
  static inline uint_4 GetMaxPrint() { return max_nprt_; }
  //! Maximum number of printed elements arint level
  static inline uint_4 GetPrintLevel() { return prt_lev_; }

  //! Set Default Memory Mapping
  static short   SetDefaultMemoryMapping(short mm=CMemoryMapping);
  //! Get Default Memory Mapping
  static inline short GetDefaultMemoryMapping() { return default_memory_mapping; }
  //! Set Default Vector Type
  static short   SetDefaultVectorType(short vt=ColumnVector);
  //! Get Default Vector Type
  static inline short GetDefaultVectorType() { return default_vector_type; }
  
  // Creator / destructor
  BaseArray();
  virtual ~BaseArray();

  // Returns true if ndim and sizes are equal
  virtual bool CompareSizes(const BaseArray& a);

  // Compacts \b size=1 array dimensions
  virtual void CompactAllDim(); // suppresses all size==1 dimensions 
  virtual void CompactTrailingDim(); // suppresses size==1 dimensions after the last size>1 dimension

  // Array dimensions
  //! Return number of dimensions
  inline uint_4 NbDimensions() const { return( ndim_ ); }

  //! Return total size of the array
  inline uint_8 Size() const { return(totsize_); }
  //! Return size along the first dimension
  inline uint_4 SizeX() const { return(size_[0]); }
  //! Return size along the second dimension
  inline uint_4 SizeY() const { return(size_[1]); }
  //! Return size along the third dimension
  inline uint_4 SizeZ() const { return(size_[2]); }
  //! Return size along the \b ka th dimension
  inline uint_4 Size(int ka) const { return(size_[CheckDI(ka,1)]); }

  uint_4 MaxSizeKA() const ;

  //! Get memory organization
  inline short GetMemoryMapping() const 
         { return ( (marowi_ == 1) ? CMemoryMapping : FortranMemoryMapping) ; }
  //! line index dimension
  inline uint_4 RowsKA() const {return marowi_; }
  //! line column dimension
  inline uint_4 ColsKA() const {return macoli_; }
  //! Index dimension of the elements of a vector
  inline uint_4 VectKA() const {return veceli_; }
  void          SetMemoryMapping(short mm=AutoMemoryMapping);

  //! Get Vector type ( \b Line or \b Column vector )
  inline short  GetVectorType() const 
                { return((marowi_ == veceli_) ? ColumnVector : RowVector); }
  void SetVectorType(short vt=AutoVectorType);

  // memory organisation - packing information
  //! return true if array is packed in memory
  inline bool   IsPacked() const { return(moystep_ == 1); }
  //! return true if array is packed along the first dimension
  inline bool   IsPackedX() const { return(step_[0] == 1); }
  //! return true if array is packed along the second dimension
  inline bool   IsPackedY() const { return(step_[1] == 1); }
  //! return true if array is packed along the third dimension
  inline bool   IsPackedZ() const { return(step_[2] == 1); }
  //! return true if array is packed along the \b ka th dimension
  inline bool   IsPacked(int ka) const { return(step_[CheckDI(ka,2)] == 1); }

  //! return the minimum step value along all the dimensions
  inline uint_4 MinStep() const  { return(minstep_); }
  //! return the average step value along all the dimensions
  inline uint_4 AvgStep() const  { return(moystep_); }
  //! return the step along the first dimension
  inline uint_4 StepX() const { return(step_[0]); }
  //! return the step along the second dimension
  inline uint_4 StepY() const { return(step_[1]); }
  //! return the step along the third dimension
  inline uint_4 StepZ() const { return(step_[2]); }
  //! return the step along the \b ka th dimension
  inline uint_4 Step(int ka) const { return(step_[CheckDI(ka,3)]); }

  uint_4 MinStepKA() const ;

  // Offset of element ip
  uint_8 Offset(uint_8 ip=0) const ;
  // Offset of the i'th vector along axe ka 
  uint_8  Offset(uint_4 ka, uint_8 i) const ;
  inline uint_8  Offset(uint_4 ix, uint_4 iy, uint_4 iz, uint_4 it=0, uint_4 iu=0) const;  

  // an abstract element acces methode
  virtual double ValueAtPosition(uint_8 ip) const = 0;

 // Impression, I/O, ...
  void           Show(ostream& os, bool si=false) const;
  //! Show information on \b cout
  inline void    Show() const { Show(cout); }
  virtual string InfoString() const;

  //  Objet DVList info 
  DVList&       Info();

protected:
  inline int CheckDI(int ka, int msg) const ;
  inline void CheckBound(int ix, uint_4 iy, uint_4 iz, uint_4 it, uint_4 iu, int msg) const ;
  // Changing Sizes/NDim ... return true if OK
  bool UpdateSizes(uint_4 ndim, const uint_4 * siz, uint_4 step, uint_8 offset, string & exmsg);
  bool UpdateSizes(uint_4 ndim, const uint_4 * siz, const uint_4 * step, uint_8 offset, string & exmsg);
  bool UpdateSizes(const BaseArray& a, string & exmsg);
  static uint_8 ComputeTotalSize(uint_4 ndim, const uint_4 * siz, uint_4 step, uint_8 offset) ;
  //  Organisation memoire
  static short SelectMemoryMapping(short mm);
  static short SelectVectorType(short vt);
  void UpdateMemoryMapping(short mm);
  void UpdateMemoryMapping(BaseArray const & a, short mm);

  // Pour Extraction de sous-tableau 
  virtual void UpdateSubArraySizes(BaseArray & ra, uint_4 ndim, uint_4 * siz, uint_4 * pos, uint_4 * step) const;

  uint_4 ndim_; //! number of dimensions of array
  uint_4 size_[BASEARRAY_MAXNDIMS]; //! array of the size in each dimension
  uint_8 totsize_; //! Total number of elements
  //! two consecutive elements distance in a given dimension 
  uint_4 step_[BASEARRAY_MAXNDIMS];
  uint_4 minstep_; //! minimal step (in any axes)
  uint_4 moystep_; //! mean step, if == 0 --\> non regular steps
  uint_8 offset_;  //! global offset -\> position of elem[0] in DataBlock
  uint_4 marowi_; //! For matrices, Row index in dimensions
  uint_4 macoli_; //! For matrices, Column index in dimensions
  uint_4 veceli_; //! For vectors, dimension index = marowi_/macoli_ (Row/Col vectors)
  bool ck_memo_vt_; //! if true, check MemoryOrg./VectorType for CompareSize
  DVList* mInfo;    //! Infos (variables) attached to the array

  static char * ck_op_msg_[6]; //! Operation messages for CheckDI() CheckBound()
  static uint_4 max_nprt_; //! Nb maximum number of printed elements
  static uint_4 prt_lev_;  //! Print level
  static short default_memory_mapping; //! Default memory mapping
  static short default_vector_type;    //! Default vector type Row/Column
  static uint_8 openmp_size_threshold; //! Size limit for parallel routine calls
};

// --------------------------------------------------
//        Methodes inline de verification
// --------------------------------------------------
//! to verify the compatibility of the dimension index
inline int BaseArray::CheckDI(int ka, int msg)  const
{
  if ( (ka < 0) || (ka >= ndim_) ) {
    string txt = "BaseArray::CheckDimensionIndex/Error ";   txt += ck_op_msg_[msg];
    throw(RangeCheckError(txt));
  }
  return(ka);
}

//! to verify the compatibility of the indexes in all dimensions
inline void BaseArray::CheckBound(int ix, uint_4 iy, uint_4 iz, uint_4 it, uint_4 iu, int msg)  const
{
  if ( (ix >= size_[0]) || (iy >= size_[1]) || (iz > size_[2]) ||
       (it >= size_[3]) || (iu >= size_[4]) ) {
    string txt = "BaseArray::CheckArrayBound/Error ";   txt += ck_op_msg_[msg];
    throw(RangeCheckError(txt));
  }
  return;
}



// --------------------------------------------------
//        Position d'un element
// --------------------------------------------------
//! Offset of element (ix,iy,iz,it,iu)
inline uint_8 BaseArray::Offset(uint_4 ix, uint_4 iy, uint_4 iz, uint_4 it, uint_4 iu) const
{
#ifdef SO_BOUNDCHECKING
  CheckBound(ix, iy, iz, it, iu, 4);
#endif
  return ( offset_+  ix*step_[0] + iy*step_[1] + iz*step_[2] +
              it*step_[3] + iu*step_[4] );
}


} // Fin du namespace

#endif