• 求Q(x)模X^n + 1的余数多项式的FFT算法


    贴上IFFT算法的C语言代码:https://blog.csdn.net/great978/article/details/84306827


    Qleft ( x 
ight )X^{n}+1的余数多项式,最高次数不超过n,令其余数多项式的形式为fleft ( x 
ight )=a_{0}+a_{1}x+a_{2}x^{2}+......+a_{n-1}x^{n-1}

    1、因为X^{n}+1有n个解,对应余数多项式有n个,分别是left { fleft ( x_{j} 
ight ) 
ight }其中,0<=j<=n-1,而x_{_{j}}则是X^{n}+1的n个解,而x_{_{j}}是容易求的,

    x_{j} = e^{ileft ( pi /n
ight )left ( 2j+1 
ight )},但是如果每个都往left { fleft ( x_{j} 
ight ) 
ight }里面带入的话,会花费大量的计算量。

    这里提一种基于FFT的简单运算,这里要求n=2^{k}。(这里对为什么是基于FFT的我也不是很明白,有一说是

    left ( a_{0} ,a_{1},...a_{n-1}
ight )
ightarrow left ( fleft ( x_{0} 
ight ) ,fleft ( x_{1} 
ight ),....,fleft ( x_{n-1} 
ight )
ight )的多项式的离散傅里叶变换,用快速傅里叶变换FFT计算得来的,DFS的公式这儿就不列出了)

    2、计算算法,因为n=2^{k},把fleft ( x 
ight )=a_{0}+a_{1}x+a_{2}x^{2}+......+a_{n-1}x^{n-1}按奇偶项分开,每次分开一半,直到最后只剩两项的时候进行计算,然后一直递归向上,这样就把大量的乘幂运算转化为了加法运算,下面给出C语言的代码

    
    //
    
    #include "stdafx.h"
    
    #define _CRT_SECURE_NO_WARNINGS
    #include "stdlib.h"
    #include "math.h"
    #include "stdio.h"
    
    
    #define K 3
    
    #define Pi  3.1415927   //定义圆周率Pi
    #define LEN sizeof(struct Compx)  //定义复数结构体大小
    
    //-----定义复数结构体-----------------------
    struct Compx
    {
    	double real;
    	double imag;
    }Compx;
    
    //-----复数乘法运算函数---------------------
    struct Compx mult(struct Compx b1, struct Compx b2)
    {
    	struct Compx b3;
    	b3.real = b1.real*b2.real - b1.imag*b2.imag;
    	b3.imag = b1.real*b2.imag + b1.imag*b2.real;
    	return(b3);
    }
    
    //-----复数加法运算函数---------------------
    struct Compx add(struct Compx a, struct Compx b)
    {
    	struct Compx c;
    	c.real = a.real + b.real;
    	c.imag = a.imag + b.imag;
    	return(c);
    }
    
    struct Compx FFT(struct Compx *t , int n , struct Compx root , struct Compx result ) ;
    
    int main(  ) 
    {
    	int N,i;
    	N =8;
    
        int x[8];
        for( i = 0 ; i < N; i++ )
    	{
        	scanf("%d",&x[i]);
    	}
    	
        struct  Compx * Source = (struct Compx *)malloc(N*LEN);			//为结构体分配存储空间
        struct  Compx * Result = (struct Compx *)malloc(N*LEN);
        struct  Compx * Root = (struct Compx *)malloc(N*LEN);
    
        
        //初始化======================================= 
        printf("
    Source初始化:
    ");
    	for (i = 0; i<N; i++)
    	{
    		Source[i].real = x[i];
    		Source[i].imag = 0;
    		printf("%.4f ", Source[i].real);
    		printf("+%.4fj  ", Source[i].imag);
    		printf("
    ");
    	}
    	printf("
    Result初始化:
    ");
    	for (i = 0; i<N; i++)
    	{
    		Result[i].real = 0;
    		Result[i].imag = 0;
    		printf("%.4f ", Result[i].real);
    		printf("+%.4fj  ", Result[i].imag);
    		printf("
    ");
    	}
    	printf("
    Root初始化:
    ");
    	for (i = 0; i<N; i++)
    	{
    		Root[i].real = cos( 2 * Pi * ( 2 * i + 1) / (2 * N));
    		Root[i].imag = sin( 2 * Pi * ( 2 * i + 1) / (2 * N));
    		printf("%.4f ", Root[i].real);
    		printf("+%.4fj  ", Root[i].imag);
    		printf("
    ");
    	}
    	
    
        for( i = 0 ; i < N ; i++) 
        {
        	Result[i] = FFT(Source , N , Root[i] , Result[i]);
    	}
    	
        //结果表示
    	printf("
    Result计算结果:
    ");
    	for (i = 0; i<N; i++)
    	{
    		printf("%.4f ", Result[i].real);
    		printf("+%.4fj  ", Result[i].imag);
    		printf("
    ");
    	}
        
    	return 0;
    }
    
    struct Compx FFT(struct Compx *t , int n , struct Compx root , struct Compx result ) 
    {
    	int i,j;
    	struct  Compx * even = (struct Compx *)malloc((n / 2) * LEN);	
    	struct  Compx * odd = (struct Compx *)malloc((n / 2) * LEN);
    	
    	//划分奇偶项 
    	for (i = 0 , j = 0 ; i < n; i += 2 , j++ )
    	{
    		even[j].real = t[i].real;
    		even[j].imag = t[i].imag;
    	}
    	for (i = 1 , j = 0 ; i < n; i += 2 , j++)
    	{
    		odd[j].real = t[i].real;
    		odd[j].imag = t[i].imag;
    	}
    	
    	if(n == 2)
    	{
    		struct Compx s = add( even[0] , mult( root , odd[0]) );
    		return add(result , s);
    	}
    	else
    	{
    		return add( FFT(even , n / 2 , mult(root , root) , result ) , mult( root , FFT(odd , n / 2 , mult(root , root) , result ) ) );
    	}
    }
    

    2018年11月20日更新:代码中root根值计算有错误,上述代码中已修改 :

    新的代码中是:

    Root[i].real = cos( 2 * Pi * ( 2 * i + 1) / (2 * N));
    Root[i].imag = sin( 2 * Pi * ( 2 * i + 1) / (2 * N));

    之前代码是:

    Root[i].real = cos( Pi * ( 2 * i + 1) / (2 * N));
    Root[i].imag = sin( Pi * ( 2 * i + 1) / (2 * N));

    欢迎指正!

  • 相关阅读:
    Bellman-Ford 单源最短路径算法
    Prim 最小生成树算法
    Kruskal 最小生成树算法
    Kosaraju 算法检测有向图的强连通性
    Kosaraju 算法查找强连通分支
    不相交集合森林的启发式策略
    Union-Find 检测无向图有无环路算法
    redis的持久化方式RDB和AOF的区别
    Docker -v 对挂载的目录没有权限 Permission denied
    postgresql如何让主键自增
  • 原文地址:https://www.cnblogs.com/cg-bestwishes/p/10681152.html
Copyright © 2020-2023  润新知