NTT数论变换(三)

最新推荐文章于 2024-04-02 09:49:36 发布
最新推荐文章于 2024-04-02 09:49:36 发布
阅读量1k 收藏 1
点赞数 1
分类专栏： 密码
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/wxkhturfun/article/details/111937428
版权

上一篇：https://blog.csdn.net/wxkhturfun/article/details/111186205
下一篇：https://blog.csdn.net/wxkhturfun/article/details/114596494

前言

前面两篇分析的都是“库里－图基”类型的蝶形运算，蝶形运算操作是先乘后加：

u = a[k] % P
t = w * (a[k + int(h / 2)] % P) % P
a[k] = (u + t) % P
a[k + int(h / 2)] = ((u - t) % P + P) % P

还存在另一种蝶形运算：“桑德－图基”，是先加后乘：

u = (a[k] + a[k + int(h / 2)])%P
t = w * ((a[k]-a[k + int(h / 2)])%P+P)%P
a[k] = u
a[k + int(h / 2)] = t

两种操作复杂度一致，但是考虑硬件开销，桑德－图基只用了一次乘法，而库里－图基用了两次，所以决定还将其实现一遍。

１.变化

桑德－图基相当于库里－图基的“逆”，也就是说整个系统的蝶形运算的输入改为输出，各参数的顺序不变，常量的值不变，箭头方向相反。
所以我们只需要改变两个地方：一个是wn的计算，还有就是Fntt()
首先给出原来的库里－图基：

#!/usr/bin/env python
# -*- coding: utf-8 -*-
# author:ys
import random
def quick_mod(a: int, b: int,m: int):
     ans = 1
     a %= m
     while(b):
          if(b & 1):
               ans = ans * a % m
               b-=1
          b >>= 1
          a = a * a % m
     return ans
def GetWn(NUM:int ,wn:list):
    for i in range(NUM):
        t = 1 << i
        wn[i] = quick_mod(G, int((P - 1) / t), P)
def Prepare(A:list, B:list, a:list ,b:list, length=1 ):
     length_A = len(A)
     length_B = len(B)
     while((length <= 2 * length_A) or (length <= 2 * length_B)):
          length <<= 1
     A.extend( [0]*(length - length_A) ) 
     B.extend( [0]*(length - length_B) ) 
     a.extend([0]*length)
     b.extend([0]*length)
     for i in range(length_A):
          A[length - 1 - i] = A[length_A - 1 - i]
     for i in range(length - length_A):
        A[i] = 0
     for i in range(length_B):
          B[length - 1 - i] = B[length_B - 1 - i]
     for i in range(length-length_B):
        B[i] = 0
     for i in range(length):
          a[length - 1 - i] = A[i]
     for i in range(length):
          b[length - 1 - i] = B[i]
     return length
def  Rader(a:list, length:int):
     j = length >> 1
     for i in range(1,length-1):
          if(i < j): 
               a[i], a[j] = a[j], a[i] 
          k = length >> 1
          while(j >= k):
               j -= k
               k >>= 1
          if(j < k):
               j += k
     #return a
def Fntt(a:list ,length:int ,on:int ):

     #a=Rader(a, length)
     Rader(a, length)
     id = 0
     h = 2
     while (h <= length):
          print(h) 
          id += 1
          for j in range(0,length,h):
               w = 1
               for k in range(j,j+int(h/2)):
                    u = a[k] % P
                    t = w * (a[k + int(h / 2)] % P) % P
                    a[k] = (u + t) % P
                    a[k + int(h / 2)] = ((u - t) % P + P) % P
                    w = w * wn[id] % P 
          h <<= 1             
     if(on == -1):
          for i in range(1,int(length/2)):
               a[i], a[length - i] = a[length - i],a[i]
          Inv = quick_mod(length, P - 2, P)
          for i in range(length):
               a[i] = a[i] % P * Inv % P
 
