數論概論中的費馬降階法求一個質數分解成2個數的平方和與窮舉法的效率比較

原始碼改寫自《HDU 3542 （費馬降階法）》原文連結：https://blog.csdn.net/Feynman1999/article/details/81903899

import random as rd
import time

def ran(n):
    if (n%4!=1):return 0
    while True:
        a=rd.randint(1,10**16)%(n-1)+1 #1~ n-1
        b=pow(a,(n-1)//4,n);
        if(b*b%n==n-1):
          if b<=(n-1)/2 :return b
          else: return n-b

def solve(p):
    if(p==2): return [1,1];
    if(p%4!=1): return mp[-1,-1];
    A=ran(p); B=1
    tmp=1;
    M=(tmp*A*A+tmp*B*B)//p
    while M>1:
        u=(A%M+M)%M;v=(B%M+M)%M
        if(2*u>M): u-=M
        if(2*v>M): v-=M
        ta=A
        A=(tmp*u*A+tmp*v*B)//M
        B=(tmp*v*ta-tmp*u*B)//M
        M=(tmp*u*u+tmp*v*v)//M
    return [min(abs(A),abs(B)),max(abs(A),abs(B))]

def prime(n):
 sieve = [True] * (n//2+1)
 for number in range(3, int(n ** 0.5) + 1,2):
  if sieve[(number-1)//2]:
   for index in range(number * number, n, number*2):
    sieve[(index-1)//2] = False
 return [x for x in range(1,n,2) if sieve[(x-1)//2]]

def test(n):
 s=0
 p=[i for i in prime(n) if i%4==1 and i>1]
 for i in p:
  pr=solve(i)
  s+=pr[0]+pr[1]
 return s

def test2(n):
 s=0
 p=[i for i in prime(n) if i%4==1 and i>1]
 for i in p:
  ii=i*i
  for x in range(1,int(((i/2)**0.5))+1):
   y0=(i -x*x)**0.5
   if y0%1==0:s+=(x+int(y0));break
 return s

t=time.time();print(test(10**6));print(time.time()-t)

t=time.time();print(test2(10**6));print(time.time()-t)

>>> 
>>> t=time.time();print(test(10**6));print(time.time()-t)
32266295
2.1802237033843994
>>> t=time.time();print(test2(10**6));print(time.time()-t)
32266295
9.938273668289185
>>> t=time.time();print(test(10**7));print(time.time()-t)
869849116
21.672968864440918
>>> t=time.time();print(test2(10**7));print(time.time()-t)
869849116
263.66866731643677