Dixon's factorization method implemented, but not yet tested.

Shit if that's slow.

book/dixon.tex

@@ -224,7 +224,7 @@ $e^{\sqrt{\ln N \ln \ln N}}$.
     \Require $\factorBase$, the factor base
     \Function{dixon}{\PKArg}
     \State $i \gets 0$
-    \State $r \gets |\factorBase| + 5$
+    \State $r \getsRandom \naturalN_{ > |\factorBase|}$
     \Comment finding linearity requires redundance
     \While{$i < r$}
     \Comment search for suitable pairs
@@ -236,12 +236,13 @@ $e^{\sqrt{\ln N \ln \ln N}}$.
   \State $M \gets \texttt{matrix}(v_0, \ldots, v_f)$
   \For{$\lambda = \{\mu_0, \ldots, \mu_k\}
     \strong{ in } \textsc{ker}(M)$}
-  \Comment get relations
-    \State $x \gets \prod\limits_{\mu \in \lambda} x_\mu \pmod{N}$
-    \State $y, r \gets \dsqrt{\prod\limits_{\mu \in \lambda} y_\mu \pmod{N}}$
-    \If{$\gcd(x+y, N) > 1$}
-      \State $p \gets \gcd(x+y, N)$
-      \State $q \gets \gcd(x-y, N)$
+    \Comment get relations
+    \State $x \gets \prod_{\mu \in \lambda} x_\mu \pmod{N}$
+    \State $y, r \gets \dsqrt{\prod_{\mu \in \lambda} y_\mu \pmod{N}}$
+    \State $g  \gets \gcd(x+y, N)$
+    \If{$1 < \gcd < N$}
+      \State $p \gets g $
+      \State $q \gets N//p$
       \State \Return $p, q$
     \EndIf
   \EndFor

src/questions/Makefile.am

@@ -1,3 +1,4 @@
+
 SUBDIRS = tests/
 
 # following the fucking manual, I am putting headers into the _SOURCES variable
@@ -5,7 +6,7 @@ SUBDIRS = tests/
 EXAMPLE_QUESTION = example.c
 WIENER_QUESTION = wiener.c    include/qwiener.h
 POLLARD_QUESTION = pollard.c    include/qpollard.h
-DIXON_QUESTION = dixon.c
+DIXON_QUESTION = dixon.c    include/qdixon.h
 FERMAT_QUESTION = fermat.c
 METADATA_QUESTION = metadata.c
 POLLARDRHO_QUESTION = pollardrho.c

src/questions/allquestions.c

@@ -46,4 +46,5 @@ void select_all_questions(void)
   REGISTER_QUESTION(MetadataQuestion);
   REGISTER_QUESTION(PollardRhoQuestion);
   REGISTER_QUESTION(WilliamsQuestion);
+  REGISTER_QUESTION(DixonQuestion);
 }

src/questions/dixon.c

@@ -15,13 +15,14 @@
  * while e, d . ed ≡ 1 (mod φ(N))  are respectively the public and the private
  * exponent.
  */
-
-#include <stdint.h>
+#include <assert.h>
 #include <strings.h>
 
 #include <openssl/bn.h>
 
 #include "qa/questions/questions.h"
+#include "qa/questions/primes.h"
+#include "qa/questions/qarith.h"
 #include "qa/questions/qstrings.h"
 #include "qa/questions/qdixon.h"
 
@@ -98,6 +99,7 @@ kernel(matrix_t *m)
   return h;
 }
 