def Conv(a:list ,b:list ,n:int ):
     Fntt(a, n, 1)
     Fntt(b, n, 1)
     for i in range(n):
          a[i] = a[i] * b[i] % P
     Fntt(a, n, -1)
     
def Transfer(a:list ,n:list ):
     t = 0
     for i in range(n):
          a[i] += t
          if(a[i] > 9):
               t = int(a[i] / 10)
               a[i] %= 10
          else:
               t = 0
def Print (a:list, n:int):
     flag = True
     for i in range(n-1,-1,-1):
          if(a[i] != 0  and  flag):
               print( a[i],end="")
               flag = False
          elif( not flag):
               print(a[i],end="")
     print("")#这里的print是用于换行的需要，因为前面都是end=""


N = 256
P = 7681
G = 17
NUM = 13
wn=[0]*NUM
GetWn(NUM,wn)
if __name__ == "__main__":
     a=[]
     b=[]
     for i in range (256):
          #a.append(i)
          #b.append(i)
          a.append(random.randint(0,256))
          b.append(random.randint(0,256))
     length=256
     with open("a.txt",'w')as f:
          for i in a:
               '''i=str(hex(i))
               i=i[2:]
               if len(i)<2:
                    i='0'*(2-len(i))+i
                    '''
               f.write(str(i))
               f.write("\n")
          f.close()
     with open("b.txt",'w')as f:
          for i in b:
               '''i=str(hex(i))
               i=i[2:]
               if len(i)<2:
                    i='0'*(2-len(i))+i'''
               f.write(str(i))
               f.write("\n")
          f.close()
     Conv(a,b,256)
     with open("mul.txt",'w')as f:
          for i in a:
               i=str(hex(i))
               i=i[2:]
               if len(i)<4:
                    i='0'*(4-len(i))+i
               f.write(i)
               f.write("\n")

下面给出桑德－图基：

#!/usr/bin/env python
# -*- coding: utf-8 -*-
# author:ys
import random
def quick_mod(a: int, b: int,m: int):
    ans = 1
    a %= m
    while(b):
        if(b & 1):
            ans = ans * a % m
            b-=1
        b >>= 1
        a = a * a % m
    return ans
def GetWn(NUM:int ,wn:list):
    for i in range(NUM):
        t = 1 << (NUM-1-i)
        wn[i] = quick_mod(G, int((P - 1) / t), P)
def Prepare(A:list, B:list, a:list ,b:list, length=1 ):
    length_A = len(A)
    length_B = len(B)
    while((length <= 2 * length_A) or (length <= 2 * length_B)):
        length <<= 1
    A.extend( [0]*(length - length_A) ) 
    B.extend( [0]*(length - length_B) ) 
    a.extend([0]*length)
    b.extend([0]*length)
    for i in range(length_A):
        A[length - 1 - i] = A[length_A - 1 - i]
    for i in range(length - length_A):
        A[i] = 0
    for i in range(length_B):
        B[length - 1 - i] = B[length_B - 1 - i]
    for i in range(length-length_B):
        B[i] = 0
    for i in range(length):
        a[length - 1 - i] = A[i]
    for i in range(length):
        b[length - 1 - i] = B[i]
    return length
def  Rader(a:list, length:int):
    j = length >> 1
    for i in range(1,length-1):
        if(i < j): 
            a[i], a[j] = a[j], a[i] 
        k = length >> 1
        while(j >= k):
            j -= k
            k >>= 1
        if(j < k):
            j += k
    #return a
def Fntt(a:list ,length:int ,on:int ):
    id = 0
    h = length
    while(h>=2):
        '''for j in range(length-h,-h,-h):
            w = 1
            for k in range(j,j+int(h/2)):
                '''
        for j in range(0,length,h):
            w = 1
            for k in range(j,j+int(h/2)):
                u = (a[k] + a[k + int(h / 2)])%P
                t = w * ((a[k]-a[k + int(h / 2)])%P+P)%P
                a[k] = u
                a[k + int(h / 2)] = t
                w = w * wn[id] % P 
        id+=1
        h>>=1  
    Rader(a,length)          
    if(on == -1):
        for i in range(1,int(length/2)):
            a[i], a[length - i] = a[length - i],a[i]
        Inv = quick_mod(length, P - 2, P)
        for i in range(length):
            a[i] = a[i] % P * Inv % P

