progs/perftest.c: Use from Glibc syscall numbers.

[catacomb] / math / strongprime.c

diff --git a/math/strongprime.c b/math/strongprime.c
index 9eab8b2..4ea6253 100644 (file)
--- a/math/strongprime.c
+++ b/math/strongprime.c
@@ -28,6 +28,7 @@
 /*----- Header files ------------------------------------------------------*/
 
 #include <mLib/dstr.h>
+#include <mLib/macros.h>
 
 #include "grand.h"
 #include "mp.h"
@@ -39,6 +40,9 @@
 
 /*----- Main code ---------------------------------------------------------*/
 
+/* Oh, just shut up. */
+CLANG_WARNING("-Wempty-body")
+
 /* --- @strongprime_setup@ --- *
  *
  * Arguments:  @const char *name@ = pointer to name root
@@ -121,7 +125,10 @@ mp *strongprime_setup(const char *name, mp *d, pfilt *f, unsigned nbits,
                rabin_iters(nb), pgen_test, &rb)) == 0)
     goto fail_t;
 
-  /* --- Choose a suitable value for %$r = 2it + 1$% for some %$i$% --- */
+  /* --- Choose a suitable value for %$r = 2it + 1$% for some %$i$% --- *
+   *
+   * Then %$r \equiv 1 \pmod{t}$%, i.e., %$r - 1$% is a multiple of %$t$%.
+   */
 
   rr = mp_lsl(rr, t, 1);
   pfilt_create(&c.f, rr);
@@ -137,20 +144,18 @@ mp *strongprime_setup(const char *name, mp *d, pfilt *f, unsigned nbits,
 
   /* --- Select a suitable congruence class for %$p$% --- *
    *
-   * This computes %$p_0 = 2 s (s^{r - 2} \bmod r) - 1$%.
+   * This computes %$p_0 = 2 s (s^{-1} \bmod r) - 1$%.  Then %$p_0 + 1$% is
+   * clearly a multiple of %$s$%, and
+   *
+   *   %$p_0 - 1 \equiv 2 s s^{-1} - 2 \equiv 0 \pmod{r}$%
+   *
+   * is a multiple of %$r$%.
    */
 
-  {
-    mpmont mm;
-
-    mpmont_create(&mm, q);
-    rr = mp_sub(rr, q, MP_TWO);
-    rr = mpmont_exp(&mm, rr, s, rr);
-    mpmont_destroy(&mm);
-    rr = mp_mul(rr, rr, s);
-    rr = mp_lsl(rr, rr, 1);
-    rr = mp_sub(rr, rr, MP_ONE);
-  }
+  rr = mp_modinv(rr, s, q);
+  rr = mp_mul(rr, rr, s);
+  rr = mp_lsl(rr, rr, 1);
+  rr = mp_sub(rr, rr, MP_ONE);
 
   /* --- Pick a starting point for the search --- *
    *