Dixon的因式分解方法及其实现

简介

Dixon的因式分解方法是一种分解大整数的算法。它可以分解任意大的整数，但是需要一些特定条件的限制。

该算法的核心是基于线性同余方程组的求解，其具体步骤如下：

1. 随机选取一组$x_0$和一组$y_0$，其中$x_0$为待分解数模素数的余数，$y_0$为$x_0$的平方模待分解数的余数。
2. 构造一个线性同余方程组，其中每一项为模素数的余数，方程组的解即为$x_0$和$y_0$。
3. 利用高斯消元将线性同余方程组化为上三角矩阵，得到线性无关的同余方程组。
4. 将同余方程组依次代入一个公式，得到$y_0$的的某些非平凡因子的平方。
5. 计算这些因子的最大公因子，得到待分解数的因子。

实现

以下是一个Python实现Dixon算法的代码片段：

import random
from math import gcd


def dixon(n):
    """Dixon's factorization method for integer n"""
    factors = []
    
    # set up variables
    B = 20
    A = []
    x0 = int(n ** 0.5) + 1
    y0 = x0 ** 2 - n
    primes = [2, 3, 5, 7, 11, 13, 17, 19]
    
    # factor x0^2 - n into primes
    for p in primes:
        if y0 % p == 0:
            A.append(p)
            while y0 % p == 0:
                y0 //= p
    
    # factor remaining residues
    while len(A) < B:
        x = random.randint(1, n - 1)
        y = (x ** 2) % n
        z = y - n
        
        # factor z into primes
        for p in primes:
            if z % p == 0:
                A.append(p)
                while z % p == 0:
                    z //= p
        
        # check if enough factors found
        if len(A) == B:
            break
    
    # build matrix
    M = []
    for a in A:
        row = []
        t = y0 % a
        for p in primes:
            count = 0
            while t % p == 0:
                t //= p
                count += 1
            row.append(count % 2)
        M.append(row)
    
    # solve matrix with Gaussian elimination
    r = 0
    c = 0
    while r < B and c < 8:
        # find maximum value in column
        max_row = max(range(r, B), key=lambda i: M[i][c])
        if M[max_row][c] == 0:
            c += 1
        else:
            # swap rows
            M[r], M[max_row] = M[max_row], M[r]
            # eliminate other rows
            for i in range(r + 1, B):
                if M[i][c] == 1:
                    M[i] = [(M[i][j] - M[r][j]) % 2 for j in range(8)]
            r += 1
    
    # check for failure
    if r < B or M[B-1][7] == 0:
        return "Failed to find factors"
    
    # extract factors
    f1 = 1
    f2 = 1
    for i in range(B):
        if M[i][7] == 1:
            t = y0
            for j in range(8):
                if M[i][j] == 1:
                    t *= primes[j]
            f1 *= t
            f2 *= (x0 + t) % n
    
    # compute gcd of factors
    g = gcd(f1 + f2, n)
    if g == n:
        return "Failed to find factors"
    else:
        return g, n // g
    

其中，dixon(n)函数接受一个整数n作为参数，并返回它的因子。

返回的结果是一个元组(a, b)，其中a和b分别是n的两个因子。如果无法找到因子，则返回字符串"Failed to find factors"。