无心人
发表于 2019-9-28 09:53:55
继续刷楼
mathematica
发表于 2019-9-29 09:19:14
无心人 发表于 2019-9-22 09:43
基于Jacobi和的分圆域素性检验算法
我觉得你可真够无聊的,其实
baillie psw算法已经够好的了,
https://bbs.emath.ac.cn/forum.php?mod=viewthread&tid=15776&fromuid=865
这儿有现成的代码!
无心人
发表于 2019-9-29 11:21:28
mathematica 发表于 2019-9-29 09:19
我觉得你可真够无聊的,其实
baillie psw算法已经够好的了,
https://bbs.emath.ac.cn/forum.php?mod=vie ...
你可能不知道BPSW的算法里,除了Miller-Rabin的另一个算法,是跟我说的Frobenius类似的算法,你觉得,对它的优化,值不值得去讨论?
你实现那个,是用F-L序列实现的,而对F-L改进的算法是我说的基于二次域的\(1+\sqrt{c}\) 的,实现起来简单多了
mathematica
发表于 2019-10-1 09:47:07
无心人 发表于 2019-9-29 11:21
你可能不知道BPSW的算法里,除了Miller-Rabin的另一个算法,是跟我说的Frobenius类似的算法,你觉得, ...
你说的我不懂,但是我知道bpsw目前来说还没发现反例,真的很优秀!
mathematica
发表于 2019-10-1 09:56:29
无心人 发表于 2019-9-29 11:21
你可能不知道BPSW的算法里,除了Miller-Rabin的另一个算法,是跟我说的Frobenius类似的算法,你觉得, ...
一个2为底的miller rabin+一个lucas动参数的测试,
两者联合起来都已经很难找到反例了,
如果多几个miller rabin测试+几个lucas动参数,更难找到反例。
两者的因数分解的模的形式不一样,所以壬就考虑联合起来,
但是没想到很难找到反例!
无心人
发表于 2019-10-6 10:14:25
mathematica 发表于 2019-10-1 09:56
一个2为底的miller rabin+一个lucas动参数的测试,
两者联合起来都已经很难找到反例了,
如果多几个mill ...
那个二次形式的Frobenius测试,就是你说的lucas在二次域的变形,原理是一样的,你可以看我给的wiki链接,但是,有个重要区别是,整数递归形式的测试,是存在很多伪素数例子的,而二次域的形式,10^8内我没找到伪素数,你可以说整数形式已经足够了,没有反例,但是这只是因为两个测试的伪素数没有发现叠加部分而已,而二次域形式的直接来自二次域的Fermat小定理,有论文分析过出现伪素数的可能,条件要远远的比整数上的fermat小定理的伪素数形式严格,甚至只要排除n是完全平方,二次域形式的Frobenius是可以独立做测试的
无心人
发表于 2019-10-7 09:46:53
素性测试的基于GMP的C源代码看大数分解程序YAFU的源代码
解压后
头文件在./include/mpz_prp_prime.h
源文件在./top/mpz_prp_prime.c
这是从头文件里扒下来的它实现的测试名字
/* ******************************************************************
* mpz_prp: (also called a Fermat pseudoprime)
* A "pseudoprime" to the base a is a composite number n such that,
* (a,n)=1 and a^(n-1) = 1 mod n
* ******************************************************************/
int mpz_prp(mpz_t n, mpz_t a);
/* *************************************************************************
* mpz_euler_prp: (also called a Solovay-Strassen pseudoprime)
* An "Euler pseudoprime" to the base a is an odd composite number n with,
* (a,n)=1 such that a^((n-1)/2)=(a/n) mod n [(a/n) is the Jacobi symbol]
* *************************************************************************/
int mpz_euler_prp(mpz_t n, mpz_t a);
/* *********************************************************************************************
* mpz_sprp: (also called a Miller-Rabin pseudoprime)
* A "strong pseudoprime" to the base a is an odd composite n = (2^r)*s+1 with s odd such that
* either a^s == 1 mod n, or a^((2^t)*s) == -1 mod n, for some integer t, with 0 <= t < r.
* *********************************************************************************************/
int mpz_sprp(mpz_t n, mpz_t a);
/* *************************************************************************
* mpz_fibonacci_prp:
* A "Fibonacci pseudoprime" with parameters (P,Q), P > 0, Q=+/-1, is a
* composite n for which V_n == P mod n
*
* *************************************************************************/
int mpz_fibonacci_prp(mpz_t n, long int p, long int q);
/* *******************************************************************************
* mpz_lucas_prp:
* A "Lucas pseudoprime" with parameters (P,Q) is a composite n with D=P^2-4Q,
* (n,2QD)=1 such that U_(n-(D/n)) == 0 mod n [(D/n) is the Jacobi symbol]
* *******************************************************************************/
int mpz_lucas_prp(mpz_t n, long int p, long int q);
/* *********************************************************************************************
* mpz_stronglucas_prp:
* A "strong Lucas pseudoprime" with parameters (P,Q) is a composite n = (2^r)*s+(D/n), where
* s is odd, D=P^2-4Q, and (n,2QD)=1 such that either U_s == 0 mod n or V_((2^t)*s) == 0 mod n
* for some t, 0 <= t < r. [(D/n) is the Jacobi symbol]
* *********************************************************************************************/
int mpz_stronglucas_prp(mpz_t n, long int p, long int q);
/* *******************************************************************************************
* mpz_extrastronglucas_prp:
* Let U_n = LucasU(p,1), V_n = LucasV(p,1), and D=p^2-4.
* An "extra strong Lucas pseudoprime" to the base p is a composite n = (2^r)*s+(D/n), where
* s is odd and (n,2D)=1, such that either U_s == 0 mod n or V_s == +/-2 mod n, or
* V_((2^t)*s) == 0 mod n for some t with 0 <= t < r-1 [(D/n) is the Jacobi symbol]
* *******************************************************************************************/
int mpz_extrastronglucas_prp(mpz_t n, long int p);
/* ***********************************************************************************************
* mpz_selfridge_prp:
* A "Lucas-Selfridge pseudoprime" n is a "Lucas pseudoprime" using Selfridge parameters of:
* Find the first element D in the sequence {5, -7, 9, -11, 13, ...} such that Jacobi(D,n) = -1
* Then use P=1 and Q=(1-D)/4 in the Lucas pseudoprime test.
* Make sure n is not a perfect square, otherwise the search for D will only stop when D=n.
* ***********************************************************************************************/
int mpz_selfridge_prp(mpz_t n);
/* *********************************************************************************************************
* mpz_strongselfridge_prp:
* A "strong Lucas-Selfridge pseudoprime" n is a "strong Lucas pseudoprime" using Selfridge parameters of:
* Find the first element D in the sequence {5, -7, 9, -11, 13, ...} such that Jacobi(D,n) = -1
* Then use P=1 and Q=(1-D)/4 in the strong Lucase pseudoprime test.
* Make sure n is not a perfect square, otherwise the search for D will only stop when D=n.
* **********************************************************************************************************/
int mpz_strongselfridge_prp(mpz_t n);
/* **********************************************************************************
* mpz_bpsw_prp:
* A "Baillie-Pomerance-Selfridge-Wagstaff pseudoprime" is a composite n such that
* n is a strong pseudoprime to the base 2 and
* n is a Lucas pseudoprime using the Selfridge parameters.
* **********************************************************************************/
int mpz_bpsw_prp(mpz_t n);
/* ****************************************************************************************
* mpz_strongbpsw_prp:
* A "strong Baillie-Pomerance-Selfridge-Wagstaff pseudoprime" is a composite n such that
* n is a strong pseudoprime to the base 2 and
* n is a strong Lucas pseudoprime using the Selfridge parameters.
* ****************************************************************************************/
int mpz_strongbpsw_prp(mpz_t n);
/* **********************************************************************************
* APR-CL (also known as APRT-CLE) is a prime proving algorithm developed by:
* L. Adleman, C. Pomerance, R. Rumely, H. Cohen, and A. K. Lenstra
* APRT-CLE = Adleman-Pomerance-Rumely Test Cohen-Lenstra Extended version
* You can find all the details of this implementation in the Cohen & Lenstra paper:
* H. Cohen and A. K. Lenstra, "Implementation of a new primality test",
* Math. Comp., 48 (1987) 103--121
*
* This C/GMP version is a conversion of Dario Alpern's Java based APRT-CLE code
* His code was based on Yuji Kida's UBASIC code
*
* Based on APRT-CLE Written by Dario Alejandro Alpern (Buenos Aires - Argentina)
* From: Last updated September 10th, 2011. See http://www.alpertron.com.ar/ECM.HTM
*
* With improvements based on Jason Moxham's APRCL v1.15 code, from 2003/01/01
* *********************************************************************************/
int mpz_aprcl(mpz_t N);
int mpz_aprtcle(mpz_t N, int verbose);
mathematica
发表于 2019-10-9 09:36:01
无心人 发表于 2019-10-7 09:46
素性测试的基于GMP的C源代码看大数分解程序YAFU的源代码
解压后
头文件在./include/mpz_prp_prime.h
在我看来,baillie psw算法足够了,不自信的话,就再加几个miller随机的测试,不纠结!
无心人
发表于 2021-11-14 18:43:56
更新了SQFT测试的源代码,打算用Julia的Hecke包重写所有代码,Hecko包会包装GMP的一些函数
无心人
发表于 2021-12-19 15:04:39
补充下 BPSW 测试
### (一)、Fibonacci 序列素性测试
定义 Fibonacci 序列,$F_0= 0, F_1 = F_2 = 1,F_{k+2} = F_{k+1} + F_k$,素数 $n > 5$, jacobi 符号 $J(\frac{n}{5})=j$,则 $F_{n-j} \equiv 0 (mod\ \ n)$。
利用上面的结果,假如 $n$ 是待测试整数,
1、若 $5 | n$ 则返回 $n$ 是合数。
2、否则若 $n \equiv \pm 1 (mod\ \ 5)$ 则 $j = 1$,若 $n \equiv \pm 2 (mod\ \ 5)$ 则 $j = -1$。
3、若 $F_{n-j} \equiv 0 (mod\ \ n)$,则返回 $n$ 可能是素数,否则返回 $n$ 是合数。
若合数 $n$ 通过上述测试,称为 Fibonacci 伪素数,简称 fpp,fpp 是无限多的。
julia示例代码
using Hecke
function fibonacciTest(n)
r = n % 5
r == 0 && return false
if r == 1 || r == 4
j = 1
end
if r == 2 || r == 3
j = -1
end
return fibonacci(n-j) % n == 0
end
for i in range(ZZ(7), stop=ZZ(9999), step=2)
if fibonacciTest(i)
if !isprime(i)
println("Fibonacci Pseudoprime: $i")
end
end
if isprime(i)
if !fibonacciTest(i)
println("Fibonacci Primality test error at: $i")
end
end
end
输出
Fibonacci Pseudoprime: 323
Fibonacci Pseudoprime: 377
Fibonacci Pseudoprime: 1891
Fibonacci Pseudoprime: 3827
Fibonacci Pseudoprime: 4181
Fibonacci Pseudoprime: 5777
Fibonacci Pseudoprime: 6601
Fibonacci Pseudoprime: 6721
Fibonacci Pseudoprime: 8149
对 10000 内大于 5 奇数测试表明,存在 9 个奇伪素数。如果不限定 $n$ 是奇数,则相应存在偶伪素数,最小的一个是 8539786。
另外,该代码是示例性质代码,并没考虑优化,序列的计算没有在模 $n$ 环上进行,当 $n$ 很大时候,可能会内存溢出、计算时间超长等。