working on tensor - outer product

t.cc

@@ -3295,7 +3295,7 @@ cout <<eu;
 }
 
 
-if(1)
+if(0)
 {
 int n=5;
 INDEXGROUP ag;
@@ -3385,4 +3385,40 @@ cout <<e;
 cout <<eu;
 }
 
+if(0)
+{
+int n=5;
+INDEXGROUP g;
+g.number=2;
+g.symmetry= 1;
+g.offset=0;
+g.range=n;
+
+Tensor<double> e(g);
+e.randomize(1.);
+INDEXLIST il(2);
+il[0]= {0,1};
+il[1]= {0,0};
+Tensor<double> eu = e.unwind_indices(il);
+
+for(int i=0; i<n; ++i)
+        for(int j=0; j<n; ++j)
+                for(int k=0; k<n; ++k)
+                        for(int l=0; l<n; ++l)
+                                {
+				if(e(i,j)!=eu(j,i)) laerror("error in unwind_indces");
+				}
+cout <<e;
+cout <<eu;
+}
+
+if(1)
+{
+NRVec<double> a({1.,10.,100.});
+NRVec<double> b({1.,2.,3.});
+Tensor<double> aa(a);
+Tensor<double> bb(b);
+cout << aa*bb;
+}
+
 }

tensor.cc

@@ -726,11 +726,16 @@ for(int i=0; i<il.size(); ++i)
 
 //then the remaining groups with one index removed, if nonempty
 int ii=il.size();
+int emptied_groups=0;
 for(int i=0; i<oldshape.size(); ++i)
 	if(oldshape[i].number>0)
 		{
 		newshape[ii++] = oldshape[i];
 		}
+	else
+		++emptied_groups;
+
+if(emptied_groups) newshape.resize(newshape.size()-emptied_groups,true);
 
 Tensor<T> r(newshape);
 if(r.rank()!=rank()) laerror("internal error 2 in unwind_indces");
@@ -1019,6 +1024,14 @@ return r;
 }
 
 
+//outer product, rhs indices will be the less significant than this
+template<typename T>
+Tensor<T> Tensor<T>::operator*(const Tensor &rhs) const
+{
+Tensor<T> r(rhs.shape.concat(shape));
+r.data= data.otimes2vec(rhs.data);
+return r;
+}
 
 
 

tensor.h

@@ -36,9 +36,9 @@
 #include "smat.h"
 #include "miscfunc.h"
 
-//@@@todo - outer product
 //@@@permutation of individual indices??? how to treat the symmetry groups
 //@@@todo - index names and contraction by named index list
+//@@@contraction inside one tensor
 
 namespace LA {
 
@@ -159,6 +159,8 @@ 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;};
 
+	Tensor operator*(const Tensor &rhs) const; //outer product
+
 	Tensor& conjugateme() {data.conjugateme(); return *this;};
         inline Tensor conjugate() const {Tensor r(*this); r.conjugateme(); return r;};