[update] tensror: add contractProd optimized, in one thread calc

tensor_t/src/tensor_t/tensor_t.c

@@ -736,6 +736,94 @@ void tensorContractnProd_##type(tensor_##type** MM, tensor_##type *M0, tensor_##
     }\
     FREE_dM_S_ \
 }\
+\
+/* M[x0,x1,x3..xn] X M[y0,y1,y3..ym] = M[z0,z1...zp] (deep = l > 0) /exists 1<= l<...<l=n /  xl = y0,x{l+1}=y1, x{n}=yl  et zi=xi i<n-l et zj=y{j-(n-l)} j>=n-l alor p=n+m-2l\
+ M[x0,x1,x3..xl x{l+1}...xn] X M[xn,x{n-1},x{n-2}...xl y{l+1} ..ym] = M[x0,x1..xly{l+1}...y{n+m-2l}] (deep = l > 0)\
+M[[i][j]]=sum_{[k]}M0[[i][k]]*M[[k][j]]*/\
+\
+void tensorContractnProdOpt0_##type(tensor_##type** MM, tensor_##type *M0, tensor_##type *M1, size_t contractionNumber) {\
+   /* if (!checkMatchProdtensor(M0->dim, M1->dim, contractionNumber)) {\
+        prsize_tf("Deep = %d\n", contractionNumber);\
+    }*/\
+    if(checkContractProdTensorDim(M0->dim, M1->dim, contractionNumber)==0){\
+      printf("checkContractProdTensorDim %ld contractionNumber\n", contractionNumber);\
+        printDebug_dimension(M0->dim, "M0 dim");\
+        printDebug_dimension(M1->dim, "M1 dim");\
+        getchar();\
+    }\
+\
+    size_t len0 = M0->dim->size - contractionNumber;\
+    size_t len1 = M1->dim->size - contractionNumber;\
+\
+    size_t* tsub0 = malloc(sizeof(size_t) *len0);\
+    size_t* tsub1 = malloc(sizeof(size_t) *len1);\
+    size_t* tDk1 = malloc(sizeof(size_t) *contractionNumber);\
+    size_t* tDk0 = malloc(sizeof(size_t) *contractionNumber);\
+    subArray(tsub0, M0->dim->perm, 0, len0, 0);\
+    subArray(tsub1, M1->dim->perm, 0, len1, contractionNumber);\
+    subArray(tDk1, M1->dim->perm, 0, contractionNumber, 0);\
+    subArray(tDk0, M0->dim->perm, 0, contractionNumber, len0);\
+    /*printArraySzt(tsub0,len0,"tsub0");\
+    printArraySzt(tsub1,len1,"tsub1");\
+    printArraySzt(tDk0,contractionNumber,"tDk0");\
+    printArraySzt(tDk1,contractionNumber,"tDk1");*/\
+    dimension *dSub0 = init_dim(tsub0, len0);\
+    dimension *dSub1 = init_dim(tsub1, len1);\
+    dimension *dM1 = init_dim(tDk1, contractionNumber);\
+    dimension *dM0 = init_dim(tDk0, contractionNumber);\
+  /*printDebug_dimension(dSub0,"dSub0");\
+  printDebug_dimension(dSub1,"dSub1");\
+  printDebug_dimension(dM0,"dM0");\
+  printDebug_dimension(dM1,"dM1");*/\
+    dimension *dM;\
+    min_copy_dimension(&dM, dM0, dM1);\
+  /*printDebug_dimension(dM,"dM");*/\
+    \
+    dimension *dd;\
+    add_dimension(&dd, dSub0, dSub1);\
+  /*printDebug_dimension(dd,"dd");*/\
+    updateRankDim(dd);\
+    _RECREATE_TENSOR_IF_NOT_THE_SAME_DIM_OR_NULL_##type(MM,dd);\
+    tensor_##type *M= *MM;\
+\
+   \
+\
+    size_t a0_id, a1_id, n0_id, n1_id, begin0, begin1;\
+    for (size_t i = 0; i < M->dim->rank; i++) {\
+        if(endian){\
+          a0_id=i/dSub1->rank;\
+          a1_id=i%dSub1->rank;\
+          begin0=a0_id*dM->rank ;\
+          begin1= a1_id ;\
+        }\
+        else{\
+          a0_id=i%dSub0->rank;\
+          a1_id=i/dSub0->rank;\
+          begin1= a1_id*dM->rank ;\
+          begin0= a0_id ;\
+        }\
+        M->x[i] = 0;\
+        for (size_t k = 0; k < dM->rank; k++) {\
+          if(endian){\
+            /*n0_id= a0_id*dM->rank + k;*/\
+            /*n1_id= a1_id + dSub1->rank * k;*/\
+            /*M->x[i] += M0->x[begin0++] * M1->x[n1_id];*/\
+            M->x[i] += M0->x[begin0++] * M1->x[begin1];\
+            begin1 +=dSub1->rank ;\
+          }\
+          else{\
+            /*n0_id= a0_id + dSub0->rank * k;*/\
+            /*n1_id= a1_id*dM->rank + k;*/\
+            /*M->x[i] += M0->x[n0_id] * M1->x[begin1++];*/\
+            M->x[i] += M0->x[begin0] * M1->x[begin1++];\
+            begin0 += dSub0->rank ;\
+          }\
+            \
+        }\
+    }\
+    FREE_dM_S_ \
+}\
+\
 struct arg_Prod_##type{\
   type *M0x;\
   type *M1x;\