def Conv(a:list ,b:list ,n:int ):
    Fntt(a, n, 1)
    Fntt(b, n, 1)
    for i in range(n):
        a[i] = a[i] * b[i] % P
    Fntt(a, n, -1)
    
def Transfer(a:list ,n:list ):
    t = 0
    for i in range(n):
        a[i] += t
        if(a[i] > 9):
            t = int(a[i] / 10)
            a[i] %= 10
        else:
            t = 0
def Print (a:list, n:int):
    flag = True
    for i in range(n-1,-1,-1):
        if(a[i] != 0  and  flag):
            print( a[i],end="")
            flag = False
        elif( not flag):
            print(a[i],end="")
    print("")#这里的print是用于换行的需要，因为前面都是end=""


N = 256
P = 7681
G = 17
NUM = 9
wn=[0]*NUM
GetWn(NUM,wn)
if __name__ == "__main__":
    length=256
    a=[]
    b=[]
    with open("a.txt","r") as f:
        for line in f.readlines():
            line = line.strip('\n')  #去掉列表中每一个元素的换行符
            nums = line.split()
            for i in nums:
                a.append(int(i))
    with open("b.txt","r") as f:
        for line in f.readlines():
            line = line.strip('\n')  #去掉列表中每一个元素的换行符
            nums = line.split()
            for i in nums:
                b.append(int(i))
    Conv(a,b,256)
    with open("GS.txt",'w')as f:
        for i in a:
            i=str(hex(i))
            i=i[2:]
            if len(i)<4:
                    i='0'*(4-len(i))+i
            f.write(i)
            f.write("\n")

上面两个代码是相辅相成的，先执行第一个，再执行第二个，第二个输入与第一个保持一致，然后对比两者的输出：diff -u GS.txt mul.txt,如果终端没有任何打印，说明两者完全相同。
另外，有一点需要注意，代码给出的是256序列长度的数相乘，也就是８级蝶形运算，当不是256长时，库里－图基给出的代码不变，但是桑德－图基给出的代码中的wn的顺序要变!!! 两者的wn要在X级下保证完全颠倒（Ｘ是２的次幂，例子中Ｘ=8）
改变：wn的顺序颠倒、id从０开始计数。Fntt内部由先码位倒置(Rader)再计算，变为先计算，再进行码位倒置。其他还有些微不足道的改变，不予赘述。

关于两种蝶形运算的具体操作，可参见：《数字信号处理》清华大学出版社　程佩青　第四版

二.硬计算（不用蝶形运算）

def quick_mod(a: int, b: int,m: int):
     ans = 1
     a %= m
     while(b):
          if(b & 1):
               ans = ans * a % m
               b-=1
          b >>= 1
          a = a * a % m
     return ans
def originNTT(a:list):
     g=quick_mod(G, int((P - 1)/N), P)
     re=[0]*N
     for i in range(N):
          for j in range(N):
               fuck=quick_mod(g, int(i*j), P)
               re[i]+= a[j]*fuck
               re[i] = re[i]%7681
     return re
G=17
P=7681
N=256
if __name__=="__main__":
     a=[]
     with open("a.txt","r") as f:
          for line in f.readlines():
               line = line.strip('\n')  #去掉列表中每一个元素的换行符
               nums = line.split()
               for i in nums:
                    a.append(int(i))
     result=originNTT(a)
     with open("originNTT.txt",'w')as f:
        for i in result:
            i=str(hex(i))
            i=i[2:]
            if len(i)<4:
                    i='0'*(4-len(i))+i
            f.write(str(i))
            f.write("\n")

diff originNTT.txt GS.txt

