lift comba limit for s_mp_sqr_comba

minad · minad · commit c3fe1be86a97 · 2019-11-10T16:46:53.000+01:00
diff --git a/demo/test.c b/demo/test.c
@@ -1929,7 +1929,7 @@ static int test_s_mp_sqr_karatsuba(void)
    for (size = MP_SQR_KARATSUBA_CUTOFF; size < (MP_SQR_KARATSUBA_CUTOFF + 20); size++) {
       DO(mp_rand(&a, size));
       DO(s_mp_sqr_karatsuba(&a, &b));
-      DO(s_mp_sqr(&a, &c));
+      DO(s_mp_sqr_comba(&a, &c));
       if (mp_cmp(&b, &c) != MP_EQ) {
          fprintf(stderr, "Karatsuba squaring failed at size %d\n", size);
          goto LBL_ERR;
@@ -2002,7 +2002,7 @@ static int test_s_mp_sqr_toom(void)
    for (size = MP_SQR_TOOM_CUTOFF; size < (MP_SQR_TOOM_CUTOFF + 20); size++) {
       DO(mp_rand(&a, size));
       DO(s_mp_sqr_toom(&a, &b));
-      DO(s_mp_sqr(&a, &c));
+      DO(s_mp_sqr_comba(&a, &c));
       if (mp_cmp(&b, &c) != MP_EQ) {
          fprintf(stderr, "Toom-Cook 3-way squaring failed at size %d\n", size);
          goto LBL_ERR;
diff --git a/etc/tune.c b/etc/tune.c
@@ -133,7 +133,7 @@ static uint64_t s_time_sqr(int size)
          goto LBL_ERR;
       }
       if (s_check_result == 1) {
-         if ((e = s_mp_sqr(&a,&c)) != MP_OKAY) {
+         if ((e = s_mp_sqr_comba(&a,&c)) != MP_OKAY) {
             t1 = UINT64_MAX;
             goto LBL_ERR;
          }
diff --git a/libtommath_VS2008.vcproj b/libtommath_VS2008.vcproj
@@ -912,10 +912,6 @@
 			RelativePath="s_mp_rand_platform.c"
 			>
 		</File>
-		<File
-			RelativePath="s_mp_sqr.c"
-			>
-		</File>
 		<File
 			RelativePath="s_mp_sqr_comba.c"
 			>
diff --git a/makefile b/makefile
@@ -48,8 +48,8 @@ s_mp_div_recursive.o s_mp_div_school.o s_mp_div_small.o s_mp_exptmod.o s_mp_expt
 s_mp_invmod.o s_mp_invmod_odd.o s_mp_log.o s_mp_log_d.o s_mp_log_pow2.o s_mp_montgomery_reduce_comba.o \
 s_mp_mul_balance.o s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o \
 s_mp_mul_toom.o s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o s_mp_rand_jenkins.o \
-s_mp_rand_platform.o s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o \
-s_mp_zero_buf.o s_mp_zero_digs.o
+s_mp_rand_platform.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o \
+s_mp_zero_digs.o
 
 #END_INS
 
diff --git a/makefile.mingw b/makefile.mingw
@@ -50,8 +50,8 @@ s_mp_div_recursive.o s_mp_div_school.o s_mp_div_small.o s_mp_exptmod.o s_mp_expt
 s_mp_invmod.o s_mp_invmod_odd.o s_mp_log.o s_mp_log_d.o s_mp_log_pow2.o s_mp_montgomery_reduce_comba.o \
 s_mp_mul_balance.o s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o \
 s_mp_mul_toom.o s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o s_mp_rand_jenkins.o \
-s_mp_rand_platform.o s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o \
-s_mp_zero_buf.o s_mp_zero_digs.o
+s_mp_rand_platform.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o \
+s_mp_zero_digs.o
 
 HEADERS_PUB=tommath.h
 HEADERS=tommath_private.h tommath_class.h tommath_superclass.h tommath_cutoffs.h $(HEADERS_PUB)
diff --git a/makefile.msvc b/makefile.msvc
@@ -43,8 +43,8 @@ s_mp_div_recursive.obj s_mp_div_school.obj s_mp_div_small.obj s_mp_exptmod.obj s
 s_mp_invmod.obj s_mp_invmod_odd.obj s_mp_log.obj s_mp_log_d.obj s_mp_log_pow2.obj s_mp_montgomery_reduce_comba.obj \
 s_mp_mul_balance.obj s_mp_mul_comba.obj s_mp_mul_high.obj s_mp_mul_high_comba.obj s_mp_mul_karatsuba.obj \
 s_mp_mul_toom.obj s_mp_prime_is_divisible.obj s_mp_prime_tab.obj s_mp_radix_map.obj s_mp_rand_jenkins.obj \
-s_mp_rand_platform.obj s_mp_sqr.obj s_mp_sqr_comba.obj s_mp_sqr_karatsuba.obj s_mp_sqr_toom.obj s_mp_sub.obj \
-s_mp_zero_buf.obj s_mp_zero_digs.obj
+s_mp_rand_platform.obj s_mp_sqr_comba.obj s_mp_sqr_karatsuba.obj s_mp_sqr_toom.obj s_mp_sub.obj s_mp_zero_buf.obj \
+s_mp_zero_digs.obj
 
 HEADERS_PUB=tommath.h
 HEADERS=tommath_private.h tommath_class.h tommath_superclass.h tommath_cutoffs.h $(HEADERS_PUB)
diff --git a/makefile.shared b/makefile.shared
@@ -45,8 +45,8 @@ s_mp_div_recursive.o s_mp_div_school.o s_mp_div_small.o s_mp_exptmod.o s_mp_expt
 s_mp_invmod.o s_mp_invmod_odd.o s_mp_log.o s_mp_log_d.o s_mp_log_pow2.o s_mp_montgomery_reduce_comba.o \
 s_mp_mul_balance.o s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o \
 s_mp_mul_toom.o s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o s_mp_rand_jenkins.o \
-s_mp_rand_platform.o s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o \
-s_mp_zero_buf.o s_mp_zero_digs.o
+s_mp_rand_platform.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o \
+s_mp_zero_digs.o
 
 #END_INS
 
diff --git a/makefile.unix b/makefile.unix
@@ -51,8 +51,8 @@ s_mp_div_recursive.o s_mp_div_school.o s_mp_div_small.o s_mp_exptmod.o s_mp_expt
 s_mp_invmod.o s_mp_invmod_odd.o s_mp_log.o s_mp_log_d.o s_mp_log_pow2.o s_mp_montgomery_reduce_comba.o \
 s_mp_mul_balance.o s_mp_mul_comba.o s_mp_mul_high.o s_mp_mul_high_comba.o s_mp_mul_karatsuba.o \
 s_mp_mul_toom.o s_mp_prime_is_divisible.o s_mp_prime_tab.o s_mp_radix_map.o s_mp_rand_jenkins.o \
-s_mp_rand_platform.o s_mp_sqr.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o \
-s_mp_zero_buf.o s_mp_zero_digs.o
+s_mp_rand_platform.o s_mp_sqr_comba.o s_mp_sqr_karatsuba.o s_mp_sqr_toom.o s_mp_sub.o s_mp_zero_buf.o \
+s_mp_zero_digs.o
 
 HEADERS_PUB=tommath.h
 HEADERS=tommath_private.h tommath_class.h tommath_superclass.h tommath_cutoffs.h $(HEADERS_PUB)
diff --git a/mp_mul.c b/mp_mul.c
@@ -21,12 +21,8 @@ mp_err mp_mul(const mp_int *a, const mp_int *b, mp_int *c)
               (a->used >= MP_SQR_KARATSUBA_CUTOFF)) {
       err = s_mp_sqr_karatsuba(a, c);
    } else if ((a == b) &&
-              MP_HAS(S_MP_SQR_COMBA) &&
-              (a->used < (MP_MAX_COMBA / 2))) {
+              MP_HAS(S_MP_SQR_COMBA)) {
       err = s_mp_sqr_comba(a, c);
-   } else if ((a == b) &&
-              MP_HAS(S_MP_SQR)) {
-      err = s_mp_sqr(a, c);
    } else if (MP_HAS(S_MP_MUL_BALANCE) &&
               /* Check sizes. The smaller one needs to be larger than the Karatsuba cut-off.
                * The bigger one needs to be at least about one MP_MUL_KARATSUBA_CUTOFF bigger
diff --git a/s_mp_sqr.c b/s_mp_sqr.c
diff --git a/s_mp_sqr_comba.c b/s_mp_sqr_comba.c
@@ -16,7 +16,7 @@ After that loop you do the squares and add them in.
 mp_err s_mp_sqr_comba(const mp_int *a, mp_int *b)
 {
    int       oldused, pa, ix;
-   mp_word   W1;
+   mp_digit c0, c1, c2;
    mp_err err;
    mp_int tmp, *b_;
 
@@ -31,10 +31,9 @@ mp_err s_mp_sqr_comba(const mp_int *a, mp_int *b)
    }
 
    /* number of output digits to produce */
-   W1 = 0;
+   c0 = c1 = c2 = 0;
    for (ix = 0; ix < pa; ix++) {
       int      tx, ty, iy, iz;
-      mp_word  W = 0;
 
       /* get offsets into the two bignums */
       ty = MP_MIN(a->used-1, ix);
@@ -53,22 +52,33 @@ mp_err s_mp_sqr_comba(const mp_int *a, mp_int *b)
 
       /* execute loop */
       for (iz = 0; iz < iy; iz++) {
-         W += (mp_word)a->dp[tx + iz] * (mp_word)a->dp[ty - iz];
+         mp_word t = (mp_word)a->dp[tx + iz] * (mp_word)a->dp[ty - iz];
+         int j;
+         for (j = 0; j < 2; ++j) {
+            mp_word w = (mp_word)c0 + t;
+            c0 = (mp_digit)(w & MP_MASK);
+            w = (mp_word)c1 + (w >> MP_DIGIT_BIT);
+            c1 = (mp_digit)(w & MP_MASK);
+            c2 += (mp_digit)(w >> MP_DIGIT_BIT);
+         }
       }
 
-      /* double the inner product and add carry */
-      W = W + W + W1;
-
       /* even columns have the square term in them */
       if (((unsigned)ix & 1u) == 0u) {
-         W += (mp_word)a->dp[ix>>1] * (mp_word)a->dp[ix>>1];
+         mp_word w = (mp_word)c0 + ((mp_word)a->dp[ix / 2] * (mp_word)a->dp[ix / 2]);
+         c0 = (mp_digit)(w & MP_MASK);
+         w = (mp_word)c1 + (w >> MP_DIGIT_BIT);
+         c1 = (mp_digit)(w & MP_MASK);
+         c2 += (mp_digit)(w >> MP_DIGIT_BIT);
       }
 
-      /* store it */
-      b_->dp[ix] = (mp_digit)W & MP_MASK;
+      /* store term */
+      b_->dp[ix] = c0;
 
       /* make next carry */
-      W1 = W >> (mp_word)MP_DIGIT_BIT;
+      c0 = c1;
+      c1 = c2;
+      c2 = 0;
    }
 
    /* setup dest */
diff --git a/tommath_class.h b/tommath_class.h
@@ -161,7 +161,6 @@
 #   define S_MP_RADIX_MAP_C
 #   define S_MP_RAND_JENKINS_C
 #   define S_MP_RAND_PLATFORM_C
-#   define S_MP_SQR_C
 #   define S_MP_SQR_COMBA_C
 #   define S_MP_SQR_KARATSUBA_C
 #   define S_MP_SQR_TOOM_C
@@ -544,7 +543,6 @@
 #   define S_MP_MUL_COMBA_C
 #   define S_MP_MUL_KARATSUBA_C
 #   define S_MP_MUL_TOOM_C
-#   define S_MP_SQR_C
 #   define S_MP_SQR_COMBA_C
 #   define S_MP_SQR_KARATSUBA_C
 #   define S_MP_SQR_TOOM_C
@@ -1204,14 +1202,6 @@
 #if defined(S_MP_RAND_PLATFORM_C)
 #endif
 
-#if defined(S_MP_SQR_C)
-#   define MP_CLAMP_C
-#   define MP_CLEAR_C
-#   define MP_GROW_C
-#   define MP_INIT_SIZE_C
-#   define S_MP_ZERO_DIGS_C
-#endif
-
 #if defined(S_MP_SQR_COMBA_C)
 #   define MP_CLAMP_C
 #   define MP_CLEAR_C
diff --git a/tommath_private.h b/tommath_private.h
@@ -129,6 +129,7 @@ extern void MP_FREE(void *mem, size_t size);
 
 #define MP_SIZEOF_BITS(type)    ((size_t)CHAR_BIT * sizeof(type))
 
+/* TODO: Remove MP_MAX_COMBA and MP_WARRAY */
 #define MP_MAX_COMBA            (int)(1uL << (MP_SIZEOF_BITS(mp_word) - (2u * (size_t)MP_DIGIT_BIT)))
 #define MP_WARRAY               (int)(1uL << ((MP_SIZEOF_BITS(mp_word) - (2u * (size_t)MP_DIGIT_BIT)) + 1u))
 
@@ -183,7 +184,6 @@ MP_PRIVATE mp_err s_mp_mul_karatsuba(const mp_int *a, const mp_int *b, mp_int *c
 MP_PRIVATE mp_err s_mp_mul_toom(const mp_int *a, const mp_int *b, mp_int *c) MP_WUR;
 MP_PRIVATE mp_err s_mp_prime_is_divisible(const mp_int *a, bool *result) MP_WUR;
 MP_PRIVATE mp_err s_mp_rand_platform(void *p, size_t n) MP_WUR;
-MP_PRIVATE mp_err s_mp_sqr(const mp_int *a, mp_int *b) MP_WUR;
 MP_PRIVATE mp_err s_mp_sqr_comba(const mp_int *a, mp_int *b) MP_WUR;
 MP_PRIVATE mp_err s_mp_sqr_karatsuba(const mp_int *a, mp_int *b) MP_WUR;
 MP_PRIVATE mp_err s_mp_sqr_toom(const mp_int *a, mp_int *b) MP_WUR;

Original file line number	Diff line number	Diff line change
`@@ -133,7 +133,7 @@ static uint64_t s_time_sqr(int size)`
`133`	`133`	`goto LBL_ERR;`
`134`	`134`	`}`
`135`	`135`	`if (s_check_result == 1) {`
`136`		`- if ((e = s_mp_sqr(&a,&c)) != MP_OKAY) {`
	`136`	`+ if ((e = s_mp_sqr_comba(&a,&c)) != MP_OKAY) {`
`137`	`137`	`t1 = UINT64_MAX;`
`138`	`138`	`goto LBL_ERR;`
`139`	`139`	`}`
Original file line number	Diff line number	Diff line change
`@@ -912,10 +912,6 @@`
`912`	`912`	`RelativePath="s_mp_rand_platform.c"`
`913`	`913`	`>`
`914`	`914`	`</File>`
`915`		`- <File`
`916`		`- RelativePath="s_mp_sqr.c"`
`917`		`- >`
`918`		`- </File>`
`919`	`915`	`<File`
`920`	`916`	`RelativePath="s_mp_sqr_comba.c"`
`921`	`917`	`>`