aboutsummaryrefslogtreecommitdiffstats
path: root/src/zptest.c
blob: a1d112c09fc42fd28c2d3a56d34f7c85f697a80a (plain) (blame)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
/* See LICENSE file for copyright and license details. */
#include "internals.h"

#define x   libzahl_tmp_ptest_x
#define a   libzahl_tmp_ptest_a
#define d   libzahl_tmp_ptest_d
#define n1  libzahl_tmp_ptest_n1
#define n4  libzahl_tmp_ptest_n4


enum zprimality
zptest(z_t witness, z_t n, int t)
{
	/*
	 * Miller–Rabin primarlity test.
	 */

	size_t i, r;

	if (zcmpu(n, 3) <= 0) {
		if (zcmpu(n, 1) <= 0) {
			if (witness)
				SET(witness, n);
			return NONPRIME;
		} else {
			return PRIME;
		}
	}
	if (zeven(n)) {
		if (witness)
			SET(witness, n);
		return NONPRIME;
	}

	zsub_unsigned(n1, n, libzahl_const_1);
	zsub_unsigned(n4, n, libzahl_const_4);

	r = zlsb(n1);
	zrsh(d, n1, r);

	while (t--) {
		zrand(a, FAST_RANDOM, UNIFORM, n4);
		zadd_unsigned(a, a, libzahl_const_2);
		zmodpow(x, a, d, n);

		if (!zcmp(x, libzahl_const_1) || !zcmp(x, n1))
			continue;

		for (i = 1; i < r; i++) {
			zsqr(x, x);
			zmod(x, x, n);
			if (!zcmp(x, libzahl_const_1)) {
				if (witness)
					zswap(witness, a);
				return NONPRIME;
			}
			if (!zcmp(x, n1))
				break;
		}
		if (i == r) {
			if (witness)
				zswap(witness, a);
			return NONPRIME;
		}
	}

	return PROBABLY_PRIME;
}