【问题标题】:Fermat Factorisation with Python使用 Python 进行费马因式分解
【发布时间】:2013-12-09 06:34:24
【问题描述】:

Python 新手,不知道为什么我的费马分解方法失败了?我认为这可能与实施大量数字的方式有关,但我对语言的了解不足,无法确定我哪里出错了。

当 n=p*q 与 p 和 q 非常接近(如在彼此之间大约 20 以内)时,下面的代码可以工作,但如果它们相距更远,似乎会永远运行。例如,n=991*997 代码可以正常工作并在 n=104729*104659 也是如此。但是,如果我将其更改为n=103591*104659,它将永远运行(好吧,我让它运行 2 小时然后停止它)。

任何正确方向的观点将不胜感激!

代码:

import math

def isqrt(n):
  x = n
  y = (x + n // x) // 2
  while y < x:
    x = y
    y = (x + n // x) // 2
  return x

n=103591*104729

a=isqrt(n) + 1
b2=a*a - n
b=isqrt(b2)

while b*b!=b2:
  a=a+1
  b2=b2+2*a+1
  b=isqrt(b2)

p=a+b
q=a-b

print('a=',a,'\n')
print('b=',b,'\n')
print('p=',p,'\n')
print('q=',q,'\n')
print('pq=',p*q,'\n')
print('n=',n,'\n')
print('diff=',n-p*q,'\n')

【问题讨论】:

    标签: python factorization


    【解决方案1】:

    我在Wikipedia 上查找了算法,这对我有用:

    #from math import ceil
    
    def isqrt(n):
      x = n
      y = (x + n // x) // 2
      while y < x:
        x = y
        y = (x + n // x) // 2
      return x
    
    def fermat(n, verbose=True):
        a = isqrt(n) # int(ceil(n**0.5))
        b2 = a*a - n
        b = isqrt(n) # int(b2**0.5)
        count = 0
        while b*b != b2:
            if verbose:
                print('Trying: a=%s b2=%s b=%s' % (a, b2, b))
            a = a + 1
            b2 = a*a - n
            b = isqrt(b2) # int(b2**0.5)
            count += 1
        p=a+b
        q=a-b
        assert n == p * q
        print('a=',a)
        print('b=',b)
        print('p=',p)
        print('q=',q)
        print('pq=',p*q)
        return p, q
    
    n=103591*104729
    fermat(n)
    

    我尝试了几个测试用例。这个来自维基百科页面:

    >>> fermat(5959)
    Trying: a=78 b2=125 b=11
    Trying: a=79 b2=282 b=16
    a= 80
    b= 21
    p= 101
    q= 59
    pq= 5959
    (101, 59)
    

    这是您的示例案例:

    >>> fermat(103591*104729)
    Trying: a=104159 b2=115442 b=339
    a= 104160
    b= 569
    p= 104729
    q= 103591
    pq= 10848981839
    (104729, 103591)
    

    查看标有“尝试”的行表明,在这两种情况下,它都很快收敛。

    更新:来自 cmets 的非常长的整数因素如下:

    n_long=316033277426326097045474758505704980910037958719395560565571239100878192955228495343184968305477308460190076404967552110644822298179716669689426595435572597197633507818204621591917460417859294285475630901332588545477552125047019022149746524843545923758425353103063134585375275638257720039414711534847429265419
    
    fermat(n_long, verbose=False)
    a= 17777324810733646969488445787976391269105128850805128551409042425916175469326288448917184096591563031034494377135896478412527365012246902424894591094668262
    b= 157517855001095328119226302991766503492827415095855495279739107269808590287074235
    p= 17777324810733646969488445787976391269105128850805128551409042425916175469483806303918279424710789334026260880628723893508382860291986009694703181381742497
    q= 17777324810733646969488445787976391269105128850805128551409042425916175469168770593916088768472336728042727873643069063316671869732507795155086000807594027
    pq= 316033277426326097045474758505704980910037958719395560565571239100878192955228495343184968305477308460190076404967552110644822298179716669689426595435572597197633507818204621591917460417859294285475630901332588545477552125047019022149746524843545923758425353103063134585375275638257720039414711534847429265419
    

    【讨论】:

    • 非常感谢!这似乎适用于小数字,但是当我使用大数字时,我得到“long int too large to convert to float”错误。我认为这是因为 n**0.5 声明?你有什么办法让这个工作在更大的数字上?例如n = 316033277426326097045474758505704980910037958719395560565571239100878192955228495343184968305477308460190076404967552110644822298179716669689426595435572597197633507818204621591917460417859294285475630901332588545477552125047019022149746524843545923758425353103063134585375275638257720039414711534847429265419 跨度>
    • @user3081739 这并不难:我使用了你的 isqrt 算法,它很容易处理你的长号码,尽管我在计算时确实出去喝咖啡休息了。
    • 太完美了,谢谢!你知道为什么我的代码失败了吗?我可以看到两者之间的唯一区别是 b2 的迭代计算。你使用 b2 = a*a - n 而我使用 b2 = b2 + 2*a + 1?
    • @user3081739 我相信你是对的,这是唯一的重大变化。
    【解决方案2】:

    错误是在增加 a 之后进行加法,因此新值不是 a 的平方。 这按预期工作:

    while b*b!=b2:
      b2+=2*a+1  
      a=a+1  
      b=isqrt(b2)
    

    对于大数字,它应该比计算具有更多位数的正方形更快。

    【讨论】:

      猜你喜欢
      • 2012-05-15
      • 1970-01-01
      • 2013-04-07
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 2018-12-12
      相关资源
      最近更新 更多