结果相同

Greate AUK
关注
  • 1
    点赞
  • 踩
  • 1
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 0
    评论
专栏目录
NTT数论变换
03-20
这是我学习数论变换时做的一个程序<br>把数论变换从理论到实践有一段很长的路，在这里我发布出我的一些体会，希望对大家有帮助。
快速数论变换
05-25
该代码是，快速数论变换NNT的C语言源代码，里面有5个函数，采用递归的思想。
EhViewer: 高级漫画阅读器的开源实现
最新发布
gitblog_00009的博客
04-02 5192
EhViewer: 高级漫画阅读器的开源实现 项目地址:https://gitcode.com/exzhawk/EhViewer EhViewer 是一个开源的Android应用，专为漫画爱好者设计，提供了便捷、个性化的漫画阅读体验。项目基于先进的网络爬虫技术和本地缓存策略，使用户能够流畅地浏览和下载来自多个来源的海量漫画资源。 技术分析 网络爬虫 EhViewer 使用智能的网络爬虫技术，能够自...
【XSY2535】整数 NTT
ez_yww的博客
10-06 321
题目描述　　问有多少个满足以下要求的kk进制数：　　　1.每个数字出现的次数不超过nn　　　2.00没有出现过　　　3.若gi,j=0g_{i,j}=0，则ii不能出现恰好jj次。　　两次询问之间会修改gg中一个位置的值（00变11或11变00）。　　输出所有询问的答案的和。　　3≤k≤10,n≤14000,m≤203\leq k\leq 10,n\leq 14000,m\leq 20　　模数p=7
大整数类代码
hehepig的博客
09-27 160
#include <cstdio> #include <vector> #include <iostream> using namespace std; namespace NTT { const long long P=998244353; // =119*2^23+1 const long long G=3; vector<int> r; long long qpow(long long x, int k) { int ret
8点流水线型FFT的Verilog实现
06-28
8点流水线型FFT的Verilog实现，本来是从网上下载的文档，我自己整理成了QuartusII工程文件，并且自己加了很多注释。可以编译通过，但是我没有验证算法是不是完全正确。供大家参考学习。
NTT数论变换（四）
wxkhturfun的博客
03-09 1023
上一篇：https://blog.csdn.net/wxkhturfun/article/details/111937428 Stockham算法 关于该算法，可以详细见链接：http://wwwa.pikara.ne.jp/okojisan/otfft-en/stockham1.html Stockham和FFT的区别在于是否需要额外进行码位变换，在基于库里-图基的FFT(或NTT)中，需要首先进行Rader变换，而基于桑德-图基的FFT（或NTT）中，最后需要进行Rader变换，而存储数组的下标是始终在
快速数论变换NTT
hefenghhhh的专栏
12-17 2861
快速数论变换(NTT) 一.相关概念: 1.剩余系: 一个整数去模n得到的结果的集合称为n的简化剩余系,即小于n的非负整数,记为ZnZ_nZn​ 2.简化剩余系: 在n的剩余系中与n互质的元素的集合,称为n的简化剩余系,记为Zn∗Z_n^*Zn∗​ 3.欧拉函数: n的简化剩余系中元素的个数,称为欧拉函数,记为ϕ(n)\phi(n)ϕ(n) 3.原根: 对于正整数x和n,如果x模n的阶为ϕ(n)\...
ntt:数论变换（NTT）
03-31
恩特数论变换（NTT） 我在NTT进行试验的游乐场。
NTT.rar_NTT_NTT c++_NTT%20c%2B%2B_nttc++_快速数论变换
09-20
用c++实现的快速数论变换(NTT)，利用原根性质，实现和FFT的功能，但是效率比FFT高。
数论变换：此代码用于求一个序列的NTT-matlab开发
05-30
在这里，我们不是将数字提高到很大的幂，而是进行递归乘法以减少计算和错误。 将结果与 DFT 进行比较。
矩形变换和数论变换
11-13
卷积的快速实现和离散傅立叶变换（discrete Fourier transform，DFT）的计算都是信号和图像处理中经常遇到的问题。在实践中，这些操作通常都是用快速傅立叶变换（fast Fourier transform，FFT）算法实现的。NTT在...
NTT快速数论变换
ruizhangj的博客
08-23 315
NTT NTTNTTNTT和FTTFTTFTT十分相似，在FTTFTTFTT中使用了单位复根ω\omegaω来取函数值，在NTTNTTNTT中则使用原根ggg来取函数值。在NTTNTTNTT中，所有的数都是模意义下的。 ggg​的一些性质 g1,g2,⋯ ,gp−1g^1,g^2,\cdots,g^{p-1}g1,g2,⋯,gp−1在模ppp意义下互不相同 根据欧拉定理可知gp−1≡1g^{p-1} \equiv 1gp−1≡1 用ggg来替换ω\omegaω实现FFTFFTFFT到NTTNTTNTT的
快速傅里叶变换
linjiayang2016的博客
05-27 1248
分类 缩写 全称 作用 时间复杂度 DFT 离散傅立叶变换 时频域转换 O(n2)O(n2)O(n^2) FFT 快速傅立叶变换 时频域转换(((带精度误差))) O(大常数+nlog2n)O(大常数+nlog2n)O(大常数+nlog_2n) NTT/FNTT 快速数论变换 模意义下的时频域转换 O(小常数+n......
快速数论变换（NTT）板子
დ♂Hany01's Blog Space♂ღ
02-23 450
原理讲解？反正我这么菜当然只能背板子了。。原理讲解？反正我这么菜当然只能背板子了。。\color{white}{原理讲解？反正我这么菜当然只能背板子了。。} /************************************************ * Au: Hany01 * Date: Feb 23th, 2018 * Prob: [UOJ34] NTT Template *...
三对角、五对角追赶法求解线性方程组
wxkhturfun的博客
05-28 5282
前言 谨以此文纪念2022-5-27 日19：00-21：30学业生涯最后一次考试《数值分析》 参考《五对角线性方程组追赶法》——王礼广等；以及《求解五对角方程组追赶法》——夏爱生等；以及《数值分析》课本第二版 一.三对角 1.1问题引入 问题引入： A={613861...........................86186}，A= \begin{Bmatrix}6& 1 & 3 &\\ 8 & 6 & 1& \\ & ... &
FFT/NTT高精度
08-22
您好！对于使用FFT（快速傅里叶变换）或NTT（快速数论变换）来实现高精度计算，可以采以下步骤： 1. 将数字表示为多项式形式：将需要进行高精度计算的数字转换为多项式的系数形式。例如，将一个整数表示为多项式的系数。 2. 补零：为了进行FFT或NTT，需要将系数数组的长度扩展为2的幂次方。如果系数数组长度不是2的幂次方，可以在末尾添加零，使其长度满足要求。 3. 进行FFT或NTT：使用选择的变换算法（FFT或NTT）对系数数组进行变换。这将把多项式从系数域转换到点值域。 4. 点值运算：在点值域中，对应位置的点值可以通过乘法或加法来进行运算。例如，如果要计算乘法，则可以对应位置相乘；如果要计算加法，则对应位置相加。 5. 逆变换：对变换后的结果进行逆变换，将其从点值域转回系数域。 6. 去除高位零：在得到逆变换结果后，需要去除多余的高位零，以获得最终的高精度计算结果。 需要注意的是，FFT和NTT是用于加速多项式乘法的算法，通过将乘法复杂度从O(n^2)降低到O(nlogn)，其中n是多项式的长度。因此，它们通常与其他算法一起使用，如Karatsuba乘法或Schönhage-Strassen乘法，以实现高精度计算。

“相关推荐”对你有帮助么？

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

Greate AUK

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付：¥1
获取中
扫码支付

您的余额不足，请更换扫码支付或充值

打赏作者

实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值