+
 /**
  * \brief Check for smoothness, incuding negative numbers.
  *
@@ -117,18 +119,42 @@ int dixon_smooth(BIGNUM *y, BN_CTX *ctx, char *v, size_t len)
   return ret;
 }
 
+/**
+ * \brief Discover a number x such that x² - n is smooth.
+ *
+ */
+inline void
+discover_smooth(BIGNUM *y, BIGNUM *x, BIGNUM *n,
+                BN_CTX *ctx, char *v, size_t len)
+{
+  do {
+    BN_pseudo_rand_range(x, n);
+    /* yᵢ = xᵢ² - N */
+    BN_sqr(y, x, ctx);
+    BN_sub(y, y, n);
+
+  } while (!dixon_smooth(y, ctx, v, len));
+}
+
 static RSA*
 dixon_question_ask_rsa(const RSA *rsa)
 {
-  size_t primes = 5;
-  size_t r = primes + 1;
+  /*
+   * take exp(sqrt(ln N ln ln N))
+   *    ≅ 1.44 * 2^(sqrt(lg N lg lg N))
+   *    ≅ ³/₂ * 2^(sqrt(lg N 10)) for keys of 1024 bits.
+   */
+  size_t primes = 3 * (1 << (BN_num_bits(rsa->n) * 5)) / 2;
   size_t f = primes + 5;
+  size_t r = primes + 1;
   size_t i, j;
   RSA *ret = NULL;
   BIGNUM
-    *x, *y,
-    *sqy, *rem,
-    *gcd;
+    *x = BN_new(),
+    *y = BN_new(),
+    *sqy = BN_new(),
+    *rem = BN_new(),
+    *gcd = BN_new();
   BN_CTX *ctx = BN_CTX_new();
   struct bnpair {
     BIGNUM *x;
@@ -138,31 +164,26 @@ dixon_question_ask_rsa(const RSA *rsa)
   matrix_t *h;
 
 
-  /** STEP 1: INITIALIZATION **/
+  /** STEP 1: initialization **/
   /* plus one for the sign */
   m = matrix_new(f, r);
   R = malloc(sizeof(struct bnpair) * f);
-  for (i=0; i != r; R[i].x = BN_new(), R[i].y = BN_new(), i++);
+  for (i=0; i!=f; i++) {
+    R[i].x = BN_new();
+    R[i].y = BN_new();
+  }
 
-  /** STEP 2: GENERATING R, THE POOL OF B-SMOOTH NUMBERS */
-  for (i=0; i < r; ) {
-    BN_pseudo_rand_range(R[i].x, rsa->n);
-    /* yᵢ = xᵢ² - N */
-    BN_sqr(R[i].y, R[i].x, ctx);
-    BN_sub(R[i].y, R[i].y, rsa->n);
-    if (dixon_smooth(R[i].y, ctx, m->M[i], r)) i++;
+  /** STEP 2 generating R */
+  for (i=0; i < m->f; i++) {
+    fprintf(stderr, "[!] Discovering %zdth smooth number\n", i);
+    discover_smooth(R[i].y, R[i].x, rsa->n,
+                    ctx, m->M[i], m->r);
   }
 
-  fprintf(stderr, "dio can\n");
-  /** STEP 3: FINDING EVEN POWERS AND ATTEMPTING FACTORIZATION **/
+  /** STEP 3: break & enter. */
   h = kernel(m);
-
-  x = BN_new(); BN_one(x);
-  sqy = BN_new(); BN_one(sqy);
-  y = BN_new();
-  rem = BN_new();
-  gcd = BN_new();
-
+  BN_one(x);
+  BN_one(sqy);
   for (i=0; i!=f; i++)
     /* if we found an even power */
     if (is_vzero(m->M[i], f)) {
@@ -173,9 +194,10 @@ dixon_question_ask_rsa(const RSA *rsa)
           BN_mul(sqy, sqy, R[j].y, ctx);
         }
       BN_sqrtmod(y, rem, sqy, ctx);
-      if (!BN_is_zero(rem)) { fprintf(stderr, "Fatal sqrt() error!\n"); exit(1);}
+      assert(!BN_is_zero(rem));
       BN_gcd(gcd, x, y, ctx);
-      if (BN_cmp(gcd, rsa->n) < 0 && BN_cmp(gcd, BN_value_one()) > 0) {
+      if (BN_cmp(gcd, rsa->n) < 0 &&
+          BN_cmp(gcd, BN_value_one()) > 0) {
         ret = RSA_new();
         ret->p = BN_dup(gcd);
         ret->q = BN_new();
@@ -184,6 +206,12 @@ dixon_question_ask_rsa(const RSA *rsa)
       }
     }
 
+  /* free all the shit */
+  for (i=0; i!=f; i++) {
+    BN_free(R[i].x);
+    BN_free(R[i].y);
+  }
+  free(R);
   BN_free(x);
   BN_free(y);
   BN_free(sqy);
@@ -191,10 +219,12 @@ dixon_question_ask_rsa(const RSA *rsa)
   BN_free(gcd);
   BN_CTX_free(ctx);
   matrix_free(m);
+
   return ret;
 }
 
 qa_question_t DixonQuestion = {
   .name = "dixon",
-  .pretty_name = "Dixon's Factorization"
+  .pretty_name = "Dixon's Factorization",
+  .ask_rsa = dixon_question_ask_rsa
 };

src/questions/include/qdixon.h

@@ -18,6 +18,9 @@ void matrix_free(matrix_t *m);
 
 matrix_t *kernel(matrix_t *m);
 
+void discover_smooth(BIGNUM *y, BIGNUM *x, BIGNUM *n,
+                     BN_CTX *ctx, char *v, size_t len);
+
 int dixon_smooth(BIGNUM *x, BN_CTX *ctx, char *v, size_t len);
 
 #endif /* _QA_DIXON_H_ */

src/questions/primes.c

@@ -68,7 +68,6 @@ smooth(BIGNUM *x, BN_CTX *ctx, char* v, size_t thresh)
   pit_t *it;
   size_t i;
 
-  i = 0;
   BN_zero(rem);
   bzero(v, thresh);
   if (BN_cmp(x, BN_value_one()) < 1) return 0;
@@ -77,8 +76,12 @@ smooth(BIGNUM *x, BN_CTX *ctx, char* v, size_t thresh)
   for (it = primes_init();
        primes_next(it, p) && i < thresh;
        i++) {
+    if (BN_cmp(p, x) > 0) break;
+
     BN_div(dv, rem, x, p, ctx);
-    for (v[i] = 0; BN_is_zero(rem); BN_div(dv, rem, x, p, ctx)) {
+    for (v[i] = 0;
+         BN_is_zero(rem);
+         BN_div(dv, rem, x, p, ctx)) {
       v[i] = (v[i] + 1) % 2;
       BN_copy(x, dv);
 

src/questions/tests/test_dixon.c

@@ -91,12 +91,30 @@ test_dixon_smooth(void)
   return;
 }
 
+void
+test_discover_smooth(void)
+{
+  BIGNUM
+    *x = BN_new(),
+    *y = BN_new(),
+    *n = BN_new();
+  BN_CTX *ctx = BN_CTX_new();
+  int len = 50;
+  char v[len + 1];
+
+  BN_dec2bn(&n, "541");
+  discover_smooth(y, x, n, ctx, v, len);
+  assert(v[0] == 1 || v[0] == 0);
+}
+
+
 int
 main(int argc, char **argv)
 {
   test_matrix();
   test_kernel();
   test_dixon_smooth();
+  test_discover_smooth();
 
   return 0;
 }