diff --git a/tensor.cc b/tensor.cc index 16ab17a..b520d8b 100644 --- a/tensor.cc +++ b/tensor.cc @@ -326,6 +326,9 @@ shape[0].symmetry=0; #ifndef LA_TENSOR_ZERO_OFFSET shape[0].offset=0; #endif +#ifdef LA_TENSOR_INDEXPOSITION +shape[0].upperindex=false; +#endif shape[0].range=x.size(); calcsize(); } @@ -342,6 +345,9 @@ if(x.nrows()==x.ncols() && !flat) shape[0].symmetry=0; #ifndef LA_TENSOR_ZERO_OFFSET shape[0].offset=0; +#endif +#ifdef LA_TENSOR_INDEXPOSITION +shape[0].upperindex=false; #endif shape[0].range=x.nrows(); } @@ -352,6 +358,9 @@ else shape[0].symmetry=0; shape[1].symmetry=0; #ifndef LA_TENSOR_ZERO_OFFSET shape[0].offset=0; shape[1].offset=0; +#endif +#ifdef LA_TENSOR_INDEXPOSITION + shape[0].upperindex=false; #endif shape[0].range=x.ncols(); shape[1].range=x.nrows(); @@ -371,6 +380,9 @@ shape[0].symmetry=1; #ifndef LA_TENSOR_ZERO_OFFSET shape[0].offset=0; #endif +#ifdef LA_TENSOR_INDEXPOSITION +shape[0].upperindex=false; +#endif shape[0].range=x.nrows(); calcsize(); } @@ -524,6 +536,9 @@ s<>x.number>>x.symmetry; #ifndef LA_TENSOR_ZERO_OFFSET s>>x.offset; #endif +#ifdef LA_TENSOR_INDEXPOSITION +s>> x.upperindex; +#endif s>>x.range; return s; } @@ -724,6 +742,9 @@ newshape[0].range=shape[group].range; #ifndef LA_TENSOR_ZERO_OFFSET newshape[0].offset = shape[group].offset; #endif +#ifdef LA_TENSOR_INDEXPOSITION +newshape[0].upperindex = shape[group].upperindex; +#endif int flatindex=0; //(group,index) in flat form for(int i=0; i=rhs1.shape[group].number) laerror("wrong index number in con if(rhsindex<0||rhsindex>=rhs.shape[rhsgroup].number) laerror("wrong index number in conntraction"); if(rhs1.shape[group].offset != rhs.shape[rhsgroup].offset) laerror("incompatible index offset in contraction"); if(rhs1.shape[group].range != rhs.shape[rhsgroup].range) laerror("incompatible index range in contraction"); +#ifdef LA_TENSOR_INDEXPOSITION +if(rhs1.shape[group].upperindex ^ rhs.shape[rhsgroup].upperindex == false) laerror("can contact only upper with lower index"); +#endif const Tensor u = conjugate1? (rhs1.unwind_index(group,index)).conjugate() : rhs1.unwind_index(group,index); const Tensor rhsu = rhs.unwind_index(rhsgroup,rhsindex); @@ -1040,6 +1070,9 @@ for(int i=0; i=rhs2.shape[il2[i].group].number) laerror("wrong index2 number in conntractions"); if(rhs1.shape[il1[i].group].offset != rhs2.shape[il2[i].group].offset) laerror("incompatible index offset in contractions"); if(rhs1.shape[il1[i].group].range != rhs2.shape[il2[i].group].range) laerror("incompatible index range in contractions"); +#ifdef LA_TENSOR_INDEXPOSITION + if(rhs1.shape[il1[i].group].upperindex ^ rhs2.shape[il2[i].group].upperindex == false) laerror("can contact only upper with lower index"); +#endif } Tensor u = rhs1.unwind_indices(il1); diff --git a/tensor.h b/tensor.h index f06e6ec..45a2c08 100644 --- a/tensor.h +++ b/tensor.h @@ -46,13 +46,16 @@ //@@@ will need to store vector of INDEX to the original tensor for the result's flatindex //@@@ will not be particularly efficient // -//@@@conversions to/from fourindex, optional negarive rande for beta spin handling -//@@@ optional distinguish covariant and contravariant check in contraction +//@@@conversions to/from fourindex, optional negarive range for beta spin handling // //@@@?general permutation of individual indices - check the indices in sym groups remain adjacent, calculate result's shape, loopover the result and permute using unwind_callback // // + +//do not distinguish covariant/contravariant indices +#undef LA_TENSOR_INDEXPOSITION + namespace LA { @@ -91,8 +94,15 @@ static const LA_index offset = 0; //compiler can optimize away some computations LA_index offset; //indices start at a general offset #endif LA_index range; //indices span this range +#ifdef LA_TENSOR_INDEXPOSITION +bool upperindex; +#endif - bool operator==(const indexgroup &rhs) const {return number==rhs.number && symmetry==rhs.symmetry && offset==rhs.offset && range==rhs.range;}; + bool operator==(const indexgroup &rhs) const {return number==rhs.number && symmetry==rhs.symmetry && offset==rhs.offset && range==rhs.range +#ifdef LA_TENSOR_INDEXPOSITION + && upperindex == rhs.upperindex +#endif + ;}; inline bool operator!=(const indexgroup &rhs) const {return !((*this)==rhs);}; } INDEXGROUP;