|
发表于 2023-9-25 16:11:34
|
显示全部楼层
EQFT的简化实现在附件的论文里,
称为\(Simplified \: Quadratic \: Frobenius \: test\),简称\(SQFT\)
下面算法中的\( \phi_3(x)=x^2+x+1, \phi_4(x)=x^2+1, \phi_8(x)=x^4+1\)
一、前置算法 \(MR2 (Miller-Rabin \: test \: with \: basis \: two \: or \: a \: small \: nonresidue)\)
输入: 奇数\(n\)
输出: \(n\)为合数,或者
\( c, (\frac{c}{n}) = -1,\epsilon \in R(n, c), \epsilon^4 = -1, \phi_8(\epsilon) = 0 \)
1、\( n \mod 4 = 3 \)
计算 \( \alpha = 2^{(n-3)/4} (\mod n) \)
如果\( 2\alpha^2 \neq \pm 1 (\mod n) \) 输出\(n\)是合数
否则,输出\(c = -1, \epsilon = \alpha + \alpha x \)
2、\( n \mod 8 = 5 \)
计算 \( \alpha = 2^{(n-1)/4} (\mod n) \)
如果\( \alpha^2 \neq - 1 (\mod n) \) 输出\(n\)是合数
否则,输出\(c = 2, \epsilon = \frac{1 + \alpha}{2} x \)
(这里的分数,在实际上的模 n 环 $R_n$ 中的计算可以转换成整数,即如果$\alpha$是偶数,
$\frac{1+\alpha}{2} \equiv \frac{1+\alpha+n}{2} (\mod n)$)
3、\( n \mod 8 = 1 \)
如果\(n\)是完全平方数,输出\(n\)是合数
否则找到小的\(c\)满足\( (\frac{c}{n}) = - 1 \) //c从3开始
计算 \( \alpha = c^{(n-1)/8} (\mod n) \)
如果\( \alpha^4 \neq - 1 (\mod n) \) 输出\(n\)是合数
否则,输出\(c = -1, \epsilon = \alpha \)
二、SQFT循环 \(SQFT_{round}\)
输入: 奇数\(n\),\( c, (\frac{c}{n}) = -1,\epsilon \in R(n, c), \epsilon^4 = -1, \phi_8(\epsilon) = 0 \)
输出:\(n\)是合数或者可能是素数
1、随机选择\( z \in R(n, c) \)满足\( (\frac{N(z)}{n}) = -1 \)
无心人备注:这个满足条件可以去掉,考虑1加上小素数的集合,\( z \in R(n, c) \) 表达成\(a + bx \),那么\( a, b \in \{ q | q = 1 \quad or \quad q \in P, q < T\} \) 这里\(T\)可以取个固定上限比如100,则满足26*25=650次测试需求了
2、如果\( z^n \neq \overline{z} \) 输出\(n\)是合数
3、如果\( z^{(n^2-1)/8} \notin \{ \pm \epsilon, \pm \epsilon^3 \} \) 输出\(n\)是合数
无心人备注:经过用29实际测试,发现这里应该是 \( \{ \pm 1, \pm \epsilon, \pm \epsilon^2, \pm \epsilon^3 \} \)
4、输出\(n\)可能是素数
三、SQFT测试 \(SQFT: (Simplified \: Quadratic \: Frobenius \: test) \)
输入: \(n, n > 20 \), 测试次数\(t\)
输出:\(n\)是合数或者可能是素数
1、如果\(n\)被小于\(20\)素数整除,输出\(n\)是合数
2、调用算法\(MR2\)
3、如果\(n\)没被判定为合数,则\(MR2\)输出\(c, \epsilon \)
4、重复\(t\)次:
用\(n, c, \epsilon \)调用算法\( SQFT_{round} \),如果算法判定\(n\)是合数,则终止
5、输出\(n\)可能是素数
附加三次根测试的\(SQFT3\)算法
四、SQFT3循环 \( SQFT3_{round} \)
输入:整数\( n, gcd(n, 6) = 1 \), 小整数\( c, (\frac{c}{n} ) = -1 \),
值\( \epsilon, \epsilon^4 = -1 \)和\( \epsilon_3 \)满足 \( \epsilon_3 =1 \) 或者 \( \phi_3(\epsilon_3) = 0\),
输出:\( n \)是合数,或者是可能的素数,同时输出下面的值
\( {\epsilon'}_3 \), 满足 \( {\epsilon'}_3 =1 \) 或者 \( \phi_3({\epsilon'}_3) = 0\),
如果, \( \epsilon_3 \neq 1\),那么 \({\epsilon'}_3 = \epsilon_3^{\pm 1} \)
1、随机选择\( z \in R(n, c) \)满足\( (\frac{N(z)}{n}) = -1 \)
无心人备注:这里应该跟上面一样,忽略掉对z的范的Jacobi符号要求
2、如果\( z^n \neq \overline{z} \) 输出\(n\)是合数
3、如果\( z^{(n^2-1)/8} \notin \{ \pm \epsilon, \pm \epsilon^3 \} \) 输出\(n\)是合数
无心人备注:这里跟上面的SQFTRound一致即可
4、置\( u = v_3(n^2 - 1), n^2 - 1 = 3^u r \) //这里的\(v_3\)看不出什么意思来~
5、令\( i = \min\{ j: 0 \leq j \leq u, z^{3^j r} = 1 \} \)
6、如果\( i = 0\),输出\(n\)可能是素数,\( {\epsilon'}_3 = \epsilon_3 \)
7、置\( {\epsilon'}_3 = z^{3^{i-1}r} \)
8、如果\( \epsilon_3 = 1 \) 且 \( \phi_3({\epsilon'}_3) \neq 0 \),输出n是合数
9、如果\( \epsilon_3 \neq 1 \) 且\( {\epsilon'}_3 \neq \epsilon_3^{\pm 1} \),输出n是合数
无心人备注:这里猜测和SQFTRound一样,是符合三次根,SQFTRound是八个八次根
10、输出\(n\)可能是素数和 \({\epsilon'}_3\)
五、SQFT3测试 \( SQFT3 \)
输入: \(n, n > 200 \), 测试次数\(t\)
输出:\(n\)是合数或者可能是素数
1、如果 \(n\)被小于\(200\)素数整除,输出\(n\) 是合数
2、调用\(MR2\),如果判定 \(n\)是合数,结束
3、从\(MR2\)获得输出的\( c, \epsilon \),并且置 \( \epsilon_3 = 1 \)
4、循环\(t\)次
用\( n, c, \epsilon, \epsilon_3 \)调用\( SQFT3_{round} \)
如果\( SQFT3_{round} \)判定\(n\)是合数,结束
置\( \epsilon_3 = {\epsilon'}_3 \),这里 \( {\epsilon'}_3 \)是算法\( SQFT3_{round} \)的输出
5、输出\(n\)可能是素数
乘法实现:
\(w, z \in R(n, c), w = A_wx + B_w, z = A_zx + B_z \)
\(z · w = (m_1 + m_2)x + (cA_z + B_z)(A_w + B_w) − cm_1 − m_2 \)
\( m_1 = A_zB_w , m_2 = B_zA_w \)
这里其实是带入\(x=\sqrt{c}\)来计算的
Julia源代码,不包含SQFT3实现- using Heckefunction qfnegmod(F, x, n) re = mod(ZZ(-coeff(x, 0)), n) ir = mod(ZZ(-coeff(x, 1)), n) return F([re, ir])endfunction qfconjmod(F, x, n) re = coeff(x, 0) ir = mod(ZZ(-coeff(x, 1)), n) return F([re, ir])endfunction qfmulmod(F, x, y, n) z = x * y re = mod(ZZ(coeff(z, 0)), n) ir = mod(ZZ(coeff(z, 1)), n) return return F([re, ir])endfunction qfpowermod(F, x, p , n) @assert p >= 0 p == 0 && return F([1, 0]) b = x t = ZZ(prevpow(BigInt(2), BigInt(p))) r = F([1, 0]) while true if p >= t r = qfmulmod(F, r, b, n) p -= t end t >>= 1 t <= 0 && break r = qfmulmod(F, r, r ,n) end return rendfunction makeEpsilonList(F, e, n) list = [] one = F([1, 0]) append!(list, [one]) append!(list, [qfnegmod(F, one, n)]) append!(list, [e]) append!(list, [qfnegmod(F, e, n)]) e2 = qfmulmod(F, e, e, n) append!(list, [e2]) append!(list, [qfnegmod(F, e2, n)]) e3 = qfmulmod(F, e2, e, n) append!(list,[e3]) append!(list, [qfnegmod(F, e3, n)]) return listendfunction MR2(n) if n % 4 == 3 alpha = powermod(ZZ(2), (n-3)>>2, n) tmp = (2*alpha*alpha) % n (tmp != 1) && (tmp != (n-1)) && return (false, 0, []) F, t = quadratic_field(-1) return (true, F, makeEpsilonList(F, F([alpha, alpha]), n)) end if n % 8 == 5 alpha = powermod(ZZ(2), (n-1)>>2, n) tmp = (alpha*alpha) % n (tmp != (n-1)) && return (false, 0, []) iseven(alpha) && (alpha+=n) alpha = (alpha+1)>>1 F, t = quadratic_field(2) return (true, F, makeEpsilonList(F, F([0, alpha]), n)) end if n % 8 == 1 issquare(n) && return (false, 0, []) c = ZZ(3) while jacobi_symbol(c, n) != -1 c = next_prime(c) end alpha = powermod(c, (n-1)>>3, n) tmp = powermod(alpha, 4, n) (tmp != (n-1)) && return (false, 0, []) F, t = quadratic_field(c) return (true, F, makeEpsilonList(F, F([alpha, 0]), n)) endendfunction SQFTRound(F, x, n, list) r = qfpowermod(F, x, n, n) tmp = qfconjmod(F, x, n) #println("1:$n: x:$x, r:$r, conj(x):$tmp") r != tmp && return false tmp = (n^2-1) >> 3 (td, tr) = divrem(tmp, n) tmp1 = qfpowermod(F, r, td, n) tmp2 = qfpowermod(F, x, tr, n) r = qfmulmod(F, tmp1, tmp2, n) #println("2: $n: $x, $r") return in(r, list)endfunction SQFT(n, t) Prm = [2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97] for p in Prm (n % p == 0) && return (n == p) end Q = append!([1], Prm) (r, F, list) = MR2(n) if r #println("MR2 True at $n: $F") #println("E: $list") i = 0 for x in Q for y in Q if x != y if SQFTRound(F, F([x, y]), n, list) i = i + 1 (i > t) && return true else return false end end end end else return false endend########################################################global b = ZZ(10)^67@time for i in range(3, step = 2, stop = 19999) test = b + i isprime(test) && !SQFT(test, 1) && println("frobenius error at $test") SQFT(test, 1) && !isprime(test) && println("frobenius error at $test")endprintln("end")
复制代码
对于可能是奇数的平方数的情况
如果n不能被3,5整除
n模240必然是[1, 49, 121, 169]之一
|
|