欧拉因子分解法——C语言实现

原创已于 2025-08-17 21:15:30 修改 · 476 阅读

1 ·

CC 4.0 BY-SA版权

文章标签：

#c语言 #算法 #开发语言

于 2024-08-31 15:21:43 首次发布

欧拉因子分解法原理详解

欧拉因子分解法是一种基于数论的整数分解算法，特别适用于可以表示为两个平方和的正整数。该方法由数学家莱昂哈德·欧拉（Leonhard Euler）提出，其核心思想是：如果一个数有多个不同的两平方和表示，则可以利用这些表示之间的差异来找到非平凡因子。以下将逐步详细解释其原理、数学基础和实现步骤，确保内容清晰易懂。

1. 基本概念

两平方和定理：一个正整数 ( n ) 可以表示为两个整数的平方和，即 ( n = a^2 + b^2 )，其中 ( a ) 和 ( b ) 是整数。例如，( 5 = 1^2 + 2^2 )，( 13 = 2^2 + 3^2 )。
欧拉因子分解法的前提：如果 ( n ) 有至少两个不同的两平方和表示（即表示方式不通过符号变换或顺序交换得到），那么 ( n ) 一定是合数（非质数），且可以通过这些表示找到其因子。这里的“不同”指 ( (a, b) \neq (c, d) ) 且 ( (a, b) \neq \pm (d, c) ) 等。

2. 原理与数学基础

假设 ( n ) 有两个不同的两平方和表示： $$ n = a^2 + b^2 $$ $$ n = c^2 + d^2 $$ 其中 ( a, b, c, d ) 是整数，且表示独立（即 ( t = a d - b c \neq 0 )）。

关键恒等式：从上述表示出发，可以推导出以下关系： $$ n^2 = (a^2 + b^2)(c^2 + d^2) $$ 利用代数恒等式，右边展开为： $$ (a^2 + b^2)(c^2 + d^2) = (a c + b d)^2 + (a d - b c)^2 $$ 因此： $$ n^2 = (a c + b d)^2 + (a d - b c)^2 $$ 令 ( t = a d - b c )，则上式简化为： $$ n^2 = (a c + b d)^2 + t^2 $$
- 由于 ( n^2 ) 和 ( (a c + b d)^2 ) 都是整数，且 ( n^2 - (a c + b d)^2 = t^2 )，这意味着： $$ n^2 - (a c + b d)^2 = t^2 $$ 即： $$ (n - (a c + b d))(n + (a c + b d)) = t^2 $$
- 由此可得 ( n ) 整除 ( t^2 )（因为左边是 ( t^2 )，而 ( n ) 是因子），即： $$ n \mid t^2 $$
因子存在性：既然 ( n ) 整除 ( t^2 )，且 ( t \neq 0 )，那么 ( n ) 和 ( t ) 必有公因子（除非 ( n ) 是质数，但前提已保证 ( n ) 有多个表示，故 ( n ) 为合数）。具体地，计算最大公约数 ( g = \gcd(t, n) )。如果 ( g > 1 ) 且 ( g < n )，则 ( g ) 是 ( n ) 的非平凡因子，另一个因子为 ( n / g )。

这一原理的核心在于：多个两平方和表示诱导出差值 ( t )，而 ( n ) 与 ( t ) 的公因子揭示了 ( n ) 的分解结构。

3. 算法步骤

欧拉因子分解法的实现步骤如下（假设已找到两个不同的表示）：

输入：一个合数 ( n )，且 ( n ) 至少有两个不同的两平方和表示 ( n = a^2 + b^2 ) 和 ( n = c^2 + d^2 )。
计算差值：计算 ( t = a d - b c )。如果 ( t = 0 )，则表示相同，需重新寻找表示。
计算最大公约数：计算 ( g = \gcd(t, n) )。这里 ( \gcd ) 是最大公约数函数。
检查因子：
- 如果 ( g > 1 ) 且 ( g < n )，则 ( g ) 是 ( n ) 的因子，另一个因子为 ( n / g )。
- 如果 ( g = 1 ) 或 ( g = n )，则需尝试其他表示或组合（如交换 ( a,b ) 或 ( c,d ) 的顺序）。
输出：( n ) 的因子对 ( (g, n/g) ).

注意：该方法依赖于高效找到多个两平方和表示。实践中，可通过遍历可能整数对或利用数论优化来实现。

4. 实例演示

以 ( n = 65 ) 为例，演示欧拉因子分解法的应用。

步骤1：找到两个不同的两平方和表示。
- ( 65 = 1^2 + 8^2 ) （因为 ( 1 + 64 = 65 ))，所以 ( a = 1 ), ( b = 8 )。
- ( 65 = 4^2 + 7^2 ) （因为 ( 16 + 49 = 65 ))，所以 ( c = 4 ), ( d = 7 )。
- 表示不同（因为 ( (1,8) \neq (4,7) ))。
步骤2：计算 ( t )： $$ t = a d - b c = 1 \cdot 7 - 8 \cdot 4 = 7 - 32 = -25 $$
- 取绝对值或直接使用 ( t = -25 )（因为 ( \gcd ) 对符号不敏感）。
步骤3：计算 ( g )： $$ g = \gcd(|t|, n) = \gcd(25, 65) $$
- 最大公约数为 5（因为 ( 25 = 5^2 ), ( 65 = 5 \times 13 ))。
步骤4：输出因子：( g = 5 > 1 ) 且 ( g < 65 )，所以因子为 5 和 ( 65 / 5 = 13 )。
- 验证：( 5 \times 13 = 65 ).