tensor_t/src/tensor_t/tensor_t.h

@@ -42,6 +42,7 @@ void tensorProdThrea2d_##type(tensor_##type **MM, tensor_##type *M0, tensor_##ty
 void tensorContractnProdThread_##type(tensor_##type **MM, tensor_##type *M0, tensor_##type *M1, size_t contractionNumber, size_t nbthread); \
 void tensorContractnPro2dThread_##type(tensor_##type **MM, tensor_##type *M0, tensor_##type *M1, size_t contractionNumber, size_t nbthread); \
 void tensorContractnProdNotOpt_##type(tensor_##type **MM, tensor_##type *M0, tensor_##type *M1, size_t contractionNumber); \
+void tensorContractnProdOpt0_##type(tensor_##type **MM, tensor_##type *M0, tensor_##type *M1, size_t contractionNumber); \
 type scalarProduct_dep_contractProd_##type(tensor_##type *M0, tensor_##type *M1, size_t nbthreads ,void (*tensorContractVar)(tensor_##type **MM, tensor_##type *M0, tensor_##type *M1, size_t contractionNumber, size_t nbthread ));\
 type scalarProduct_0_##type(tensor_##type *M0, tensor_##type *M1);\
 void init_random_x_##type(tensor_##type *M, type minR, type maxR,  int randomRange);\

tensor_t/test/is_good.c

@@ -17,7 +17,7 @@
 //#include "permutation_t/permutation_t.h"
 #include "tensor_t/tensor_t.h"
 
-#define VALGRIND_ 1
+#define VALGRIND_ 0
 
 TEST(rank){
   endian =true;
@@ -1098,13 +1098,13 @@ TEST(tensorContractnProd_TYPE_FLOAT2 ){
 #else
 
 
-  d0->perm[0]=35;
+  d0->perm[0]=335;
   d0->perm[1]=32; //3;
-  d0->perm[2]=23;
+  d0->perm[2]=43;
 
   d1->perm[0]=32;
-  d1->perm[1]=23;//3;
+  d1->perm[1]=43;//3;
-  d1->perm[2]=44;
+  d1->perm[2]=244;
 #endif
 
   updateRankDim(d0);
@@ -1135,7 +1135,68 @@ TEST(tensorContractnProd_TYPE_FLOAT2 ){
   // for(size_t i=0;i<M->dim->rank;++i)
   //    EXPECT_EQ_TYPE_FLOAT(M->x[i],MnO->x[i]);
     
-  //EXPECT_ARRAY_EQ_TYPE_FLOAT(M->x,M->dim->rank,MnO->x,MnO->dim->rank);
+  EXPECT_ARRAY_EQ_TYPE_FLOAT(M->x,M->dim->rank,MnO->x,MnO->dim->rank);
+
+  free_tensor_TYPE_FLOAT(M);
+  free_tensor_TYPE_FLOAT(MnO);
+  free_tensor_TYPE_FLOAT(M0);
+  free_tensor_TYPE_FLOAT(M1);
+
+}
+
+TEST(tensorContractnProdOpt0_TYPE_FLOAT2 ){
+  dimension *d0=create_dim(3);
+  dimension *d1=create_dim(3);
+#if VALGRIND_
+  d0->perm[0]=5;
+  d0->perm[1]=2; //3;
+  d0->perm[2]=3;
+
+  d1->perm[0]=2;
+  d1->perm[1]=3;//3;
+  d1->perm[2]=8;
+
+#else
+
+
+  d0->perm[0]=335;
+  d0->perm[1]=32; //3;
+  d0->perm[2]=43;
+
+  d1->perm[0]=32;
+  d1->perm[1]=43;//3;
+  d1->perm[2]=244;
+#endif
+
+  updateRankDim(d0);
+  updateRankDim(d1);
+
+
+  tensor_TYPE_FLOAT *M0 = CREATE_TENSOR_TYPE_FLOAT(d0);
+  tensor_TYPE_FLOAT *M1 = CREATE_TENSOR_TYPE_FLOAT(d1);
+
+  LOG("M0->dim->rank = %ld\n",M0->dim->rank);
+  LOG("M1->dim->rank = %ld\n",M1->dim->rank);
+  for(size_t i=0; i<M0->dim->rank;++i) M0->x[i]=i*0.1 +1;
+  for(size_t i=0; i<M1->dim->rank;++i) M1->x[i]=i*0.003 + 2;
+
+//  print_tensor_float(M0,"M0");
+//  print_tensor_float(M1,"M1");
+
+  tensor_TYPE_FLOAT *M=NULL;
+  tensor_TYPE_FLOAT *MnO=NULL;
+
+  tensorContractnProdOpt0_TYPE_FLOAT(&M, M0,M1,2);
+//  print_tensor_float(M,"M");
+  tensorContractnProdNotOpt_TYPE_FLOAT(&MnO, M0,M1,2);
+
+
+//  print_tensor_float(MnO,"MnO");
+ 
+  // for(size_t i=0;i<M->dim->rank;++i)
+  //    EXPECT_EQ_TYPE_FLOAT(M->x[i],MnO->x[i]);
+    
+  EXPECT_ARRAY_EQ_TYPE_FLOAT(M->x,M->dim->rank,MnO->x,MnO->dim->rank);
 
   free_tensor_TYPE_FLOAT(M);
   free_tensor_TYPE_FLOAT(MnO);