*** empty log message ***

fourindex.h

@@ -524,6 +524,9 @@ istream& operator>>(istream  &s, fourindex<I,T> &x)
 //we use a general template forward declaration, but then it has to be done differently for (almost) each case
 //by means of partial template specialization
 
+//note - loops for the twoelectronrealmullikan integral to be unique and in canonical order
+// i=1..n, j=1..i, k=1..i, l=1..(i==k?j:k)
+
 template<fourindexsymtype S, class T, class DUMMY> class fourindex_dense;
 
 //make it as a derived class in order to be able to use it in a base class context - "supermatrix" operations

smat.h

@@ -283,12 +283,38 @@ inline T SMat_index(T i, T j)
 return (i>=j) ? i*(i+1)/2+j : j*(j+1)/2+i;
 }
 
+template<typename T>
+inline T SMat_index_igej(T i, T j)
+{
+return i*(i+1)/2+j;
+}
+
+template<typename T>
+inline T SMat_index_ilej(T i, T j)
+{
+return j*(j+1)/2+i;
+}
+
+
 template<typename T>
 inline T SMat_index_1(T i, T j)
 {
 return (i>=j)? i*(i-1)/2+j-1 : j*(j-1)/2+i-1;
 }
 
+template<typename T>
+inline T SMat_index_1igej(T i, T j)
+{
+return i*(i-1)/2+j-1;
+}
+
+template<typename T>
+inline T SMat_index_1ilej(T i, T j)
+{
+return j*(j-1)/2+i-1;
+}
+
+
 // access the element, 2-dim array case
 template <typename T>
 inline T & NRSMat<T>::operator()(const int i, const int j)