另一个例子：( n = 85 )。

表示：( 85 = 2^2 + 9^2 )（( a=2, b=9 ))，( 85 = 6^2 + 7^2 )（( c=6, d=7 ))。
计算 ( t = 2 \cdot 7 - 9 \cdot 6 = 14 - 54 = -40 )。
( g = \gcd(40, 85) = 5 )。
因子：5 和 17（因为 ( 85 / 5 = 17 )）。

5. 优缺点与适用性

优点：方法简洁，基于优雅的数论原理，特别适合能快速找到多个两平方和表示的数（如形如 ( 4k+1 ) 的质数的乘积）。
缺点：
- 并非所有数都有多个两平方和表示（例如，形如 ( 4k+3 ) 的质数不能表示为两平方和）。
- 寻找多个表示可能计算量大，尤其对大 ( n )。
应用：主要用于教学和理论分析，现代因子分解算法（如Pollard rho）更高效，但欧拉方法在数论中仍有重要地位。

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

// 计算最大公约数（GCD）
long gcd(long a, long b) {
    while (b != 0) {
        long temp = b;
        b = a % b;
        a = temp;
    }
    return a;
}

// 查找n的所有平方和表示（n = a^2 + b^2 且 a >= b >= 0）
int find_square_sums(long n, long reps[][2], int max_reps) {
    int count = 0;
    long sqrt_n = (long)sqrt(n);
    
    // 从sqrt(n)向下遍历到sqrt(n/2)
    for (long a = sqrt_n; a * a >= n / 2; a--) {
        long b2 = n - a * a;  // 计算b^2
        if (b2 < 0) continue;
        
        long b = (long)sqrt(b2);  // 取b的整数部分
        // 检查b^2是否等于b2
        if (b * b == b2) {
            // 确保a >= b（避免重复）
            long a_val = (a >= b) ? a : b;
            long b_val = (a >= b) ? b : a;
            
            // 避免重复记录相同的表示
            int duplicate = 0;
            for (int i = 0; i < count; i++) {
                if (reps[i][0] == a_val && reps[i][1] == b_val) {
                    duplicate = 1;
                    break;
                }
            }
            if (!duplicate) {
                reps[count][0] = a_val;
                reps[count][1] = b_val;
                count++;
                if (count >= max_reps) break;
            }
        }
    }
    return count;
}

// 欧拉因子分解法主函数
void euler_factorization(long n) {
    const int max_reps = 10;
    long reps[max_reps][2];  // 存储平方和表示
    int count = find_square_sums(n, reps, max_reps);
    
    // 检查是否找到至少两组表示
    if (count < 2) {
        printf("Euler method failed: Only %d representations found for n = %ld\n", count, n);
        return;
    }
    
    printf("Found %d square sum representations:\n", count);
    for (int i = 0; i < count; i++) {
        printf("Representation %d: %ld^2 + %ld^2\n", i+1, reps[i][0], reps[i][1]);
    }
    
    // 尝试所有表示对组合
    for (int i = 0; i < count; i++) {
        long a1 = reps[i][0], b1 = reps[i][1];
        for (int j = i + 1; j < count; j++) {
            long a2 = reps[j][0], b2 = reps[j][1];
            printf("\nTrying representations: (%ld,%ld) and (%ld,%ld)\n", a1, b1, a2, b2);
            
            // 计算关键量
            long x = labs(a1 * b2 - b1 * a2);  // |a1*b2 - b1*a2|
            long y = labs(a1 * b2 + b1 * a2);  // |a1*b2 + b1*a2|
            
            // 计算可能的因子
            long p = gcd(x, n);
            long q = gcd(y, n);
            
            printf("x = |%ld*%ld - %ld*%ld| = %ld\n", a1, b2, b1, a2, x);
            printf("y = |%ld*%ld + %ld*%ld| = %ld\n", a1, b2, b1, a2, y);
            printf("gcd(x, n) = gcd(%ld, %ld) = %ld\n", x, n, p);
            printf("gcd(y, n) = gcd(%ld, %ld) = %ld\n", y, n, q);
            
            // 检查是否找到非平凡因子
            if (p > 1 && p < n) {
                printf("\nFactorization SUCCESS: %ld = %ld * %ld\n", n, p, n / p);
                return;
            }
            if (q > 1 && q < n) {
                printf("\nFactorization SUCCESS: %ld = %ld * %ld\n", n, q, n / q);
                return;
            }
        }
    }
    printf("\nEuler method failed: No factors found for n = %ld\n", n);
}

int main() {
    // 测试数据：1000009（70191551 无法表示为两个平方和）
    long n = 1000009;
    printf("Attempting to factorize n = %ld using Euler's method...\n", n);
    euler_factorization(n);
    return 0;
}