梅森素数(Mersenne prime)判断, FFT 大数乘法 (非递归), O(n^2 log n), c++

最新推荐文章于 2025-06-01 11:25:02 发布

kevinkitty

最新推荐文章于 2025-06-01 11:25:02 发布

阅读量2.6k

点赞数

CC 4.0 BY-SA版权

分类专栏：编程文章标签： FFT c语言梅森素数数论

本文链接：https://blog.youkuaiyun.com/kevinkitty_love/article/details/51485048

本文介绍了如何运用快速傅里叶变换(FFT)算法实现大数乘法，以此优化卢卡斯-莱默检验法来判断梅森素数。详细阐述了FFT的工作原理，并给出了一个n=8的实例来说明其计算过程。同时讨论了如何选择根w以避免浮点数运算，特别针对2的幂次寻找满足条件的模数MOD。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

原创代码，请勿转载！

梅森素数判定：
卢卡斯-莱默检验法：参考https://zh.wikipedia.org/wiki/%E5%8D%A2%E5%8D%A1%E6%96%AF-%E8%8E%B1%E9%BB%98%E6%A3%80%E9%AA%8C%E6%B3%95

卢卡斯-莱默检验法是迭代算法，需要用到高精度乘法运算。
而现有的乘法运算算法最快的要数快速傅里叶变换算法，
可以将两位n-bit的大整数乘法在O(n lg n)计算出来。

FFT简介：快速傅里叶变换, 参考 http://beige.ucs.indiana.edu/B673/node12.html

n = 8 为例子, 计算两个4位数的乘法。
寻找根w使得 w^n = 1.
A[0..7][0..7]: A[i][j] = w^{i*j}, 其中 w^n = 1.
A^{-1} = w^{-1} ^ {i*j}.

对于大整数X和Y，首先做矩阵乘法运算
X' = A*X
Y' = A*Y
然后对X' Y' 对应bit位对应相乘
Z = X' .* Y' (对应相乘)
最后还原
Z' = A^{-1}*Z

快速傅里叶变换的高效之处在于，A矩阵是具有很强的结构性，因此矩阵乘法运算过程可以优化到O(n lg n).

关于如何选择根w.
采用数论中的 root https://en.wikipedia.org/wiki/Root_of_unity_modulo_n
避免浮点数运算。

对于n=2^k, 寻找大素数MOD，使得存在 w， w^{2^k}=1 (mod MOD).

#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<ctime>
#include<cmath>

#define EXC(a,b) a ^ b ? a ^= b ^= a ^= b : 0
#define EXCB(a,b) a > b ? a ^= b ^= a ^= b : 0
#define FOR(i,a,b) for(int i=(a);i<(b);i++)
#define FORE(i,a,b) for(int i=(a);i<=(b);i++)
#define FORI(TYPE,it,s) for(TYPE::iterator it=s.begin();it!=s.end();it++) 
#define MST(a,b) memset(a,b,sizeof(a))
#define SCF(a) scan