jiri/LA_library
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

tensor: indexmatrix() and zero index detection

Browse Source
This commit is contained in:
Jiri Pittner
2026-01-27 17:41:47 +01:00
parent 1580891639
commit 1d53afd257
3 changed files with 68 additions and 11 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
t.cc
View File

@@ -4597,6 +4597,8 @@ Tensor<double> x(shape); x.randomize(1.);
x.defaultnames(); x.defaultnames();
cout <<"x= "<<x.shape << " "<<x.names<<endl; cout <<"x= "<<x.shape << " "<<x.names<<endl;
cout <<"indexmatrix of x = "<<x.indexmatrix();
NRVec<INDEXGROUP> yshape(2); NRVec<INDEXGROUP> yshape(2);
yshape[0].number=1; yshape[0].number=1;
yshape[0].symmetry=0; yshape[0].symmetry=0;

62
tensor.cc
View File

@@ -417,7 +417,7 @@ calcsize();
template<typename T> template<typename T>
void loopingroups(Tensor<T> &t, int ngroup, int igroup, T **p, SUPERINDEX &I, void (*callback)(const SUPERINDEX &, T *)) void loopingroups(Tensor<T> &t, int ngroup, int igroup, T **p, SUPERINDEX &I, void (*callback)(const SUPERINDEX &, T *), bool skipzeros)
{ {
LA_index istart,iend; LA_index istart,iend;
const INDEXGROUP *sh = &(* const_cast<const NRVec<INDEXGROUP> *>(&t.shape))[ngroup]; const INDEXGROUP *sh = &(* const_cast<const NRVec<INDEXGROUP> *>(&t.shape))[ngroup];
@@ -443,6 +443,8 @@ switch(sh->symmetry)
for(LA_index i = istart; i<=iend; ++i) for(LA_index i = istart; i<=iend; ++i)
{ {
if(skipzeros && i==0) continue;
I[ngroup][igroup]=i; I[ngroup][igroup]=i;
if(ngroup==0 && igroup==0) if(ngroup==0 && igroup==0)
{ {
@@ -460,14 +462,14 @@ for(LA_index i = istart; i<=iend; ++i)
const INDEXGROUP *sh2 = &(* const_cast<const NRVec<INDEXGROUP> *>(&t.shape))[newngroup]; const INDEXGROUP *sh2 = &(* const_cast<const NRVec<INDEXGROUP> *>(&t.shape))[newngroup];
newigroup=sh2->number-1; newigroup=sh2->number-1;
} }
loopingroups(t,newngroup,newigroup,p,I,callback); loopingroups(t,newngroup,newigroup,p,I,callback,skipzeros);
} }
} }
} }
template<typename T> template<typename T>
void Tensor<T>::loopover(void (*callback)(const SUPERINDEX &, T *)) void Tensor<T>::loopover(void (*callback)(const SUPERINDEX &, T *), bool skipzeros)
{ {
SUPERINDEX I(shape.size()); SUPERINDEX I(shape.size());
for(int i=0; i<I.size(); ++i) for(int i=0; i<I.size(); ++i)
@@ -479,12 +481,12 @@ for(int i=0; i<I.size(); ++i)
T *pp=&data[0]; T *pp=&data[0];
int ss=shape.size()-1; int ss=shape.size()-1;
const INDEXGROUP *sh = &(* const_cast<const NRVec<INDEXGROUP> *>(&shape))[ss]; const INDEXGROUP *sh = &(* const_cast<const NRVec<INDEXGROUP> *>(&shape))[ss];
loopingroups(*this,ss,sh->number-1,&pp,I,callback); loopingroups(*this,ss,sh->number-1,&pp,I,callback,skipzeros);
} }
template<typename T> template<typename T>
void constloopingroups(const Tensor<T> &t, int ngroup, int igroup, const T **p, SUPERINDEX &I, void (*callback)(const SUPERINDEX &, const T *)) void constloopingroups(const Tensor<T> &t, int ngroup, int igroup, const T **p, SUPERINDEX &I, void (*callback)(const SUPERINDEX &, const T *), bool skipzeros)
{ {
LA_index istart,iend; LA_index istart,iend;
const INDEXGROUP *sh = &t.shape[ngroup]; const INDEXGROUP *sh = &t.shape[ngroup];
@@ -510,6 +512,8 @@ switch(sh->symmetry)
for(LA_index i = istart; i<=iend; ++i) for(LA_index i = istart; i<=iend; ++i)
{ {
if(skipzeros && i==0) continue;
I[ngroup][igroup]=i; I[ngroup][igroup]=i;
if(ngroup==0 && igroup==0) if(ngroup==0 && igroup==0)
{ {
@@ -527,14 +531,14 @@ for(LA_index i = istart; i<=iend; ++i)
const INDEXGROUP *sh2 = &(* const_cast<const NRVec<INDEXGROUP> *>(&t.shape))[newngroup]; const INDEXGROUP *sh2 = &(* const_cast<const NRVec<INDEXGROUP> *>(&t.shape))[newngroup];
newigroup=sh2->number-1; newigroup=sh2->number-1;
} }
constloopingroups(t,newngroup,newigroup,p,I,callback); constloopingroups(t,newngroup,newigroup,p,I,callback,skipzeros);
} }
} }
} }
template<typename T> template<typename T>
void Tensor<T>::constloopover(void (*callback)(const SUPERINDEX &, const T *)) const void Tensor<T>::constloopover(void (*callback)(const SUPERINDEX &, const T *), bool skipzeros) const
{ {
SUPERINDEX I(shape.size()); SUPERINDEX I(shape.size());
for(int i=0; i<I.size(); ++i) for(int i=0; i<I.size(); ++i)
@@ -546,7 +550,28 @@ for(int i=0; i<I.size(); ++i)
const T *pp=&data[0]; const T *pp=&data[0];
int ss=shape.size()-1; int ss=shape.size()-1;
const INDEXGROUP *sh = &shape[ss]; const INDEXGROUP *sh = &shape[ss];
constloopingroups(*this,ss,sh->number-1,&pp,I,callback); constloopingroups(*this,ss,sh->number-1,&pp,I,callback,skipzeros);
}
static INDEXMATRIX *indexmat_p;
static LA_largeindex indexmat_row;
template<typename T>
static void indexmatrix_callback(const SUPERINDEX &I, const T *p)
{
FLATINDEX f=superindex2flat(I);
indexmat_p->rowset(f,indexmat_row++);
}
template<typename T>
INDEXMATRIX Tensor<T>::indexmatrix() const
{
INDEXMATRIX r;
r.resize(size(),rank());
indexmat_p = &r;
indexmat_row = 0;
constloopover(indexmatrix_callback<T>,false);
return r;
} }
@@ -2439,6 +2464,27 @@ return true;
} }
bool zero_in_index(const FLATINDEX &I)
{
for(int i=0; i<I.size(); ++i) if(I[i]==0) return true;
return false;
}
bool zero_in_index(const INDEXMATRIX &m, const LA_largeindex row)
{
const LA_index *p = &m(row,0);
for(int i=0; i<m.ncols(); ++i) if(p[i]==0) return true;
return false;
}
bool zero_in_index(const SUPERINDEX &I)
{
for(int i=0; i<I.size(); ++i) if(zero_in_index(I[i])) return true;
return false;
}
template class Tensor<double>; template class Tensor<double>;
template class Tensor<std::complex<double> >; template class Tensor<std::complex<double> >;

