序
对于椭圆曲线方面的算法中,会常常遇到数据量很大的多项式乘除法,譬如:给定某个整数N,如何获得一个椭圆曲线方程使得其有理点群阶为N(这里,就涉及了Hilbert modular polynomials);抑或是给定一个椭圆曲线的方程,去计算该曲线在某个有限域的有点个数(这里就涉及了Division Polynomial)。
多项式除法可以在一定情形下转变为多项式乘法。从而如何快速计算多项式乘法是个很有意义的问题。实际上,处理那两个多项式更好的算法是用NTT快速数论变换。不过现在有的地方还没有搞清楚怎么弄,所以先整理整理FFT算法的实现。
这篇文章偏向于python的实现。。。。(毕竟具体教学的文章有好多。。。。。
我对傅里叶变换写了一些不是那么“直接”的文章。。。可供参考Alepha E:复变随记(五)从傅里叶变换到四平方和,以及一些几何推广zhuanlan.zhihu.comAlepha E:聊一聊【向量】视角下的函数(例子:傅里叶变换)zhuanlan.zhihu.com
(这两篇是我比较早写的,有些理解略显幼稚。。。。)
数据结构
我这里打算是用Python来做的,Python中自带一个dictionary数据类型,基于该数据类型来对多项式类进行定义。
我们首先要搞清楚,我们希望多项式要有怎样的代数结构,即:多项式指数到多项式系数的单射。从而,对于多项式
存储形式:dictionary的key值为多项式的指数,values为多项式相对应指数的系数。
在构建时,可以直接输入符合这样多项式的dictionary。
在构建多项式之前,至少要有加法。做加法时,如
,以
的dictionary为准,遍历
的指数并判断是否在
的keys中:若在,则对应key值的value做加法,并且将和为令的对应指数消去;
若不在,则直接添加到返回结果中。
傅里叶变换
对于普通的乘法,譬如
,需要做大约
次。
进行傅里叶变换时,是为了将原先:
转化为
当然,上面这两种表示是可以互逆的。
注意,单位根有这么样的性质:
,以及
从而可以重复的处理计算
上(这里
是作为区分的两个值,与上面多项式系数无关)
具体举个例子,比方说,对这个多项式做傅里叶变换
先人工整理,可有
进一步
这样,能发现,我们再计算时:对于形式为
的运算,执行了4次;
对于形式为
的运算,执行了2次;
对于形式为
的运算,执行了1次。
用单位根来代换,即可用单位根的“折半”性质去做。
这才是傅里叶变换做乘法的精髓所在。
通过卷积定理,可以满足
这就保证了变换前后满足乘积的性质了。
附上代码
def FFT4polynomial(poly,max_deg):
ind_coe = poly.ind_coe #这里ind_coe是多项式的指数系数dictionary
if '1' in bin(max_deg)[3:]: #这里的两个判断是做初始化用的
max_deg = 2 ** (len(bin(max_deg)[2:]))
M = max_deg
ind_coe = {i: ind_coe[i] if i in ind_coe else 0 for i in range(max_deg)}
else:
M = max_deg
ind_coe = {i: ind_coe[i] if i in ind_coe else 0 for i in range(max_deg)}
ans = {}
m = len(bin(max_deg)[3:])
for item in range(max_deg):
t = 1
tmp_ind = ind_coe.copy()
max_deg = len(tmp_ind)
while len(tmp_ind)!=1:
unit_root = np.exp(2 * math.pi * Complex(0, 1)*item*2**(m-t)/(M))
tmp_ind = {i: tmp_ind[i] + tmp_ind[i + max_deg // 2] * unit_root for i in tmp_ind if i < max_deg // 2}
t += 1
max_deg = max_deg//2
ans[item] = tmp_ind[0]
return Polynomials(ans)
对于逆变换,只需要对单位根部分的指数取个复数即可。
对于多项式的数据结构,,各位就自己实现吧。。。