3 * Calculate discrete logs in prime-order groups
5 * (c) 2017 Mark Wooding
8 /*----- Licensing notice --------------------------------------------------*
10 * This file is part of Rhodes, a distributed discrete-log finder.
12 * Rhodes is free software; you can redistribute it and/or modify
13 * it under the terms of the GNU General Public License as published by
14 * the Free Software Foundation; either version 2 of the License, or
15 * (at your option) any later version.
17 * Rhodes is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 * GNU General Public License for more details.
22 * You should have received a copy of the GNU General Public License
23 * along with Rhodes; if not, write to the Free Software Foundation,
24 * Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
41 #include <sys/types.h>
48 #include <NTL/GF2XFactoring.h>
52 static const char *prog
;
53 static int stdin_fdflags
= -1;
55 static void cleanup(int sig
)
57 if (stdin_fdflags
>= 0) fcntl(0, F_SETFL
, stdin_fdflags
);
58 if (sig
== SIGTSTP
) raise(SIGSTOP
);
59 else if (sig
) { signal(sig
, SIG_DFL
); raise(sig
); }
62 __attribute__((noreturn
))
63 static void barf(const char *msg
, int err
, ...)
68 std
::fprintf(stderr
, "%s: ", prog
);
69 std
::vfprintf(stderr
, msg
, ap
);
70 if (err
) std
::fprintf(stderr
, ": %s", std
::strerror(err
));
71 std
::fputc('\n', stderr
);
77 static void wakeup(int sig
)
79 if (fcntl(0, F_SETFL
, stdin_fdflags
| O_NONBLOCK
))
80 barf("fcntl(stdin)", errno
);
83 static int intarg(const char *p
, const char *what
)
90 n
= std
::strtol(p
, &q
, 0);
91 if (*q
|| errno
|| n
< 0 || n
> INT_MAX
)
92 barf("bad int %s `%s'", 0, what
, p
);
97 static NTL
::GF2X
polyarg(const char *p
, const char *what
)
101 // Oh, this is wretched: NTL reads and writes the nibbles backwards. Hack
103 if (s
.size() < 2 || s
[0] != '0' || s
[1] != 'x') goto bad
;
104 std
::reverse(s
.begin() + 2, s
.end());
107 std
::istringstream ss
{s
};
113 barf("bad poly %s `%s'", 0, what
, p
);
116 static unsigned long long hash(const unsigned char *p
, size_t n
)
119 unsigned long long h
= 0x6a078c523ae42f68ull
;
121 for (i
= 0; i
< n
; i
++) { h
+= p
[i
]; h
*= 0xbeaa913866b8d5a3ull
; }
125 static std
::string
showpoly(NTL
::GF2X f
)
127 std
::ostringstream ss
;
129 std
::string s
{ss
.str()};
130 assert(s
.size() >= 2 && s
[0] == '0' && s
[1] == 'x');
131 std
::reverse(s
.begin() + 2, s
.end());
135 static NTL
::ZZ
zzarg(const char *p
, const char *what
)
138 std
::istringstream ss
{s
};
141 if (!ss
) barf("bad int %s `%s'", 0, what
, p
);
147 unsigned char b
[256];
150 if ((fp
= std
::fopen("/dev/urandom", "rb")) == 0)
151 barf("open(/dev/urandom)", errno
);
152 if (std
::fread(b
, 1, sizeof(b
), fp
) != sizeof(b
)) {
153 barf("read(/dev/urandom)%s",
154 ferror(fp
) ? errno
: 0,
155 ferror(fp
) ?
"" : ": unexpected short read");
158 NTL
::ZZ t
{NTL
::ZZFromBytes(b
, sizeof(b
))};
168 unsigned long long k
;
177 #define CHECK_NITER 65536
179 int main(int argc
, char *argv
[])
183 unsigned char buf
[4096];
184 unsigned long long niter
, dpmask
;
191 NTL
::GF2X
::HexOutput
= 1;
193 // Collect the arguments and check that they make sense.
195 std
::fprintf(stderr
, "usage: %s DPBITS gf2x P A B L\n",
200 // Distinguished-point bits, or zero to find a full cycle.
201 dpbits
= intarg(argv
[1], "dpbits");
202 dpmask
= (1ull << dpbits
) - 1;
204 // Group specification. Currently only subgroups of (GF(2)[x]/(p(x)))^*
206 if (std
::strcmp(argv
[2], "gf2x") != 0)
207 barf("unknown group kind `%s'", 0, argv
[1]);
208 NTL
::GF2X p
= polyarg(argv
[3], "p");
209 if (!IterIrredTest(p
)) barf("not irreducible", 0);
212 // The base and operand for the logarithm.
213 NTL
::GF2E a
{to_GF2E(polyarg(argv
[4], "a"))};
214 if (a
== 0 || a
== 1) barf("a is trivial", 0);
215 NTL
::GF2E b
{to_GF2E(polyarg(argv
[5], "b"))};
216 if (b
== 0) barf("b isn't a unit", 0);
218 // The group order, which we expect to have been calculated. This must be
219 // prime, or we risk failing later.
220 NTL
::ZZ l
{zzarg(argv
[6], "l")};
221 if (!ProbPrime(l
)) barf("order isn't prime", 0);
224 // Check the points at least have the same order.
225 if (power(a
, l
) != 1) barf("a has wrong order", 0);
226 if (power(b
, l
) != 1) barf("b has wrong order", 0);
228 // Build the table of steps. Must do this deterministically! The basic
229 // structure of the walk is taken from Edlyn Teske, 1998, `Better Random
230 // Walks for Pollard's Rho Method': about twenty random steps, together
231 // with a couple of slots' worth of squaring. Deviating slightly from
232 // Teske, I'm using a prime number of slots to improve the hash
234 step tab
[NSTEP
- NSQR
];
235 SetSeed(NTL
::ZZ
::zero());
236 for (i
= 0; i
< NSTEP
- NSQR
; i
++) {
237 tab
[i
].du
= NTL
::random_ZZ_p();
238 tab
[i
].dv
= NTL
::random_ZZ_p();
239 tab
[i
].m
= power(a
, rep(tab
[i
].du
))*power(b
, rep(tab
[i
].dv
));
242 // Prepare signal handlers. We're going to make stdin be non-blocking,
243 // which will have a persistent effect on whatever file is being used
244 // (e.g., a terminal), so we should be careful to undo it when we exit.
245 stdin_fdflags
= fcntl(0, F_GETFL
);
246 if (stdin_fdflags
< 0) barf("fcntl stdin", errno
);
247 sa
.sa_handler
= cleanup
;
248 sigemptyset(&sa
.sa_mask
);
249 sigaddset(&sa
.sa_mask
, SIGTSTP
);
250 sigaddset(&sa
.sa_mask
, SIGCONT
);
252 if (sigaction(SIGINT
, &sa
, 0) ||
253 sigaction(SIGTERM
, &sa
, 0) ||
254 sigaction(SIGQUIT
, &sa
, 0) ||
255 sigaction(SIGTSTP
, &sa
, 0))
256 barf("sigaction", errno
);
257 sa
.sa_handler
= wakeup
;
258 if (sigaction(SIGCONT
, &sa
, 0))
259 barf("sigaction", errno
);
260 if (fcntl(0, F_SETFL
, stdin_fdflags
| O_NONBLOCK
))
261 barf("fcntl stdin", errno
);
263 // Now we start the random walk...
265 niter
= 8ull << (dpbits ? dpbits
: (NumBits(l
) + 1)/2);
267 NTL
::ZZ_p u
{NTL
::random_ZZ_p()}, v
{NTL
::random_ZZ_p()};
268 NTL
::GF2E t
= power(a
, rep(u
))*power(b
, rep(v
));
272 unsigned long long k
= 0;
274 // Initialize the table of previous steps in the walk, if we're going to
275 // finding a cycle. Deviating from Teske's algorithm, initialize the table
276 // with zeros. If we don't do this then, in the case where b = 1, the
277 // algorithm tends to find the cycle but can't deduce anything useful from
280 for (i
= 0; i
< NHIST
; i
++)
281 { h
[i
].k
= 0; h
[i
].y
= a
; h
[i
].u
= 1; h
[i
].v
= 0; }
284 if (k
>= niter
) goto again
;
286 // Check for input on stdin. If stdin has closed, that means that our
287 // caller has gone away, presumably because they no longer care about
288 // what we have to compute. (Maybe some other job found the answer
289 // quicker than we did.) If stdin has data in, then read it: this will
290 // keep a network connection alive (or notice that it's gone dead).
291 if (!(k
%CHECK_NITER
)) {
292 FD_ZERO(&fds_in
); FD_SET(0, &fds_in
);
293 tv
.tv_sec
= 0; tv
.tv_usec
= 0;
294 if (select(1, &fds_in
, 0, 0, &tv
) < 0)
295 { if (errno
!= EINTR
) barf("select", errno
); }
296 else if (FD_ISSET(0, &fds_in
)) {
298 n
= read(0, buf
, sizeof(buf
));
299 if (!n
) { cleanup(0); std
::exit(1); }
300 else if (n
> 0) continue;
301 else if (errno
== EAGAIN
) break;
302 else if (errno
!= EINTR
) barf("read stdin", errno
);
308 size_t nb
= NumBytes(rep(t
));
309 BytesFromGF2X(buf
, rep(t
), nb
);
310 unsigned long long hh
= hash(buf
, nb
);
313 // If we're walking to a distinguished point (the algorithm of Paul van
314 // Oorschot and Michael Wiener, 1994, `Parallel Collision Search with
315 // Application to Hash Functions and Discrete Logarithms') then check
316 // the hash to see if enough low-order bits are zero.
318 std
::cout
<< u
<< " " << v
<< " " << showpoly(rep(t
)) << " "
323 // Otherwise, check for a cycle. We use the algorithm is from Edlyn
324 // Teske, 1998, `A Space Efficient Algorithm for Group Structure
325 // Computation', rather than the usual one by Floyd, which requires
326 // about three times as much computation.
328 // Check for a match.
329 for (i
= 0; i
< NHIST
; i
++) {
331 if (h
[i
].u
== u
|| h
[i
].v
== v
) goto again
;
332 std
::cout
<< (h
[i
].u
- u
)/(v
- h
[i
].v
) << " " << k
<< std
::endl
;
337 // Update the table of earlier points if necessary.
339 h
[o
].u
= u
; h
[o
].v
= v
; h
[o
].y
= t
; h
[o
].k
= k
;
344 // Take a step in the random walk.
346 if (i
>= NSTEP
- NSQR
) { mul(u
, u
, 2); mul(v
, v
, 2); sqr(t
, t
);; }
347 else { add(u
, u
, tab
[i
].du
); add(v
, v
, tab
[i
].dv
); mul(t
, t
, tab
[i
].m
); }