15
tensor.h
View File

@@ -190,8 +190,9 @@ class LA_traits<INDEXGROUP> {
typedef NRVec<LA_index> FLATINDEX; //all indices but in a single vector typedef NRVec<LA_index> FLATINDEX; //all indices but in a single vector
typedef NRVec<NRVec<LA_index> > SUPERINDEX; //all indices in the INDEXGROUP structure typedef NRVec<FLATINDEX> SUPERINDEX; //all indices in the INDEXGROUP structure
typedef NRVec<LA_largeindex> GROUPINDEX; //set of indices in the symmetry groups typedef NRVec<LA_largeindex> GROUPINDEX; //set of indices in the symmetry groups
typedef NRMat<LA_index> INDEXMATRIX; //list of FLATINDEXes (rows of the matrix) of all tensor elements - convenient to be able to run over the whole tensor in a for loop rather than via recursive loopovers with a callback
struct INDEX struct INDEX
{ {
int group; int group;
@@ -230,6 +231,11 @@ int flatposition(int group, int index, const NRVec<INDEXGROUP> &shape);
int flatposition(const INDEX &i, const NRVec<INDEXGROUP> &shape); //position of that index in FLATINDEX int flatposition(const INDEX &i, const NRVec<INDEXGROUP> &shape); //position of that index in FLATINDEX
INDEX indexposition(int flatindex, const NRVec<INDEXGROUP> &shape); //inverse to flatposition INDEX indexposition(int flatindex, const NRVec<INDEXGROUP> &shape); //inverse to flatposition
//useful for negative offsets and 0 index excluded
bool zero_in_index(const FLATINDEX &);
bool zero_in_index(const SUPERINDEX &);
bool zero_in_index(const INDEXMATRIX &, const LA_largeindex row);
FLATINDEX superindex2flat(const SUPERINDEX &I); FLATINDEX superindex2flat(const SUPERINDEX &I);
template<typename T> template<typename T>
@@ -361,11 +367,14 @@ public:
bool fulfills_hermiticity() const; //check it is so bool fulfills_hermiticity() const; //check it is so
inline void randomize(const typename LA_traits<T>::normtype &x) {data.randomize(x); enforce_hermiticity();}; inline void randomize(const typename LA_traits<T>::normtype &x) {data.randomize(x); enforce_hermiticity();};
void loopover(void (*callback)(const SUPERINDEX &, T *)); //loop over all elements
void constloopover(void (*callback)(const SUPERINDEX &, const T *)) const; //loop over all elements void loopover(void (*callback)(const SUPERINDEX &, T *), bool skipzeros=false); //loop over all elements, optionally skip zero indices (i.e. run over ...-2,-1,1,2...) which is useful for special applications
void constloopover(void (*callback)(const SUPERINDEX &, const T *), bool skipzeros=false) const; //loop over all elements
void grouploopover(void (*callback)(const GROUPINDEX &, T *)); //loop over all elements disregarding the internal structure of index groups void grouploopover(void (*callback)(const GROUPINDEX &, T *)); //loop over all elements disregarding the internal structure of index groups
void constgrouploopover(void (*callback)(const GROUPINDEX &, const T *)) const; //loop over all elements disregarding the internal structure of index groups void constgrouploopover(void (*callback)(const GROUPINDEX &, const T *)) const; //loop over all elements disregarding the internal structure of index groups
INDEXMATRIX indexmatrix() const; //get indexmatrix - rows store FLATINDEXes matching data[]
Tensor permute_index_groups(const NRPerm<int> &p) const; //rearrange the tensor storage permuting index groups as a whole: result_i = source_p_i Tensor permute_index_groups(const NRPerm<int> &p) const; //rearrange the tensor storage permuting index groups as a whole: result_i = source_p_i
Tensor permute_index_groups(const NRVec<INDEXNAME> &names) const; //permute to requested order of group's first indices (or permute individual indices of a flat tensor) Tensor permute_index_groups(const NRVec<INDEXNAME> &names) const; //permute to requested order of group's first indices (or permute individual indices of a flat tensor)