Fixed prestore mapping.

[phys/su2clebsch.git] / su2clebsch.cpp

diff --git a/su2clebsch.cpp b/su2clebsch.cpp
index 11d506e5b1667863d24a21c147b84deb5a7c1dd7..831c9dadc6e4ee799193a477b719c85c6005a5af 100644 (file)
--- a/su2clebsch.cpp
+++ b/su2clebsch.cpp
@@ -21,24 +21,32 @@ namespace su2clebsch {
     return result;
   }
 
     return result;
   }
 

+  prestore::~prestore ()
+  {
+    delete[] cgcdata;
+  }
+

   prestore::prestore(unsigned int _max2j) 
   {
     max2j = _max2j;
   prestore::prestore(unsigned int _max2j) 
   {
     max2j = _max2j;

-    cgcdata = new double[ipow(max2j+1,6)];
+    cgcdata = new double[ipow(max2j+1,3)*ipow(2*max2j+1,3)];

     
     for (int factor1_2j=0; factor1_2j<=max2j; factor1_2j++)
     
     for (int factor1_2j=0; factor1_2j<=max2j; factor1_2j++)

-      for (int factor1_mp=0; factor1_mp<=max2j; factor1_mp++)
+      for (int factor1_mp=0; factor1_mp<=2*max2j; factor1_mp++)

        for (int factor2_2j=0; factor2_2j<=max2j; factor2_2j++)
        for (int factor2_2j=0; factor2_2j<=max2j; factor2_2j++)

-         for (int factor2_mp=0; factor2_mp<=max2j; factor2_mp++)
+         for (int factor2_mp=0; factor2_mp<=2*max2j; factor2_mp++)

            for (int irrep_2j=0; irrep_2j<=max2j; irrep_2j++)
            for (int irrep_2j=0; irrep_2j<=max2j; irrep_2j++)

-             for (int irrep_mp=0; irrep_mp<=max2j; irrep_mp++)
+             for (int irrep_mp=0; irrep_mp<=2*max2j; irrep_mp++)

                {
                {

-                 cgcdata[ factor1_2j*ipow(max2j+1,5) + factor1_mp*ipow(max2j+1,4) +
-                          factor2_2j*ipow(max2j+1,3) + factor2_mp*ipow(max2j+1,2) +
-                          irrep_2j*(max2j+1)         + irrep_mp ] 
+                 cgcdata[ factor1_2j*ipow(max2j+1,2)*ipow(2*max2j+1,3) + 
+                          factor1_mp*ipow(max2j+1,2)*ipow(2*max2j+1,2) +
+                          factor2_2j*(max2j+1)*ipow(2*max2j+1,2) + 
+                          factor2_mp*(max2j+1)*(2*max2j+1) +
+                          irrep_2j*(2*max2j+1) + 
+                          irrep_mp ] 

                    =
                    =

-                   cgc( factor1_2j, 2*factor1_mp-max2j, factor2_2j, 2*factor2_mp-max2j,
-                        irrep_2j, 2*irrep_mp-max2j );
+                   cgc( factor1_2j, factor1_mp-max2j, factor2_2j, factor2_mp-max2j,
+                        irrep_2j, irrep_mp-max2j );

                }
   }
 
                }
   }
 


@@ -59,12 +67,12 @@ namespace su2clebsch {
     if ( factor1_2j > Prestore.max2j || factor2_2j > Prestore.max2j || irrep_2j > Prestore.max2j )
       return cgc(factor1_2j, factor1_2m, factor2_2j, factor2_2m, irrep_2j, irrep_2m);
     
     if ( factor1_2j > Prestore.max2j || factor2_2j > Prestore.max2j || irrep_2j > Prestore.max2j )
       return cgc(factor1_2j, factor1_2m, factor2_2j, factor2_2m, irrep_2j, irrep_2m);
     

-    return Prestore.cgcdata[ factor1_2j*ipow(Prestore.max2j+1,5) + 
-                            (factor1_2m+Prestore.max2j)/2 *ipow(Prestore.max2j+1,4) +
-                            factor2_2j*ipow(Prestore.max2j+1,3) + 
-                            (factor2_2m+Prestore.max2j)/2 *ipow(Prestore.max2j+1,2) +
-                            irrep_2j*(Prestore.max2j+1) + 
-                            (irrep_2m+Prestore.max2j)/2 ];
+    return Prestore.cgcdata[ factor1_2j*ipow(Prestore.max2j+1,2)*ipow(2*Prestore.max2j+1,3) + 
+                            (factor1_2m+Prestore.max2j) * ipow(Prestore.max2j+1,2)*ipow(2*Prestore.max2j+1,2) +
+                            factor2_2j*(Prestore.max2j+1)*ipow(2*Prestore.max2j+1,2) + 
+                            (factor2_2m+Prestore.max2j) * (Prestore.max2j+1)*(2*Prestore.max2j+1) +
+                            irrep_2j*(2*Prestore.max2j+1) + 
+                            (irrep_2m+Prestore.max2j) ];

   } 
 
   double cgc (int factor1_2j, int factor1_2m, int factor2_2j, int factor2_2m,
   } 
 
   double cgc (int factor1_2j, int factor1_2m, int factor2_2j, int factor2_2m,