素性测试中的平方数判定问题

无心人

不论是BSPW中的Lucas测试，还是我发的Quadratic Frobenius测试，都要求测试整数n不能是完全平方数，否则用Jacobi symbol找不到参数

我想到个办法，就是用n mod m余数的方法，当选择特定m时候，余数是很少的
比如m = 120, 当n不含有2，3，5素因子时候，如果n是完全平方数，那么n模m的剩余只有1，49两个
这是暴力搜索m到65536得到的结果
其中length是满足与number互素的整数的平方剩余的个数，list里是具体的值

(length = 2, number = 120, list = Any[1, 49])
(length = 6, number = 420, list = Any[1, 109, 121, 169, 289, 361])
(length = 6, number = 840, list = Any[1, 121, 169, 289, 361, 529])
(length = 30, number = 4620, list = Any[1, 169, 289, 361, 421, 529, 709, 841, 949, 961, 1369, 1549, 1621, 1681, 1849, 2209, 2269, 2641, 2689, 2809, 2941, 3061, 3301, 3469, 3481, 3529, 3721, 4141, 4321, 4489])
(length = 30, number = 9240, list = Any[1, 169, 289, 361, 529, 841, 961, 1369, 1681, 1849, 2209, 2641, 2689, 2809, 3481, 3529, 3721, 4321, 4489, 5041, 5329, 5569, 6169, 6241, 6889, 7561, 7681, 7921, 8089, 8761])
(length = 180, number = 60060, list = Any[1, 289, 361, 529, 841, 961, 1369, 1621, 1681, 1849, 2209, 2809, 2941, 3301, 3481, 3721, 4489, 5041, 5149, 5329, 5461, 5581, 5989, 6241, 6301, 6829, 6889, 7309, 7921, 8089, 8761, 8941, 9109, 9409, 9949, 10189, 10201, 10609, 10789, 10921, 11449, 11509, 11881, 12301, 12541, 12769, 13381, 13729, 13861, 14221, 14569, 14821, 15409, 16069, 16129, 16501, 16669, 17161, 17341,
18001, 18769, 18841, 18901, 19009, 19189, 19321, 20029, 20101, 20329, 20749, 21121, 21421, 21961, 22201, 22621, 22801, 23269, 23461, 23521, 24049, 24469, 24649, 24781, 25741, 25789, 26569, 26581, 27421, 27589, 27889, 28081, 28141, 28669, 28681, 29929, 30781, 31021, 31201, 31249, 32041, 32341, 32509, 32629, 32761, 33049, 33289, 34189, 34609, 35149, 35281, 36061, 36481, 36661, 37129, 37249, 37489, 37801, 37909, 38509, 38581, 38809, 39601, 39649, 39901, 40429, 40681, 41869, 41941, 42949, 43261, 43429, 44269, 44521, 44641, 45049, 45061, 45109, 45301, 46069, 46201, 46489, 46621, 47161, 48049, 48409, 48889, 49009, 49141, 49261, 49501, 49669, 49729, 49921, 50521, 50821, 50989, 51349, 51529, 51661, 51769, 52441, 53089, 53509, 53629, 53881, 54289, 54349, 54961, 55441, 55969, 56281, 56809, 56989, 57061, 57121, 58081, 58249, 58501, 58969, 59929])

无心人

上面的数据里，有价值的是
(length = 2, number = 120, list = Any[1, 49]) 120=2^3*3*5
(length = 6, number = 840, list = Any[1, 121, 169, 289, 361, 529]) 840=2^3*3*5*7
最后一组，候选太多了，虽然数据比较好，待测试数除非很大，否则不建议选择
(length = 180, number = 60060, list = Any[1, 289, 361, 529, 841, 961, 1369, 1621, 1681, 1849, 2209, 2809, 2941, 3301, 3481, 3721, 4489, 5041, 5149, 5329, 5461, 5581, 5989, 6241, 6301, 6829, 6889, 7309, 7921, 8089, 8761, 8941, 9109, 9409, 9949, 10189, 10201, 10609, 10789, 10921, 11449, 11509, 11881, 12301, 12541, 12769, 13381, 13729, 13861, 14221, 14569, 14821, 15409, 16069, 16129, 16501, 16669, 17161, 17341,
18001, 18769, 18841, 18901, 19009, 19189, 19321, 20029, 20101, 20329, 20749, 21121, 21421, 21961, 22201, 22621, 22801, 23269, 23461, 23521, 24049, 24469, 24649, 24781, 25741, 25789, 26569, 26581, 27421, 27589, 27889, 28081, 28141, 28669, 28681, 29929, 30781, 31021, 31201, 31249, 32041, 32341, 32509, 32629, 32761, 33049, 33289, 34189, 34609, 35149, 35281, 36061, 36481, 36661, 37129, 37249, 37489, 37801, 37909, 38509, 38581, 38809, 39601, 39649, 39901, 40429, 40681, 41869, 41941, 42949, 43261, 43429, 44269, 44521, 44641, 45049, 45061, 45109, 45301, 46069, 46201, 46489, 46621, 47161, 48049, 48409, 48889, 49009, 49141, 49261, 49501, 49669, 49729, 49921, 50521, 50821, 50989, 51349, 51529, 51661, 51769, 52441, 53089, 53509, 53629, 53881, 54289, 54349, 54961, 55441, 55969, 56281, 56809, 56989, 57061, 57121, 58081, 58249, 58501, 58969, 59929])
60060=2^2*3*5*7*11*13
这组可以用
(length = 30, number = 9240, list = Any[1, 169, 289, 361, 529, 841, 961, 1369, 1681, 1849, 2209, 2641, 2689, 2809, 3481, 3529, 3721, 4321, 4489, 5041, 5329, 5569, 6169, 6241, 6889, 7561, 7681, 7921, 8089, 8761])  9240=2^3*3*5*7*11
代替，180/60060=0.002997，30/9240=0.0003247，差距很小

