VFYPR is uses a combination of methods to verify primes: Factorization of N - 1, factorization of N + 1, and the APRCL method. The implementation of the APRCL algorithm is a more-or-less direct conversion of APRT-CLE from UBASIC to C, and somewhat extended to handle numbers of up to about 3317 digits.
An outline of the algorithm is given below.
The program runs on PC systems under MSDOS. Download it from here.
Email me if you have any comments or bug reports. Email address: tonyforbes@ltkz.demon.co.uk.
Preliminary tests
We assume N is not too small. We also assume that N is not divisible by any of the relevant parameters of the algorithm. This is achieved by trial-dividing N by anything that could cause trouble.
Let N be given and suppose d is a prime divisor of N.
Pocklington's Theorem
Suppose N - 1 = F1 R1, 2|F1, gcd(F1, R1) = 1 and F1 is completely factorized into verified primes. Suppose also that for each prime factor f of F1, there exists a number b such that 1 < b < N - 1, b^(N - 1) = 1 (mod N) and gcd(b^((N - 1)/f) - 1, N) = 1. Then, by Pocklington's Theorem [BLS1975, Theorem 4],
(1) d == 1 (mod F1).
Morrison's Theorem
Suppose N + 1 = F2 R2, 2|F2, gcd(F2, R2) = 1 and F2 is completely factorized into verified primes. Suppose also that there is an integer D such that (D/N) = -1 and for each prime factor f of F, there exists a Lucas sequence {U_k} with discriminant D = P^2 - 4Q, where gcd(N, Q) = 1, for which U_(N + 1) = 0 (mod N) and gcd(U_((N + 1)/f), N) = 1. Then, by Morrison's Theorem [BLS1975, Theorem 16],
(2) d == +/-1 (mod F2).
The APRCL method
This part follows closely Section 9.1 of Henri Cohen's book [C1995]. See also the paper of Cohen & Lenstra [CL1987].
Denote by v_p(x) the exponent of the largest power of p that divides x and by zeta_n the n-th root of unity e^(2 pi i / n). Let chi be a character modulo q. The Gauss sum is defined by
tau(chi) = sum{x in (Z/qZ)*: chi(x) (zeta_q)^x}.
Let chi and psi be characters modulo q. The Jacobi sum is defined by
j(chi, psi) = sum{x in (Z/qZ)*: chi(x) psi(1 - x)}.
Let T be an even integer. Write
e(T) = 2 product{q prime, q - 1|T: q^(v_q(T) + 1)}.
Let S be a factor of e(T) such that gcd(S, e(T)/S) = 1.
For each pair of primes q|S and p|q - 1, let k = v_p(q - 1) and let chi be the character modulo q of order p^k defined by chi(x) = (zeta_p^k)^ind_g(x), where g is a primitive root modulo q. Assume gcd(N, ST) = 1. Assume that for each pair (p, q) the character chi satisfies condition (*beta) for some beta for which (zeta_p)^beta is not equal to 1. Suppose also that for all primes p|T, condition L_p is satisfied. Then
(3) d == N^i (mod S) for some i, 0 < i < T.
Condition (*beta)
Let q, p, k and chi be as above; q|S, p^k||q - 1, and chi is a character modulo q of order p^k. Let sigma_x denote the homomorphism of Z[zeta_p^k] that sends zeta_p^k to (zeta_p^k)^x. If h = A sigma_a + B sigma_b + ... is a linear combination of such homomorphims and if f(zeta_p^k) is a polynomial in zeta_p^k, we use the expression f(zeta_p^k)^h to denote the product f((zeta_p^k)^a)^A f((zeta_p^k)^b)^B ....
The condition (*beta) is defined by
(*beta) tau^(beta (N - sigma_N)) == eta(chi)^(-beta N) (mod N)
for some eta(chi) in <zeta_p^k>.
When applying the (*beta) test it is sufficient to show that tau^(beta (N - sigma_N)) is congruent to a p^k-th root of unity modulo N. If also condition L_p is satisfied, then in fact eta(chi) = chi(N) [C1995, Proposition 9.1.11].
For a practical primality-proving algorithm we replace (*beta) by equivalent conditions involving Jacobi sums.
Suppose p > 2 (and p < 1093) [C1995, Proposition 9.1.20]. Let
beta = sum{p^k/2 < x < p^k, gcd(x, p) = 1: (sigma_x)^(-1)}.
Then (zeta_p)^beta is not equal to 1, and condition (*beta) is equivalent to the congruence
j(chi, chi)^alpha == eta(chi)^(-cN) (mod N),
where
alpha = sum{0 < x < p^k, gcd(x, p) = 1: [Nx/p^k] (sigma_x)^(-1)}
and
c = 2 (2^((p - 1) p^(k - 1)) - 1) / p^k.
Suppose p = 2 and k > 2 [C1995, Corollary 9.1.23]. Let
beta = sum{0 < x < 2^k, x == 1 or 3 (mod 8): [3x/2^k] (sigma_x)^(-1)}.
Then (zeta_p)^beta is not equal to 1, and condition (*beta) is equivalent to the congruence
(j(chi, chi) j(chi, chi^2))^alpha j(chi^(2^(k - 3)), chi^(3*2^(k - 3)))^(2 delta) == (-1)^delta eta(chi)^(-cN) (mod N),
where
alpha = sum{0 < x < 2^k, x == 1 or 3 (mod 8): [Nx/2^k] (sigma_x)^(-1)},
c = 3*(3^(2^(k - 2)) - 1)/2^k
and
delta = 0 if N == 1 or 3 (mod 8),
delta = 1 if N == 5 or 7 (mod 8).
Suppose p = 2 and k = 2. Then condition (*1) is equivalent to the congruence
j(chi, chi)^((N - 1)/2) q^((N - 1)/4) == eta(chi)^(-1) (mod N)
if N == 1 (mod 4), and to the congruence
j(chi, chi)^((N + 1)/2) q^((N - 3)/4) == - eta(chi) (mod N)
if N == 3 (mod 4).
Suppose p = 2 and k = 1. Then condition (*1) is equivalent to the congruence
q^((N - 1)/2) == eta(chi) (mod N)
if N == 1 (mod 4), and to the congruence
q^((N - 1)/2) == - eta(chi) (mod N)
if N == 3 (mod 4).
Condition L_p
For prime p we say that condition L_p is satisfied if for all prime divisors r of N and all integers a > 0 we can find an integer l_p(r, a) such that
r^(p - 1) == N^((p - 1) l_p(r, a)) (mod p^a).
Assume that we can find a character chi modulo q of order p^k and a beta, with (zeta_p)^beta not equal to 1, for which (*beta) is satisfied. Then condition L_p is satisfied if the following additional condition is true:
p > 2 and N^(p - 1) is not congruent to 1 (mod p^2).
However, if (*beta) is satisfied and if eta(chi) is a primitive p^k-th root of unity, then condition L_p is satisfied if one of the following additional conditions is true [C1995, Proposition 9.1.17]:
(i) p > 2
(ii) p = 2, k = 1 and N == 1 (mod 4)
(iii) p = 2, k > 1 and q^((N - 1)/2) == -1 (mod N).
Checking the divisors d of N
Assume gcd(S, F1 F2) = 1, Let G = F1 F2 S/2. We require G^3 > N.
From the conditions (1), (2) and (3) above, by the Chinese Remainder Theorem, d must belong to one of at most 2T residue classes mod G. We compute the representatives of these residue classes which are less than G and try them one by one as possible non-trivial divisors of N. If there are none, we can conclude that d > G.
If G^2 > N, N is prime.
Otherwise, let F be the greater of F1 and F2. Let H = G F^2. In addition to G^3 > N, the final condition we require is that mH > N, where m is not too large. (m < 1000000000, say.) For then we can complete a primality proof by means of one of the following theorems.
Theorem 1. Suppose N is composite, suppose N - 1 = FR, where F is even, R is odd, R > 1 and gcd(F, R) = 1. Suppose also that for each prime factor f of F, there is an integer b, 1 < b < N-1, such that b^(N - 1) == 1 (mod N) and gcd(b^((N - 1)/f) - 1, N) = 1. If every prime factor of N is greater than G, then N is the product of two primes, cF + 1 and dF + 1. Furthermore, if
N < G (2aF^2 - G + rF + 2)
for some a > 0, where r and s are defined by R = 2Fs + r, 0 < r < 2F, then 2(s - a) < cd <= 2s.
Theorem 2. Suppose N is composite, suppose N + 1 = FR, where F is even, R is odd, R > 1 and gcd(F, R) = 1. Suppose there exists an integer D for which (D/N) = -1 and also for each prime factor f of F a Lucas sequence {U_i} with discriminant D such that N divides U_(N+1) and gcd(U_((N+1)/f), N) = 1. If every prime factor of N is greater than G, then N is the product of two primes, cF + e and dF - e, where e = 1 or -1. Furthermore, if
N < G (2aF ^2 + G - |rF - 2|)
for some a > 0, where r and s are defined by R = 2Fs + r, |r| < F, then 2(s - a) < cd <= 2(s + a).
References
[BLS1975] J. Brillhart, D.H. Lehmer and J.L. Selfridge, 'New Primality Criteria and Factorizations of 2^m +/- 1', Mathematics of Computation, 29 (1975), 620-647.
[CL1987] H. Cohen & A. K. Lenstra, 'Implementation of a new primality test', Math. Comp., 48 (1987), 103-121.
[C1995] Henri Cohen, A Course in Computational Algebraic Number Theory, Springer-Verlag, Berlin, 1995.