From 74a96d4eb62467e695d4e43ca574304b0d36790c Mon Sep 17 00:00:00 2001 From: Jiri Pittner Date: Tue, 9 Apr 2024 16:08:15 +0200 Subject: [PATCH] working on tensor class --- tensor.cc | 47 ++++++++++++++++++++++++++++++++++ tensor.h | 75 ++++++++++++++++--------------------------------------- 2 files changed, 69 insertions(+), 53 deletions(-) diff --git a/tensor.cc b/tensor.cc index baeb778..459e57b 100644 --- a/tensor.cc +++ b/tensor.cc @@ -26,6 +26,53 @@ namespace LA { +template +int Tensor:: calcrank() +{ +int r=0; +for(int i=0; i +LA_largeindex Tensor::calcsize() +{ +groupsizes.resize(shape.size()); +cumsizes.resize(shape.size()); +LA_largeindex s=1; +for(int i=0; i &I) //index of one subgroup { #ifdef DEBUG diff --git a/tensor.h b/tensor.h index 4dc5504..5ebd25d 100644 --- a/tensor.h +++ b/tensor.h @@ -65,6 +65,9 @@ static const LA_index offset = 0; //compiler can optimiza away some computations LA_index offset; //indices start at a general offset #endif LA_index range; //indices span this range + + bool operator==(const indexgroup &rhs) const {return number==rhs.number && symmetry==rhs.symmetry && offset==rhs.offset && range==rhs.range;}; + inline bool operator!=(const indexgroup &rhs) const {return !((*this)==rhs);}; } INDEXGROUP; template<> @@ -73,6 +76,7 @@ class LA_traits { static bool is_plaindata() {return true;}; static void copyonwrite(indexgroup& x) {}; typedef INDEXGROUP normtype; + static inline int gencmp(const indexgroup *a, const indexgroup *b, int n) {return memcmp(a,b,n*sizeof(indexgroup));}; static inline void put(int fd, const indexgroup &x, bool dimensions=1) {if(sizeof(indexgroup)!=write(fd,&x,sizeof(indexgroup))) laerror("write error 1 in indexgroup put"); } static inline void multiput(int nn, int fd, const indexgroup *x, bool dimensions=1) {if(nn*sizeof(indexgroup)!=write(fd,x,nn*sizeof(indexgroup))) laerror("write error 1 in indexgroup multiiput"); } static inline void get(int fd, indexgroup &x, bool dimensions=1) {if(sizeof(indexgroup)!=read(fd,&x,sizeof(indexgroup))) laerror("read error 1 in indexgroup get");} @@ -104,13 +108,14 @@ public: Tensor(const NRVec &s) : shape(s), data((int)calcsize()), myrank(calcrank()) {}; //general tensor Tensor(const indexgroup &g) {shape.resize(1); shape[0]=g; data.resize(calcsize()); myrank=calcrank();}; //tensor with a single index group Tensor(const Tensor &rhs): myrank(rhs.myrank), shape(rhs.shape), groupsizes(rhs.groupsizes), cumsizes(rhs.cumsizes), data(rhs.data) {}; + Tensor(int xrank, const NRVec &xshape, const NRVec &xgroupsizes, const NRVec xcumsizes, const NRVec &xdata) : myrank(xrank), shape(xshape), groupsizes(xgroupsizes), cumsizes(xcumsizes), data(xdata) {}; void clear() {data.clear();}; int rank() const {return myrank;}; int calcrank(); //is computed from shape LA_largeindex calcsize(); //set redundant data and return total size LA_largeindex size() const {return data.size();}; - void copyonwrite() {shape.copyonwrite(); data.copyonwrite();}; + void copyonwrite() {shape.copyonwrite(); groupsizes.copyonwrite(); cumsizes.copyonwrite(); data.copyonwrite();}; inline Signedpointer lhs(const SUPERINDEX &I) {int sign; LA_largeindex i=index(&sign,I); return Signedpointer(&data[i],sign);}; inline T operator()(const SUPERINDEX &I) {int sign; LA_largeindex i=index(&sign,I); if(sign==0) return 0; return sign>0 ?data[i] : -data[i];}; inline Signedpointer lhs(const FLATINDEX &I) {int sign; LA_largeindex i=index(&sign,I); return Signedpointer(&data[i],sign);}; @@ -125,69 +130,33 @@ public: inline Tensor& operator/=(const T &a) {data/=a; return *this;}; inline Tensor operator/(const T &a) const {Tensor r(*this); r /=a; return r;}; + + inline Tensor& operator+=(const Tensor &rhs) {if(shape!=rhs.shape) laerror("incompatible tensors for operation"); data+=rhs.data; return *this;} + inline Tensor& operator-=(const Tensor &rhs) {if(shape!=rhs.shape) laerror("incompatible tensors for operation"); data-=rhs.data; return *this;} + inline Tensor operator+(const Tensor &rhs) const {Tensor r(*this); r+=rhs; return r;}; + inline Tensor operator-(const Tensor &rhs) const {Tensor r(*this); r-=rhs; return r;}; + + Tensor operator-() const {return Tensor(myrank,shape,groupsizes,cumsizes,-data);}; //unary- + void put(int fd) const; void get(int fd); + inline void randomize(const typename LA_traits::normtype &x) {data.randomize(x);}; + //@@@TODO - unwinding to full size in a specified index - //@@@contraction by a whole index group + //@@@contraction by a whole index group or by individual single index //@@@TODO - contractions - basic and efficient? first contraction in a single index; between a given group+index in group at each tensor - //@@@ dvojite rekurzivni loopover s callbackem - nebo iterator s funkci next??? - //@@@ stream i/o na zaklade tohoto - //@@@permuteindexgroups //@@@symmetrize a group, antisymmetrize a group, expand a (anti)symmetric grtoup - obecne symmetry change krome +1 na -1 vse mozne - //@@@outer product - //@@@explicit constructors from vec mat smat and dense fourindex - //@@@@@@+= -= + - on same shape - //@@@@@@ randomize + //@@@outer product and product with a contraction + //@@@@permuteindexgroups + //@@@@@@explicit constructors from vec mat smat and dense fourindex + //@@@@@@ dvojite rekurzivni loopover s callbackem - nebo iterator s funkci next??? + //@@@@@@ stream i/o na zaklade tohoto }; -template -int Tensor:: calcrank() -{ -int r=0; -for(int i=0; i -LA_largeindex Tensor::calcsize() -{ -groupsizes.resize(shape.size()); -cumsizes.resize(shape.size()); -LA_largeindex s=1; -for(int i=0; i