无心人

1. global cnt = (length = 65536, number = 65537, list = [])
   
2. println("Begin")
   
3. for m in range(120, step = 2, stop = 65536)
   
4.     global s = []
   
5.     for i in range(1, stop = m-1)
   
6.         if gcd(i, m) == 1
   
7.             t = i * i % m
   
8.             push!(s, t)
   
9.         end
   
10.     end
    
11.     sort!(s)
    
12.     unique!(s)
    
13.     if cnt.length//cnt.number > length(s)//m
    
14.         global cnt = (length = length(s), number = m, list = s)
    
15.         println(cnt)
    
16.     end
    
17. end
    
18. println("End")
    
19.         
    

复制代码

mathe

求平方根是很快的
https://gmplib.org/manual/Integer-Roots.html
函数mpz_sqrt
得到平方根以后，在自相乘，再判断结果是否等于原整数即可

无心人

1. function isSquareOrDivBy2357(n)
   
2.     if gcd(n, 840) != 1
   
3.         return true
   
4.     end
   
5.     t = n % 840
   
6.     if t in [1, 121, 169, 289, 361, 529]
   
7.         sq = isqrt(n)
   
8.         return sq * sq == n
   
9.     else
   
10.         return false
    
11.     end
    
12. end
    
13. 
    
14. for i in range(1, stop = 65535)
    
15.     t = i * i
    
16.     if !isSquareOrDivBy2357(t)
    
17.         print("error at $i * $i = $t")
    
18.     end
    
19. end
    
20. 
    
21. t = big(10)^67 - 1
    
22. println("$t: ", isSquareOrDivBy2357(t))
    
23. 
    
24. t = t * t
    
25. println("$t: ", isSquareOrDivBy2357(t))
    

复制代码

无心人

mathe 发表于 2019-10-8 10:05
求平方根是很快的
https://gmplib.org/manual/Integer-Roots.html
函数mpz_sqrt

我知道，但是我想先筛选掉绝大部分不是平方的数

无心人

1. function isSquareOrDivBy2357(n)
   
2.     if gcd(n, 840) != 1
   
3.         return true
   
4.     end
   
5.     t = n % 840
   
6.     if t in [1, 121, 169, 289, 361, 529]
   
7.         sq = isqrt(n)
   
8.         return sq * sq == n
   
9.     else
   
10.         return false
    
11.     end
    
12. end
    
13. 
    
14. function isSquare1(n)
    
15.     sq = isqrt(n)
    
16.     return sq * sq == n
    
17. end
    
18. 
    
19. @time for i in range(1, stop = 65535)
    
20.     t = i * i
    
21.     if !isSquareOrDivBy2357(t)
    
22.         print("error at $i * $i = $t")
    
23.     end
    
24. end
    
25. 
    
26. @time for i in range(1, stop = 65535)
    
27.     t = i * i
    
28.     if !isSquare1(t)
    
29.         print("error at $i * $i = $t")
    
30.     end
    
31. end
    
32. 
    
33. t = t * t
    
34. @time println("$t: ", isSquareOrDivBy2357(t))
    
35. @time println("$t: ", isSquare1(t))
    
36. 
    

复制代码

好吧，结果果然不如直接开平方~
0.004386 seconds (14.98 k allocations: 1.828 MiB)
  0.000710 seconds
9999999999999999999999999999999999999999999999999999999999999999996000000000000000000000000000000000000000000000000000000000000000000599999999999999999999999999999999999999999999999999999999999999999960000000000000000000000000000000000000000000000000000000000000000001: true
  0.072349 seconds (197.96 k allocations: 10.581 MiB, 15.72% gc time)
9999999999999999999999999999999999999999999999999999999999999999996000000000000000000000000000000000000000000000000000000000000000000599999999999999999999999999999999999999999999999999999999999999999960000000000000000000000000000000000000000000000000000000000000000001: true
  0.018368 seconds (2.21 k allocations: 106.926 KiB)

gxqcn

借题，问一下，可有此类算法：可快速判定一个大整数，有/无 含非 1 的平方因子？

.·.·.

验证n是完全平方数
sqrt(n)求一下看看是不是整数，不就出结果了吗？

zeroieme

gxqcn 发表于 2019-10-8 17:12
借题，问一下，可有此类算法：可快速判定一个大整数，有/无 含非 1 的平方因子？

3年前我问过，也再问一下。

https://bbs.emath.ac.cn/thread-9242-1-1.html