精华内容
下载资源
问答
  • 有限域的多项式
    2020-12-22 01:41:26

    版权申明:本文为博主窗户(Colin Cai)原创,欢迎转帖。如要转贴,必须注明原文网址

    http://www.cnblogs.com/Colin-Cai/p/9489225.html

    作者:窗户

    QQ/微信:6679072

    E-mail:6679072@qq.com

    接着上两章内容,我们还是得继续寻找有限域的构造方法。上章证明矩阵环是个单环,自然是没戏了,但我们还可以考虑多项式环。

    多项式环

    多项式是我们大家熟知的概念,以下都是一元多项式:

    1

    2x+4

    x2+2x+3

    3x2+5x2+9

    ...

    所谓的一元就是只有一个未知数,在这里我就不对于一元多项式给出一个严格的定义了,直接解释多项式环。

    所谓一个环A的多项式环B,指的是如下:

    (1) B的每个元是一个一元多项式

    (2) B的每个元(一元多项式)的每一个系数都是A上的元

    (3) 系数全是A上的元的一元多项式都是B的元

    多项式的加法、减法就是合并同类项,因为系数取自一个环,所以系数间的加法、减法是合法的,会得到别的多项式。

    同样,多项式的乘法要麻烦一点,不过也是得到多项式的。多项式的乘法是利用分配律,展开各项。以下整系数多项式的例子可以让我们回忆起多项式的乘法:

    (x2+2x+3) * (3x+5)

    = (x2+2x+3) * 3x + (x2+2x+3) * 5

    = x2 * 3x + 2x * 3x + 3 * 3x + x2 * 5 + 2x * 5 + 3 *5

    = 3x3 + 6x2 + 9x + 5x2 + 10x + 15

    = 3x3 + 11x2 + 19x + 15

    因为系数都在一个环里,所有乘法、加法都是封闭的,所以多项式乘法也是一样合法的。

    所以多项式环当然是环。

    可能有的人想问,为什么这里非要用一元多项式?

    其实我们在刚才的多项式环定义那里为多项式引入任意多个未知数(甚至无穷多个未知数),其组成的代数系统依然为环,只是多元的多项式环挺复杂,这里不研究。

    不可分多项式

    我们知道质数是2、3、5、7、11、13、17...这样除了1和本身外没有其他正约数的正整数。

    我们在这里对质数做一个引申。

    一个多项式环上的任意多项式,当然可以表示为1和自身的乘积,当然也可以表示为-1(1元的相反元)和自身的相反元的乘积,这两者都是很平凡的。

    比如:

    x2+x+1 = 1 * (x2+x+1)

    = -1 * (-x2-x-1)

    这都是平凡的,没什么意义。

    如果是域上的多项式环,里面任何多项式表示成域上任何一个非0元和一个多项式的乘积。从而,这些也都是平凡的。

    而所谓真正意义上的分解,就是要求两个乘积项都不是常数,也就是次数是大于0的。

    比如,

    x2+2x+1 = (x+1) * (x+1)

    不可分解的多项式我们称之为不可分多项式。

    比如整数系数下的x2+x+1就是不可分多项式,实际上,即使是2元域(0/1两个元组成的特征2的域)上,这个多项式也是不可分多项式。

    但在7元域(0/1/2/3/4/5/6组成的特征7的域)上,

    x2+x+1 = (x+3) * (x+5)

    多项式的带余除法

    我们从小就知道自然数的带余除法,

    比如

    7÷3 = 2 ... 1

    换个写法,7 = 3*2+1

    其实,域的多项式环里的多项式也存在这样的带余除法。

    对于多项式f和g(g为非0多项式),一定存在唯一的多项式a和b,满足

    f = g*a+b

    并且b的次数小于g的次数。

    其实证明起来很简单,就如同儿时的竖式除法计算那样,一步步的把高次的项减掉。

    比如我们以2阶素域下的多项式 x5+x4+x+1 和 x2+x+1为例

    x3    + x      +1

    _________________________

    x2 + x + 1    |     x5 + x4                 + x    + 1

    x5 + x4 + x3

    ______________

    x3         + x

    x3 + x2 + x

    ___________

    x2          + 1

    x2   + x   + 1

    ____________

    x

    所以,

    x5+x4+x+1 = (x2+x+1) * (x3+x+1)  + x

    带余除法对于后面的理解有至关重要的作用。

    有限域

    既然想通过商环的方法构造域,那么当然要先考虑多项式环的理想。

    我们依然使用生成元的方法去研究。

    我们以 p阶素域作为原本的环 A, 那么A的多项式环称为 B,

    我们考虑由多项式 f 生成的理想,我们假设 f 是可以分解的,f = f1 * f2。

    f1、f2并不在理想里,很明显,f1、f2的次数比f都低,不存在f乘以一个多项式得到f1或 f2。

    我们再回忆一下商环的运算,根据f = f1 * f2,我们有

    商集(f1) * 商集(f2) = 商集(f)

    商集(f)其实就是理想,也就是商环里的0元,

    从而左边两个非0元乘法得到右边的0元,于是这个商环不是整环,当然更不可能是域了。

    于是我们考虑由单个不可分多项式生成的理想。

    考虑其下一个m次不可分多项式 f(最高未知数次数为m)生成的理想 C ,

    然后我们考虑商环B/C长什么样。

    理想C其实是所有以多项式 f为因子的多项式的集合。

    我们考虑所有次数小于m的多项式,根据排列组合的乘法原理,这样的多项式一共有pm个。

    对于任意两个不同的次数小于m的多项式,假设为g和h。

    g-h为非0的次数小于m的多项式,从而g-h不可能以f为因子,从而g-h不在理想里,从而g和h一定属于不同的商集。

    因为g和h选择的随意性,从而这pm个多项式分属于pm个不同的商集。

    上面介绍过带余除法,考虑次数大于等于m的多项式,假设有一个这样的多项式h,

    一定存在一个多项式a和一个次数小于m的多项式b,使得

    h = g*a+b

    h-b = g*a

    g*a在理想C里,于是h和b在同一个商集里。

    由于h选择的随意性,从而任何一个次数大于等于m的多项式都落在那pm个不同的商集里。

    所以,我们最终的这个商环也就有pm个元。

    这里多项式乘法的可交换性遗传自域乘法的可交换,从而这个商环可交换是必然的。

    另外,f的不可分特性导致了如果任意g、h不以f为因子,则g*h也不以f为因子。从而,这个商环是一个整环。

    有限的可交换整环,因为其有限性,那么当然是除环,从而当然就是域啦(其实,并不存在有限的不可交换整环,不过这个定理证明有那么点麻烦)。

    OK,我们终于找到了构造任意阶有限域的方法。

    我们可以用这pm个次数小于m的多项式来代表这个域的各个元素。加法、减法就是合并同类项。

    乘法就是多项式乘法结果再利用带余除法除以多项式 f 得到的余数。

    我们来举个例子。

    x2+x+1 是 2阶素域下的不可分多项式。

    利用刚才的手段得到了一个4阶域,我们可以记该域下的4个元为

    [0]

    [1]

    [x]

    [x+1]

    其四则运算为

    [0] + [0] = [0]  [0] - [0] = [0]  [0] * [0] = [0]

    [0] + [1] = [1]  [0] - [1] = [1]  [0] * [1] = [0]  [0] / [1] = [0]

    [0] + [x] = [x]  [0] - [x] = [x]  [0] * [x] = [0]  [0] / [x] = [0]

    [0] + [x+1] = [x+1]  [0] - [x+1] = [x+1]  [0] * [x+1] = [0]  [0] / [x+1] = [0]

    [1] + [0] = [1]  [1] - [0] = [1]  [1] * [0] = [0]

    [1] + [1] = [0]  [1] - [1] = [0]  [1] * [1] = [1]  [1] / [1] = [1]

    [1] + [x] = [x+1]  [1] - [x] = [x+1]  [1] * [x] = [x]  [1] / [x] = [x+1]

    [1] + [x+1] = [x]  [1] - [x+1] = [x]  [1] * [x+1] = [x+1]  [1] / [x+1] = [x]

    [x] + [0] = [x]  [x] - [0] = [x]  [x] * [0] = [0]

    [x] + [1] = [x+1]  [x] - [1] = [x+1]  [x] * [1] = [x]  [x] / [1] = [x]

    [x] + [x] = [0]  [x] - [x] = [0]  [x] * [x] = [x+1]  [x] / [x] = [1]

    [x] + [x+1] = [1]  [x] - [x+1] = [1]  [x] * [x+1] = [1]  [x] / [x+1] = [x+1]

    [x+1] + [0] = [x+1]  [x+1] - [0] = [x+1]  [x+1] * [0] = [0]

    [x+1] + [1] = [x]  [x+1] - [1] = [x]  [x+1] * [1] = [x+1]  [x+1] / [1] = [x+1]

    [x+1] + [x] = [1]  [x+1] - [x] = [1]  [x+1] * [x] = [1]  [x+1] / [x] = [x]

    [x+1] + [x+1] = [0]  [x+1] - [x+1] = [0]  [x+1] * [x+1] = [x]  [x+1] / [x+1] = [1]

    附:上一章有网友提议用LaTeX,嗯,本人确实有点懒,不过后面会考虑的。

    另外,很多网上文章这里都写本原多项式,人云亦云,悲哀,一帮只愿去抄书不愿去理解的人啊。

    建议还是先去了解一下什么叫本原多项式吧。

    更多相关内容
  • 多项式环R=GF[pd][x]/(xn−1)R=GF[p^d][x]/(x^n-1)R=GF[pd][x]/(xn−1),其中有限域GF(pd)≅Zp[x]/(m(x))GF(p^d) \cong Z_p[x]/(m(x))GF(pd)≅Zp​[x]/(m(x)),m(x)m(x)m(x)是ZpZ_pZp​上ddd次不可约多项式。...

    问题

    多项式环 R = G F [ p d ] [ x ] / ( x n − 1 ) R=GF[p^d][x]/(x^n-1) R=GF[pd][x]/(xn1),其中有限域 G F ( p d ) ≅ Z p [ x ] / ( m ( x ) ) GF(p^d) \cong Z_p[x]/(m(x)) GF(pd)Zp[x]/(m(x)) m ( x ) m(x) m(x) Z p Z_p Zp d d d次不可约多项式。对于任意 f , g ∈ R f,g \in R f,gR,如何快速计算 f ⋅ g f\cdot g fg

    NTL

    代码

    #include <iostream>
    #include "tools.h"
    
    #include <NTL/ZZ_p.h> // integers mod p
    #include <NTL/ZZ_pX.h> // polynomials over ZZ_p
    #include <NTL/ZZ_pE.h> // ring/field extension of ZZ_p
    #include <NTL/ZZ_pEX.h> // polynomials over ZZ_pE
    #include <NTL/ZZ_pXFactoring.h>
    #include <NTL/ZZ_pEXFactoring.h>
    
    using namespace std;
    using namespace NTL;
    
    #pragma comment(lib, "NTL")
    
    int main()
    {
    	ZZ p(2); //初始化为2
    
    	//群Z_p
    	ZZ_p::init(p);
    
    	//随机生成Z_p[x]中的d次不可约多项式
    	int d = 149;
    	ZZ_pX m;
    	//BuildIrred(m, d);
    	SetCoeff(m, d);
    	SetCoeff(m, 10);
    	SetCoeff(m, 9);
    	SetCoeff(m, 7);
    	SetCoeff(m, 0);
    
    	//域GF(p^d) = Z_p[x]/m(x)
    	ZZ_pE::init(m);
    
    	//GF(p^d)[x]中的多项式
    	ZZ_pEX f, g, h, P;
    	int n = 53;
    
    	//P(x)=x^n-1
    	SetCoeff(P, 0, -1);
    	SetCoeff(P, n);
    
    	// f(x) = x^8 - 1
    	SetCoeff(f, 8); //将 x^8 系数置为 1
    	SetCoeff(f, 0, -1); //将 x^0 系数置为 -1
    
    	//随机生成n长多项式g(x)
    	random(g, n-1);
    
    	// 环上多项式的运算
    	h = g * g;
    
    	cout << "p = " << p << endl;
    	cout << "d = " << d << endl;
    	cout << "n = " << n << endl;
    	cout << "m(x) = " << m << endl;
    	cout << "P(x) = " << P << endl;
    
    	//cout << "f = " << f << endl;
    	//cout << "g = " << g << endl;
    	//cout << "h = " << h << endl;
    
    	int loop = 10;
    
    	printf("%d rounds, NTL mul: ", loop);
    	Loop(10, f*g);
    
    	printf("%d rounds, NTL mod: ", loop);
    	Loop(10, h%P);
    
    
    	return 0;
    }
    

    结果

    p = 2
    d = 149
    n = 53
    m(x) = [1 0 0 0 0 0 0 1 0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1]
    P(x) = [[1] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [] [1]]
    10 rounds, NTL mul: 5.20162 s
    10 rounds, NTL mod: 17.3201 s
    

    NTT+Barrett

    代码

    IncompleteNTT
    #include "IncompleteNTT.h"
    
    
    int32 brv(int32 b, int32 l)
    {
    	int32 bb = 0;
    
    	for (int32 i = 0; i < l; i++)
    	{
    		bb <<= 1;
    		bb |= (b & 1);
    		b >>= 1;
    	}
    
    	return bb;
    }
    
    //int64 exgcd(int64* x, int64* y, int64 a, int64 b)
    //{
    //	if (b == 0)
    //	{
    //		*x = 1;
    //		*y = 0;
    //		return a;
    //	}
    //	int64 ret = exgcd(x, y, b, a%b);
    //	int64 tmp = *x;
    //	*x = *y;
    //	*y = tmp - (a / b) * (*y);
    //	return ret;
    //}
    
    //int32 deg(Poly& f)
    //{
    //	int32 d = f.size() - 1;
    //	int64* pf = f.data();
    //
    //	while (pf[d] == 0)
    //		d--;
    //	return d;
    //}
    
    ostream& operator<<(ostream& cout, Poly& f)
    {
    	int32 d = deg(f);
    	int64* pf = f.data();
    
    	printf("[ ");
    	for (int32 i = 0; i < d; i++)
    		printf("%lld, ", pf[i]);
    	printf("%lld ]", pf[d]);
    
    	return cout;
    }
    
    ostream& print_vector(int64* f, int32 len)
    {
    	len--;
    
    	printf("[ ");
    	for (int32 i = 0; i < len; i++)
    		printf("%lld, ", f[i]);
    	printf("%lld ]", f[len]);
    
    	return std::cout;
    }
    
    
    //############################# 分割线 ###############################
    
    
    void IncompleteNTT::init(int64 p, int64 w, int32 k, int32 h)
    {
    	this->p = p;
    	this->w = w;
    	this->k = k;
    	this->h = h;
    	this->m = (1L << k);
    	this->n = this->h * this->m;
    
    	int64 minv, pinv;
    	int64 gcd = exgcd(&minv, &pinv, m, p);
    	if (gcd != 1)
    		ErrorInfo("%s\n", "gcd(2^k, p) != 1");
    	this->factor = minv < 0 ? minv + p : minv;
    	this->factor_pre = PreCompute(this->factor, p);
    	this->one_pre = PreComputeMod(p);
    
    	this->W.resize(this->m + 1); //重复存储:w^{0} = w^{m}
    	this->W_pre.resize(this->m + 1);
    
    	int64 *pW = W.data();
    	int64 *pW_pre = W_pre.data();
    
    	int64 wi = 1;	//单位根的幂次
    	for (int32 i = 0; i <= this->m; i++)
    	{
    		pW[i] = wi;
    		pW_pre[i] = PreCompute(wi, p);
    		wi = (wi * w) % p;
    	}
    
    	brvlist.resize(k + 1);
    	for (int j = 0; j <= k; j++)
    	{
    		auto& list = brvlist[j];
    		int32 len = (1L << j);
    		for (int i = 0; i < len; i++)
    			list.push_back(brv(i, j));
    	}
    }
    
    
    
    
    void IncompleteNTT::FwdNTT(NTTRep& F, Poly& f)
    {
    	int32 d = deg(f);
    	if (d >= this->n)
    		ErrorInfo("%s\n", "deg(f) > n - 1");
    
    	if (F.size() < this->n)
    		F.resize(this->n);
    
    	int64* pF = F.data();
    	int64* pf = f.data();
    
    	FwdNTT(pF, pf, f.size());
    }
    
    void IncompleteNTT::FwdNTT(int64* pF, int64* pf, int32 flen)
    {
    	if (pF != pf)
    	{
    		int32 nn = min(flen, n);
    		for (int32 i = 0; i < nn; i++)
    			pF[i] = BarrettMod(pf[i], one_pre, p);
    	}
    
    	int64* pW = W.data();
    	int64* pW_pre = W_pre.data();
    
    	int64 P = 2 * p;
    	for (int32 j = 0; j < k; j++)
    	{
    		int32 blocknum = 1 << j; //对j层分块
    		int32 blocksize = this->n >> j; //对j层块大小
    		int32 offset = blocksize >> 1; //计算第j+1层时蝴蝶的偏移量(块大小的一半)
    		int32* list = brvlist[j].data();
    
    		if (offset >= 4 * h)
    			for (int32 i = 0; i < blocknum; i++)
    			{
    				int64 X, WY;
    				int32 s = i * blocksize; //块的起始位置
    				int64 ww = pW[list[i] << (k - j - 1)]; //第j层第i个分块使用的单位根,w_{2^{j+1}}^{brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[list[i] << (k - j - 1)];
    
    				int64* pX = pF + (s);
    				int64* pY = pF + (s + offset);
    
    				for (int32 k = 0; k < offset; k += 4)
    				{
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    				}
    			}
    		else if (offset == 2 * h)
    			for (int32 i = 0; i < blocknum; i++)
    			{
    				int64 X, WY;
    				int32 s = i * blocksize; //块的起始位置
    				int64 ww = pW[list[i] << (k - j - 1)]; //第j层第i个分块使用的单位根,w_{2^{j+1}}^{brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[list[i] << (k - j - 1)];
    
    				int64* pX = pF + (s);
    				int64* pY = pF + (s + offset);
    
    				for (int32 k = 0; k < h; k++)
    				{
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    				}
    			}
    		else
    			for (int32 i = 0; i < blocknum; i++)
    			{
    				int64 X, WY;
    				int32 s = i * blocksize; //块的起始位置
    				int64 ww = pW[list[i] << (k - j - 1)]; //第j层第i个分块使用的单位根,w_{2^{j+1}}^{brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[list[i] << (k - j - 1)];
    
    				int64* pX = pF + (s);
    				int64* pY = pF + (s + offset);
    
    				for (int32 k = 0; k < h; k++)
    				{
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    				}
    			}
    	}
    
    	//迭代结束,结果是m=2^k个h长多项式
    	//约束范围到[0,2p)
    	if (k >= 2)
    		for (int32 i = 0; i < n; i += 4)
    		{
    			*pF -= ((P - *pF) >> 63)&P;
    			pF++;
    
    			*pF -= ((P - *pF) >> 63)&P;
    			pF++;
    
    			*pF -= ((P - *pF) >> 63)&P;
    			pF++;
    
    			*pF -= ((P - *pF) >> 63)&P;
    			pF++;
    		}
    	else
    		for (int32 i = 0; i < n; i++)
    		{
    			*pF -= ((P - *pF) >> 63)&P;
    			pF++;
    		}
    }
    
    
    
    void IncompleteNTT::InvNTT(Poly& f, NTTRep& F)
    {
    	int32 d = deg(F);
    	if (d >= this->n)
    		ErrorInfo("%s\n", "deg(F) > n - 1");
    
    	if (f.size() < this->n)
    		f.resize(this->n);
    
    	int64* pF = F.data();
    	int64* pf = f.data();
    
    	InvNTT(pf, pF);
    }
    
    void IncompleteNTT::InvNTT(int64* pf, int64* pF)
    {
    	if (pF != pf)
    	{
    		for (int32 i = 0; i < n; i++)
    			pf[i] = pF[i];
    	}
    
    	int64* pW = W.data();
    	int64* pW_pre = W_pre.data();
    
    	int64 P = 2 * p;
    	for (int32 j = k - 1; j >= 0; j--)
    	{
    		int32 blocknum = 1 << j; //对j+1层分块
    		int32 blocksize = this->n >> j; //第j+1层块大小
    		int32 offset = blocksize >> 1; //计算第j层时蝴蝶的偏移量(块大小的一半)
    		int32* list = brvlist[j].data();
    
    		if (offset >= 4 * h)
    			for (int32 i = 0; i < blocknum; i++)
    			{
    				int64 XY, T;
    				int32 s = i * blocksize; //块的起始位置
    				int64 ww = pW[m - (list[i] << (k - j - 1))]; //第j+1层第i个分块使用的单位根,w_{2^{j+1}}^{-brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[m - (list[i] << (k - j - 1))];
    
    				int64* pX = pf + (s);
    				int64* pY = pf + (s + offset);
    
    				for (int32 k = 0; k < offset; k += 4)
    				{
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    				}
    			}
    		else if (offset == 2 * h)
    			for (int32 i = 0; i < blocknum; i++)
    			{
    				int64 XY, T;
    				int32 s = i * blocksize; //块的起始位置
    				int64 ww = pW[m - (list[i] << (k - j - 1))]; //第j+1层第i个分块使用的单位根,w_{2^{j+1}}^{-brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[m - (list[i] << (k - j - 1))];
    
    				int64* pX = pf + (s);
    				int64* pY = pf + (s + offset);
    
    				for (int32 k = 0; k < h; k++)
    				{
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    				}
    			}
    		else
    			for (int32 i = 0; i < blocknum; i++)
    			{
    				int64 XY, T;
    				int32 s = i * blocksize; //块的起始位置
    				int64 ww = pW[m - (list[i] << (k - j - 1))]; //第j+1层第i个分块使用的单位根,w_{2^{j+1}}^{-brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[m - (list[i] << (k - j - 1))];
    
    				int64* pX = pf + (s);
    				int64* pY = pf + (s + offset);
    
    				for (int32 k = 0; k < h; k++)
    				{
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    				}
    			}
    	}
    
    	//迭代结束,再整理一下
    	int64 tmp;
    	if (k >= 2)
    		for (int32 i = 0; i < n; i += 4)
    		{
    			tmp = BarrettMulMod(*pf, factor, factor_pre, p);
    			*pf = tmp - (tmp >= p)*p;
    			pf++;
    
    			tmp = BarrettMulMod(*pf, factor, factor_pre, p);
    			*pf = tmp - (tmp >= p)*p;
    			pf++;
    
    			tmp = BarrettMulMod(*pf, factor, factor_pre, p);
    			*pf = tmp - (tmp >= p)*p;
    			pf++;
    
    			tmp = BarrettMulMod(*pf, factor, factor_pre, p);
    			*pf = tmp - (tmp >= p)*p;
    			pf++;
    		}
    	else
    		for (int32 i = 0; i < n; i++)
    		{
    			tmp = BarrettMulMod(*pf, factor, factor_pre, p);
    			*pf = tmp - (tmp >= p)*p;
    			pf++;
    		}
    }
    
    
    void IncompleteNTT::convolution(int64* FG, int64* F, int64* G, int64 w)
    {
    	for (int32 i = 0; i < h; i++)
    	{
    		FG[i] = 0;
    		int32 s = h + i;
    		for (int32 j = 0; j <= i; j++)
    			FG[i] = BarrettMod(FG[i] + G[j] * F[i - j], one_pre, p); //冗余,[0,2p)
    		for (int32 j = i + 1; j < h; j++)
    		{
    			int64 tmp = FG[i] + BarrettMod(w * G[j], one_pre, p) * F[s - j];
    			FG[i] = BarrettMod(tmp, one_pre, p); //冗余,[0,2p)
    		}
    	}
    }
    
    
    void IncompleteNTT::NTTMul(NTTRep& FG, NTTRep& F, NTTRep& G)
    {
    	if (F.size() < n || G.size() < n)
    		ErrorInfo("%s\n", "F.size < n or G.size < n");
    
    	if (FG.size() < n)
    		FG.resize(n);
    
    	int64* pFG = FG.data();
    	int64* pF = F.data();
    	int64* pG = G.data();
    
    	NTTMul(pFG, pF, pG);
    }
    
    void IncompleteNTT::NTTMul(int64* pFG, int64* pF, int64* pG)
    {
    	int64* pW = W.data();
    
    	int64 ww;
    	int64 r0, r1;
    	int64 P = 2 * p;
    	int32* list = brvlist[k].data();
    	if (h == 1)
    		if (k >= 2)
    			for (int32 i = 0; i < m; i += 4) //m个h长小多项式
    			{
    				//使用冗余表示,[0,2p)
    				*pFG = BarrettMod(*pF * *pG, one_pre, p);
    
    				pFG++;
    				pF++;
    				pG++;
    
    				//使用冗余表示,[0,2p)
    				*pFG = BarrettMod(*pF * *pG, one_pre, p);
    
    				pFG++;
    				pF++;
    				pG++;
    
    				//使用冗余表示,[0,2p)
    				*pFG = BarrettMod(*pF * *pG, one_pre, p);
    
    				pFG++;
    				pF++;
    				pG++;
    
    				//使用冗余表示,[0,2p)
    				*pFG = BarrettMod(*pF * *pG, one_pre, p);
    
    				pFG++;
    				pF++;
    				pG++;
    			}
    		else
    			for (int32 i = 0; i < m; i++) //m个h长小多项式
    			{
    				//使用冗余表示,[0,2p)
    				*pFG = BarrettMod(*pF * *pG, one_pre, p);
    
    				pFG++;
    				pF++;
    				pG++;
    			}
    	else if (h == 2)
    		if (k >= 2)
    			for (int32 i = 0; i < m; i += 4) //m个h长小多项式
    			{
    				int64 tmp;
    
    				//使用冗余表示,[0,2p)
    				ww = pW[list[i]]; //第k层第i个多项式使用的单位根,w_{2^k}^{brv_k(i)}
    
    				tmp = ww * BarrettMod(pF[1] * pG[1], one_pre, p);
    				r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(tmp, one_pre, p);
    				//r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(ww * BarrettMod(pF[1] * pG[1], one_pre, p), one_pre, p);
    				r1 = BarrettMod(pF[0] * pG[1], one_pre, p) + BarrettMod(pF[1] * pG[0], one_pre, p);
    				pFG[0] = r0 - (((P - r0) >> 63)&P);
    				pFG[1] = r1 - (((P - r1) >> 63)&P);
    
    				pFG += 2;
    				pF += 2;
    				pG += 2;
    
    				//使用冗余表示,[0,2p)
    				ww = pW[list[i + 1]]; //第k层第i个多项式使用的单位根,w_{2^k}^{brv_k(i)}
    
    				tmp = ww * BarrettMod(pF[1] * pG[1], one_pre, p);
    				r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(tmp, one_pre, p);
    				//r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(ww * BarrettMod(pF[1] * pG[1], one_pre, p), one_pre, p);
    				r1 = BarrettMod(pF[0] * pG[1], one_pre, p) + BarrettMod(pF[1] * pG[0], one_pre, p);
    				pFG[0] = r0 - (((P - r0) >> 63)&P);
    				pFG[1] = r1 - (((P - r1) >> 63)&P);
    
    				pFG += 2;
    				pF += 2;
    				pG += 2;
    
    				//使用冗余表示,[0,2p)
    				ww = pW[list[i + 2]]; //第k层第i个多项式使用的单位根,w_{2^k}^{brv_k(i)}
    
    				tmp = ww * BarrettMod(pF[1] * pG[1], one_pre, p);
    				r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(tmp, one_pre, p);
    				//r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(ww * BarrettMod(pF[1] * pG[1], one_pre, p), one_pre, p);
    				r1 = BarrettMod(pF[0] * pG[1], one_pre, p) + BarrettMod(pF[1] * pG[0], one_pre, p);
    				pFG[0] = r0 - (((P - r0) >> 63)&P);
    				pFG[1] = r1 - (((P - r1) >> 63)&P);
    
    				pFG += 2;
    				pF += 2;
    				pG += 2;
    
    				//使用冗余表示,[0,2p)
    				ww = pW[list[i + 3]]; //第k层第i个多项式使用的单位根,w_{2^k}^{brv_k(i)}
    
    				tmp = ww * BarrettMod(pF[1] * pG[1], one_pre, p);
    				r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(tmp, one_pre, p);
    				//r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(ww * BarrettMod(pF[1] * pG[1], one_pre, p), one_pre, p);
    				r1 = BarrettMod(pF[0] * pG[1], one_pre, p) + BarrettMod(pF[1] * pG[0], one_pre, p);
    				pFG[0] = r0 - (((P - r0) >> 63)&P);
    				pFG[1] = r1 - (((P - r1) >> 63)&P);
    
    				pFG += 2;
    				pF += 2;
    				pG += 2;
    			}
    		else
    			for (int32 i = 0; i < m; i++) //m个h长小多项式
    			{
    				//使用冗余表示,[0,2p)
    				ww = pW[list[i]]; //第k层第i个多项式使用的单位根,w_{2^k}^{brv_k(i)}
    
    				int64 tmp = ww * BarrettMod(pF[1] * pG[1], one_pre, p);
    				r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(tmp, one_pre, p);
    				//r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(ww * BarrettMod(pF[1] * pG[1], one_pre, p), one_pre, p);
    				r1 = BarrettMod(pF[0] * pG[1], one_pre, p) + BarrettMod(pF[1] * pG[0], one_pre, p);
    				pFG[0] = r0 - (((P - r0) >> 63)&P);
    				pFG[1] = r1 - (((P - r1) >> 63)&P);
    
    				pFG += 2;
    				pF += 2;
    				pG += 2;
    			}
    	else
    		for (int32 i = 0; i < m; i++) //m个h长小多项式
    		{
    			ww = pW[list[i]]; //第k层第i个多项式使用的单位根,w_{2^k}^{brv_k(i)}
    			convolution(pFG, pF, pG, ww);
    
    			pFG += h;
    			pF += h;
    			pG += h;
    		}
    }
    
    void IncompleteNTT::Param()
    {
    	cout << "p = " << p; pn;
    	cout << "w = " << w; pn;
    	cout << "k = " << k; pn;
    	cout << "h = " << h; pn;
    	cout << "n = " << n; pn;
    	cout << "factor = " << factor; pn;
    }
    
    
    //############################# 分割线 ###############################
    
    void NegacyclicNTT::init(int64 p, int64 w, int32 k, int32 h, bool tag)
    {
    	if (tag == 0)
    		IncompleteNTT::init(p, w, k, 2 * h); //构建父类IncompleteNTT,并初始化
    	else
    		IncompleteNTT::init(p, w, k + 1, h); //构建父类IncompleteNTT,并初始化
    
    	hh = h;
    	hm = m >> 1;
    	hn = n >> 1;
    
    	int64 hminv, pinv;
    	int64 gcd = exgcd(&hminv, &pinv, hm, p);
    	if (gcd != 1)
    		ErrorInfo("%s\n", "gcd(2^{k-1}, p) != 1");
    	this->factor = hminv < 0 ? hminv + p : hminv;
    	this->factor_pre = PreCompute(this->factor, p);
    }
    
    
    void NegacyclicNTT::FwdNTT(NTTRep& F, Poly& f)
    {
    	int32 d = deg(f);
    	if (d >= this->hn)
    		ErrorInfo("%s\n", "deg(f) > hn - 1");
    
    	if (F.size() < this->hn)
    		F.resize(this->hn);
    
    	int64* pF = F.data();
    	int64* pf = f.data();
    
    	FwdNTT(pF, pf, f.size());
    }
    
    void NegacyclicNTT::FwdNTT(int64* pF, int64* pf, int32 flen)
    {
    	if (pF != pf)
    	{
    		int32 nn = min(flen, hn);
    		for (int32 i = 0; i < nn; i++)
    			pF[i] = BarrettMod(pf[i], one_pre, p);
    	}
    
    	int64* pW = W.data();
    	int64* pW_pre = W_pre.data();
    
    	int64 P = 2 * p;
    	for (int32 j = 1; j < k; j++) //从j=1层开始迭代
    	{
    		int32 blocknum = 1 << j; //对j层分块
    		int32 blocksize = this->n >> j; //第j层块大小
    		int32 offset = blocksize >> 1; //计算第j+1层时蝴蝶的偏移量(块大小的一半)
    		int32* list = brvlist[j].data();
    
    		if (offset >= 4 * hh)
    			for (int32 i = blocknum >> 1; i < blocknum; i++) //只处理下半部分
    			{
    				int64 X, WY;
    				int32 s = (i - (blocknum >> 1)) * blocksize; //块的起始位置
    				int64 ww = pW[list[i] << (k - j - 1)]; //第j层第i个分块使用的单位根,w_{2^{j+1}}^{brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[list[i] << (k - j - 1)];
    
    				int64* pX = pF + (s);
    				int64* pY = pF + (s + offset);
    
    				for (int32 k = 0; k < offset; k += 4)
    				{
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    				}
    			}
    		else if (offset == 2 * hh)
    			for (int32 i = blocknum >> 1; i < blocknum; i++) //只处理下半部分
    			{
    				int64 X, WY;
    				int32 s = (i - (blocknum >> 1)) * blocksize; //块的起始位置
    				int64 ww = pW[list[i] << (k - j - 1)]; //第j层第i个分块使用的单位根,w_{2^{j+1}}^{brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[list[i] << (k - j - 1)];
    
    				int64* pX = pF + (s);
    				int64* pY = pF + (s + offset);
    
    				for (int32 k = 0; k < hh; k++)
    				{
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    				}
    			}
    		else
    			for (int32 i = blocknum >> 1; i < blocknum; i++) //只处理下半部分
    			{
    				int64 X, WY;
    				int32 s = (i - (blocknum >> 1)) * blocksize; //块的起始位置
    				int64 ww = pW[list[i] << (k - j - 1)]; //第j层第i个分块使用的单位根,w_{2^{j+1}}^{brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[list[i] << (k - j - 1)];
    
    				int64* pX = pF + (s);
    				int64* pY = pF + (s + offset);
    
    				for (int32 k = 0; k < hh; k++)
    				{
    					//CT蝴蝶(Harvey,输入输出范围[0,4p))
    					X = *pX;
    					X -= ((P - X) >> 63)&P;
    					WY = BarrettMulMod(*pY, ww, ww_pre, p); //Barrett算法
    					*pX = X + WY;
    					*pY = X + (P - WY);
    					pX++;
    					pY++;
    				}
    			}
    	}
    
    	//迭代结束,结果是hm=2^{k-1}个h长多项式
    	//约束范围到[0,2p)
    	if (k >= 2)
    		for (int32 i = 0; i < hn; i += 4)
    		{
    			*pF -= ((P - *pF) >> 63)&P;
    			pF++;
    
    			*pF -= ((P - *pF) >> 63)&P;
    			pF++;
    
    			*pF -= ((P - *pF) >> 63)&P;
    			pF++;
    
    			*pF -= ((P - *pF) >> 63)&P;
    			pF++;
    		}
    	else
    		for (int32 i = 0; i < hn; i++)
    		{
    			*pF -= ((P - *pF) >> 63)&P;
    			pF++;
    		}
    }
    
    
    
    void NegacyclicNTT::InvNTT(Poly& f, NTTRep& F)
    {
    	int32 d = deg(F);
    	if (d >= this->hn)
    		ErrorInfo("%s\n", "deg(F) > hn - 1");
    
    	if (f.size() < this->hn)
    		f.resize(this->hn);
    
    	int64* pF = F.data();
    	int64* pf = f.data();
    
    	InvNTT(pf, pF);
    }
    
    void NegacyclicNTT::InvNTT(int64* pf, int64* pF)
    {
    	if (pF != pf)
    	{
    		for (int32 i = 0; i < hn; i++)
    			pf[i] = pF[i];
    	}
    
    	int64* pW = W.data();
    	int64* pW_pre = W_pre.data();
    
    	int64 P = 2 * p;
    	for (int32 j = k - 1; j >= 1; j--) //回到第j=1层
    	{
    		int32 blocknum = 1 << j; //对j+1层分块
    		int32 blocksize = this->n >> j; //第j+1层块大小
    		int32 offset = blocksize >> 1; //计算第j层时蝴蝶的偏移量(块大小的一半)
    		int32* list = brvlist[j].data();
    
    		if (offset >= 4 * hh)
    			for (int32 i = blocknum >> 1; i < blocknum; i++) //只处理下半部分
    			{
    				int64 XY, T;
    				int32 s = (i - (blocknum >> 1)) * blocksize; //块的起始位置
    				int64 ww = pW[m - (list[i] << (k - j - 1))]; //第j+1层第i个分块使用的单位根,w_{2^{j+1}}^{-brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[m - (list[i] << (k - j - 1))];
    
    				int64* pX = pf + (s);
    				int64* pY = pf + (s + offset);
    
    				for (int32 k = 0; k < offset; k += 4)
    				{
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    				}
    			}
    		else if (offset == 2 * hh)
    			for (int32 i = blocknum >> 1; i < blocknum; i++) //只处理下半部分
    			{
    				int64 XY, T;
    				int32 s = (i - (blocknum >> 1)) * blocksize; //块的起始位置
    				int64 ww = pW[m - (list[i] << (k - j - 1))]; //第j+1层第i个分块使用的单位根,w_{2^{j+1}}^{-brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[m - (list[i] << (k - j - 1))];
    
    				int64* pX = pf + (s);
    				int64* pY = pf + (s + offset);
    
    				for (int32 k = 0; k < hh; k++)
    				{
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    				}
    			}
    		else
    			for (int32 i = blocknum >> 1; i < blocknum; i++) //只处理下半部分
    			{
    				int64 XY, T;
    				int32 s = (i - (blocknum >> 1)) * blocksize; //块的起始位置
    				int64 ww = pW[m - (list[i] << (k - j - 1))]; //第j+1层第i个分块使用的单位根,w_{2^{j+1}}^{-brv_{j+1}(2i)}
    				int64 ww_pre = pW_pre[m - (list[i] << (k - j - 1))];
    
    				int64* pX = pf + (s);
    				int64* pY = pf + (s + offset);
    				for (int32 k = 0; k < hh; k++)
    				{
    					//GS蝴蝶(Harvey,输入输出范围[0,2p))
    					XY = *pX + *pY;
    					XY -= ((P - XY) >> 63)&P;
    					T = *pX + (P - *pY);
    					*pX = XY;
    					*pY = BarrettMulMod(T, ww, ww_pre, p); //Barrett算法
    					pX++;
    					pY++;
    				}
    			}
    	}
    
    	//迭代结束,再整理一下
    	int64 tmp;
    	if (k >= 2)
    		for (int32 i = 0; i < hn; i += 4)
    		{
    			tmp = BarrettMulMod(*pf, factor, factor_pre, p);
    			*pf = tmp - (tmp >= p)*p;
    			pf++;
    
    			tmp = BarrettMulMod(*pf, factor, factor_pre, p);
    			*pf = tmp - (tmp >= p)*p;
    			pf++;
    
    			tmp = BarrettMulMod(*pf, factor, factor_pre, p);
    			*pf = tmp - (tmp >= p)*p;
    			pf++;
    
    			tmp = BarrettMulMod(*pf, factor, factor_pre, p);
    			*pf = tmp - (tmp >= p)*p;
    			pf++;
    		}
    	else
    		for (int32 i = 0; i < hn; i++)
    		{
    			tmp = BarrettMulMod(*pf, factor, factor_pre, p);
    			*pf = tmp - (tmp >= p)*p;
    			pf++;
    		}
    }
    
    
    
    void NegacyclicNTT::NTTMul(NTTRep& FG, NTTRep& F, NTTRep& G)
    {
    	if (F.size() < hn || G.size() < hn)
    		ErrorInfo("%s\n", "F.size < n or G.size < n");
    
    	if (FG.size() < hn)
    		FG.resize(hn);
    
    	int64* pFG = FG.data();
    	int64* pF = F.data();
    	int64* pG = G.data();
    
    	NTTMul(pFG, pF, pG);
    }
    
    void NegacyclicNTT::NTTMul(int64* pFG, int64* pF, int64* pG)
    {
    	int64* pW = W.data();
    
    	int64 ww;
    	int64 r0, r1;
    	int64 P = 2 * p;
    	int32* list = brvlist[k].data() + hm;
    	if (h == 1)
    		if (k >= 2)
    			for (int32 i = 0; i < hm; i += 4) //m个h长小多项式
    			{
    				//使用冗余表示,[0,2p)
    				*pFG = BarrettMod(*pF * *pG, one_pre, p);
    
    				pFG++;
    				pF++;
    				pG++;
    
    				//使用冗余表示,[0,2p)
    				*pFG = BarrettMod(*pF * *pG, one_pre, p);
    
    				pFG++;
    				pF++;
    				pG++;
    
    				//使用冗余表示,[0,2p)
    				*pFG = BarrettMod(*pF * *pG, one_pre, p);
    
    				pFG++;
    				pF++;
    				pG++;
    
    				//使用冗余表示,[0,2p)
    				*pFG = BarrettMod(*pF * *pG, one_pre, p);
    
    				pFG++;
    				pF++;
    				pG++;
    			}
    		else
    			for (int32 i = 0; i < hm; i++) //m个h长小多项式
    			{
    				//使用冗余表示,[0,2p)
    				*pFG = BarrettMod(*pF * *pG, one_pre, p);
    
    				pFG++;
    				pF++;
    				pG++;
    			}
    	else if (h == 2)
    		if (k >= 3)
    			for (int32 i = 0; i < hm; i += 4) //hm个h长小多项式
    			{
    				int64 tmp;
    
    				//使用冗余表示,[0,2p)
    				ww = pW[list[i]]; //第k层第2^{k-1}+i个多项式使用的单位根,w_{2^k}^{brv_k(2^{k-1}+i)}
    
    				tmp = ww * BarrettMod(pF[1] * pG[1], one_pre, p);
    				r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(tmp, one_pre, p);
    				//r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(ww * BarrettMod(pF[1] * pG[1], one_pre, p), one_pre, p);
    				r1 = BarrettMod(pF[0] * pG[1], one_pre, p) + BarrettMod(pF[1] * pG[0], one_pre, p);
    				pFG[0] = r0 - (((P - r0) >> 63)&P);
    				pFG[1] = r1 - (((P - r1) >> 63)&P);
    
    				pFG += 2;
    				pF += 2;
    				pG += 2;
    
    				//使用冗余表示,[0,2p)
    				ww = pW[list[i + 1]]; //第k层第2^{k-1}+i个多项式使用的单位根,w_{2^k}^{brv_k(2^{k-1}+i)}
    
    				tmp = ww * BarrettMod(pF[1] * pG[1], one_pre, p);
    				r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(tmp, one_pre, p);
    				//r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(ww * BarrettMod(pF[1] * pG[1], one_pre, p), one_pre, p);
    				r1 = BarrettMod(pF[0] * pG[1], one_pre, p) + BarrettMod(pF[1] * pG[0], one_pre, p);
    				pFG[0] = r0 - (((P - r0) >> 63)&P);
    				pFG[1] = r1 - (((P - r1) >> 63)&P);
    
    				pFG += 2;
    				pF += 2;
    				pG += 2;
    
    				//使用冗余表示,[0,2p)
    				ww = pW[list[i + 2]]; //第k层第2^{k-1}+i个多项式使用的单位根,w_{2^k}^{brv_k(2^{k-1}+i)}
    
    				tmp = ww * BarrettMod(pF[1] * pG[1], one_pre, p);
    				r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(tmp, one_pre, p);
    				//r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(ww * BarrettMod(pF[1] * pG[1], one_pre, p), one_pre, p);
    				r1 = BarrettMod(pF[0] * pG[1], one_pre, p) + BarrettMod(pF[1] * pG[0], one_pre, p);
    				pFG[0] = r0 - (((P - r0) >> 63)&P);
    				pFG[1] = r1 - (((P - r1) >> 63)&P);
    
    				pFG += 2;
    				pF += 2;
    				pG += 2;
    
    				//使用冗余表示,[0,2p)
    				ww = pW[list[i + 3]]; //第k层第2^{k-1}+i个多项式使用的单位根,w_{2^k}^{brv_k(2^{k-1}+i)}
    
    				tmp = ww * BarrettMod(pF[1] * pG[1], one_pre, p);
    				r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(tmp, one_pre, p);
    				//r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(ww * BarrettMod(pF[1] * pG[1], one_pre, p), one_pre, p);
    				r1 = BarrettMod(pF[0] * pG[1], one_pre, p) + BarrettMod(pF[1] * pG[0], one_pre, p);
    				pFG[0] = r0 - (((P - r0) >> 63)&P);
    				pFG[1] = r1 - (((P - r1) >> 63)&P);
    
    				pFG += 2;
    				pF += 2;
    				pG += 2;
    			}
    		else
    			for (int32 i = 0; i < hm; i++) //hm个h长小多项式
    			{
    				//使用冗余表示,[0,2p)
    				ww = pW[list[i]]; //第k层第2^{k-1}+i个多项式使用的单位根,w_{2^k}^{brv_k(2^{k-1}+i)}
    
    				int64 tmp = ww * BarrettMod(pF[1] * pG[1], one_pre, p);
    				r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(tmp, one_pre, p);
    				//r0 = BarrettMod(pF[0] * pG[0], one_pre, p) + BarrettMod(ww * BarrettMod(pF[1] * pG[1], one_pre, p), one_pre, p);
    				r1 = BarrettMod(pF[0] * pG[1], one_pre, p) + BarrettMod(pF[1] * pG[0], one_pre, p);
    				pFG[0] = r0 - (((P - r0) >> 63)&P);
    				pFG[1] = r1 - (((P - r1) >> 63)&P);
    
    				pFG += 2;
    				pF += 2;
    				pG += 2;
    			}
    	else
    		for (int32 i = 0; i < hm; i++) //hm个h长小多项式
    		{
    			ww = pW[list[i]]; //第k层第2^{k-1}+i个多项式使用的单位根,w_{2^k}^{brv_k(2^{k-1}+i)}
    			convolution(pFG, pF, pG, ww);
    
    			pFG += h;
    			pF += h;
    			pG += h;
    		}
    }
    
    PolyBarrett
    #include "PolyBarrett.h"
    
    
    // h = f / g mod p
    void Div(Poly& h, Poly f, Poly g, int64 p)
    {
    	int df = deg(f);
    	int dg = deg(g);
    
    	int64 a = g[dg];
    	int64 ainv = inv(a, p);
    
    	int64* pf = f.data();
    	int64* pg = g.data();
    
    	for (int i = df; i >= dg; i--)
    	{
    		int64 b = pf[i];
    		int64 c = (b * ainv) % p;
    		h[i - dg] = c;
    		pf[i] = 0;
    
    		if (c == 0)
    			continue;
    
    		int64* ptr = pf + (i - dg);
    		for (int j = 0; j < dg; j++)
    		{
    			ptr[j] -= c * pg[j];
    			ptr[j] %= p;
    			ptr[j] = ptr[j] < 0 ? ptr[j] + p : ptr[j];
    		}
    	}
    }
    
    
    void PolyBarrett::init(uint64 p, Poly& g, int32 acc, int32 h)
    {
    	this->p = p;
    	this->dg = deg(g);
    	this->acc = acc;
    	this->k = dg + acc;
    	this->g = g;
    
    	this->N = 2 * dg + acc;	//因为需要计算f*m
    	this->H = h;
    	this->K = ceil(log2(N / h));
    	this->n = pow(2, K)*H;
    
    	uint64 x = (this->n * p*p) >> K;	// P > n*p^2
    	do {
    		x += 1;
    		P = (x << K) + 1;
    	} while (!Miller_Rabin(P, 10));	// 2^k | P-1
    
    	int HN = pow(2, K) / 2;
    	for (W = 2; W < P; W++)
    		if (pow_mod(W, HN, P) == P - 1)	// order(w)=2^K
    			break;
    
    	Poly xk(k + 1);
    	xk[k] = 1;
    	this->m.resize(this->acc + 1);
    	Div(this->m, xk, g, p);
    
    	printf("Barrett of Polynomials: acc=%d, dg=%d, p=%lld -> ", acc, dg, p);
    	printf("P=%lld, W=%lld, K=%d, H=%d.\n", P, W, K, H);
    
    	this->NTT.init(P, W, K, H);
    	this->NTT.FwdNTT(this->M, this->m);
    	this->NTT.FwdNTT(this->G, this->g);
    }
    
    
    void PolyBarrett::Mod(int64* h, int64* f, int flen)
    {
    	Poly t;
    	NTTRep T, F(this->n);
    
    	NTT.FwdNTT(F.data(), f, flen);
    	NTT.NTTMul(T, F, M);
    	NTT.InvNTT(t, T);
    
    	int64* pt = t.data();
    	for (int i = 0; i < this->n - this->k; i++)
    	{
    		pt[i] = pt[i + k] % p;
    	}
    	for (int i = this->n - this->k; i < this->n; i++)
    	{
    		pt[i] = 0;
    	}
    
    	NTT.FwdNTT(T, t);
    	NTT.NTTMul(T, T, G);
    	NTT.InvNTT(t, T);
    
    	for (int i = 0; i < this->dg; i++)
    	{
    		h[i] = (f[i] - pt[i]) % p;
    		h[i] += h[i] < 0 ? p : 0;
    	}
    }
    
    GFX
    #include "GFX.h"
    
    
    //g = f mod (x^n - 1)
    void ModXnSubOne(GFX& g, GFX& f, int n)
    {
    	int fn = f.n;
    	int fd = f.d;
    	int fp = f.p;
    
    	if (&g != &f)
    		for (int i = 0; i < fn; i++)
    		{
    			int64* p1 = g.rep[i%n].data();
    			int64* p2 = f.rep[i].data();
    			for (int j = 0; j < fd; j++)
    			{
    				p1[j] += p2[j];
    				p1[j] %= fp;
    			}
    		}
    	else
    		for (int i = n; i < fn; i++)
    		{
    			int64* p1 = g.rep[i%n].data();
    			int64* p2 = f.rep[i].data();
    			for (int j = 0; j < fd; j++)
    			{
    				p1[j] += p2[j];
    				p1[j] %= fp;
    				p2[j] = 0;
    			}
    		}
    }
    
    //g = f mod m(x),p
    void Modm(int64* g, int64* f, int flen, Poly& m, uint64 p)
    {
    	int df;
    	for (df = flen - 1; df >= 0; df--)
    		if (f[df] != 0)
    			break;
    	int dm = deg(m);
    
    	if (df < dm)
    	{
    		for (int i = 0; i <= df; i++)
    			g[i] = f[i] % p;
    		return;
    	}
    
    	int64* pf = new int64[df + 1];
    	for (int i = 0; i < df + 1; i++)
    		pf[i] = f[i];
    
    	int64 a = m[dm];
    	int64 ainv = inv(a, p);
    	int64* pm = m.data();
    
    	for (int i = df; i >= dm; i--)
    	{
    		int64 b = pf[i];
    		int64 c = (b * ainv) % p;
    		g[i - dm] = c;
    		pf[i] = 0;
    
    		if (c == 0)
    			continue;
    
    		int64* ptr = pf + (i - dm);
    		for (int j = 0; j < dm; j++)
    		{
    			ptr[j] -= c * pm[j];
    			ptr[j] %= p;
    			ptr[j] = ptr[j] < 0 ? ptr[j] + p : ptr[j];
    		}
    	}
    
    	for (int i = 0; i < dm; i++)
    		g[i] = pf[i];
    	delete pf;
    }
    
    ostream& operator<<(ostream& cout, GFX& f)
    {
    	int deg1 = f.deg();
    	printf("[");
    	for (int i = 0; i <= deg1; i++)
    	{
    		int deg2 = deg(f.rep[i]);
    		int64* p = f.rep[i].data();
    
    		printf(" [");
    		for (int j = 0; j <= deg2; j++)
    		{
    			printf(" %lld ", p[j]);
    		}
    		printf("] ");
    	}
    	printf("]");
    	return cout;
    }
    
    //c = a+b
    void GFX_Calculator::Add(GFX&c, GFX& a, GFX&b)
    {
    	int d = a.d;
    	int n = a.n;
    	uint64 p = a.p;
    
    	for (int i = 0; i < n; i++)
    	{
    		int64* pa = a.rep[i].data();
    		int64* pb = b.rep[i].data();
    		int64* pc = c.rep[i].data();
    
    		for (int j = 0; j < d; j++)
    		{
    			pc[j] = pa[j] + pb[j];
    			pc[j] -= ((p - pc[j]) >> 63)&p; //if c>p then c = c-p
    		}
    	}
    }
    
    //c = a-b
    void GFX_Calculator::Sub(GFX&c, GFX& a, GFX&b)
    {
    	int d = a.d;
    	int n = a.n;
    	uint64 p = a.p;
    
    	for (int i = 0; i < n; i++)
    	{
    		int64* pa = a.rep[i].data();
    		int64* pb = b.rep[i].data();
    		int64* pc = c.rep[i].data();
    
    		for (int j = 0; j < d; j++)
    		{
    			pc[j] = pa[j] - pb[j];
    			pc[j] += (pc[j] >> 63)&p; //if c<0 then c = c+p
    		}
    	}
    }
    
    
    //寻找合适的Z_P,以及2^k次单位根w
    void GFX_Calculator::findP(int d, int n, uint64 p, uint64& P, uint64& w, int& k, int h)
    {
    	int N = 2 * ((2 * d)*n);
    	k = ceil(log2(N / h));
    	int N2 = pow(2, k);	// N2 * h = h*2^k >= N
    
    	uint64 x = (N2*h * p*p) >> k;	// P > n*p^2
    	do {
    		x += 1;
    		P = (x << k) + 1;
    	} while (!Miller_Rabin(P, 10));	// 2^k | P-1
    
    	int HN2 = N2 >> 1;
    	for (w = 2; w < P; w++)
    		if (pow_mod(w, HN2, P) == P - 1)	// order(w)=N2
    			break;
    }
    
    
    //c = a*b,GF(p^d)=Z_p[x]/(m(x))
    void GFX_Calculator::Mul(GFX&c, GFX& a, GFX&b)
    {
    	int d = a.d;
    	int n = a.n;
    	uint64 p = a.p;
    	int BlockSize = 2 * d;
    	int BlockNum = 2 * n;
    	int Len = BlockSize * BlockNum;
    
    	Poly A(Len);
    	Poly B(Len);
    	Poly C(Len);
    
    	for (int i = 0; i < n; i++)
    	{
    		int64* pa = a.rep[i].data();
    		int64* pb = b.rep[i].data();
    		int64* pA = A.data() + (2 * d)*i;
    		int64* pB = B.data() + (2 * d)*i;
    
    		for (int j = 0; j < d; j++)
    		{
    			pA[j] = pa[j];
    			pB[j] = pb[j];
    		}
    	}
    
    	NTTRep A2, B2, C2;
    
    	intt.FwdNTT(A2, A);
    	intt.FwdNTT(B2, B);
    	intt.NTTMul(C2, A2, B2);
    	intt.InvNTT(C, C2);	//C=A*B
    
    	int64* pC = C.data();
    	for (int i = 0; i < Len; i++)
    		pC[i] %= p;
    
    	/*cout << "A = " << A; pn;
    	cout << "B = " << B; pn;
    	cout << "C = " << C; pn;*/
    
    	c.zero();
    	for (int i = 0; i < BlockNum; i++)
    	{
    		int64* pC = C.data() + BlockSize * i;
    		PB.Mod(c.rep[i].data(), pC, BlockSize);
    	}
    }
    
    
    //c = a*b,GF(p^d)=Z_p[x]/(m(x))
    void GFX_Calculator::Mul2(GFX&c, GFX& a, GFX&b)
    {
    	int d = a.d;
    	int n = a.n;
    	uint64 p = a.p;
    	int BlockSize = 2 * d;
    	int BlockNum = 2 * n;
    	int Len = BlockSize * BlockNum;
    
    	Poly A(Len);
    	Poly B(Len);
    	Poly C(Len);
    
    	for (int i = 0; i < n; i++)
    	{
    		int64* pa = a.rep[i].data();
    		int64* pb = b.rep[i].data();
    		int64* pA = A.data() + (2 * d)*i;
    		int64* pB = B.data() + (2 * d)*i;
    
    		for (int j = 0; j < d; j++)
    		{
    			pA[j] = pa[j];
    			pB[j] = pb[j];
    		}
    	}
    
    	NTTRep A2, B2, C2;
    
    	intt.FwdNTT(A2, A);
    	intt.FwdNTT(B2, B);
    	intt.NTTMul(C2, A2, B2);
    	intt.InvNTT(C, C2);	//C=A*B
    
    	int64* pC = C.data();
    	for (int i = 0; i < Len; i++)
    		pC[i] %= p;
    
    	/*cout << "A = " << A; pn;
    	cout << "B = " << B; pn;
    	cout << "C = " << C; pn;*/
    
    	c.zero();
    	for (int i = 0; i < BlockNum; i++)
    	{
    		int64* pC = C.data() + BlockSize * i;
    		Modm(c.rep[i].data(), pC, BlockSize, m, p);
    	}
    }
    

    结果

    GFX Mul
    m = [ 1, 0, 0, 0, 0, 0, 0, 1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1 ]
    
    GFX Mul: d=149,n=53,p=2 -> P=147457, W=22, K=14, H=2.
    Barrett of Polynomials: acc=149, dg=149, p=2 -> P=3329, W=17, K=8, H=2.
    a = [ [ 1 ]  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  [ 1 ] ]
    
    b = [ [ 1  0  1 ]  [ 1 ] ]
    
    10 rounds
    
    使用Barrett:a*b = [ [ 1  0  1 ]  [ 1 ]  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  [ 1  0  1 ]  [ 1 ] ]
    0.155162 s
    
    使用长除法:a*b = [ [ 1  0  1 ]  [ 1 ]  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  [ 1  0  1 ]  [ 1 ] ]
    0.291891 s
    
    用长除法计算:a*b mod (x^n - 1) = [ [ 0  0  1 ]  [ 1 ]  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  []  [ 1  0  1 ] ]
    0.0015902 s
    
    展开全文
  • 二元多项式除法电路原理 例:g(x)=x4+x2+x+1g(x)=x^4 + x^2+x+1g(x)=x4+x2+x+1,a(x)=x5+x3+1a(x) = x^5+x^3+1a(x)=x5+x3+1,寄存器个数等于g(x)g(x)g(x)次数为4,a(x)a(x)a(x)除以g(x)g(x)g(x)的运算过程如下: g...

    关注公号【逆向通信猿】试读更多内容!!!

    二元多项式除法电路原理

    在这里插入图片描述
    例: g ( x ) = x 4 + x 2 + x + 1 g(x)=x^4 + x^2+x+1

    展开全文
  • AES:有限域多项式乘法逆元求解

    千次阅读 2017-11-07 17:23:00
    cout多项式(";polynomialtostring(b);cout)mod(";polynomialtostring(m);cout)的乘法逆元是(";polynomialtostring(x);cout)"; system("pause"); return 0; } 运行结果如下图

    题目:


    求解算法,扩展的欧几里得算法

    /*
    @author tilltheendwjx
    @blog http://blog.csdn.net/wjh200821或者http://www.cnblogs.com/tilltheendwjx/ 
    */ 
    #include<iostream>
    using namespace std;
    int indexofmax1(int value)
    {
        int tmp=1;
        int count=0;
        for(int i=0;i<sizeof(int)*8;++i)
        {
              if((value&tmp))
                 count=i;
              tmp=tmp*2;
        }
        return count;
    }
    void polynomialtostring(int value)
    {
        int tmp=1;
        int flag=0;
        int c=indexofmax1(value);
        for(int i=0;i<sizeof(int)*8;++i)
        {
              if((value&tmp))
              {
                 if(i==0)
                 {
                   cout<<"1";
                 }else if(i==1)
                 {
                   cout<<"x";
                 }else
                 {
                   cout<<"x^"<<i;
                 }
                 flag=1;
                 if(i<c)
                   cout<<"+";
              }   
              tmp=tmp*2;
        }
        if(flag==0)
          cout<<"0";
    }
    int powofvalue(int value)
    {
        return 1<<(value);
    }
    int divide(int m,int b,int &remainvalue)
    {
        int mindex=indexofmax1(m);
        int vindex=indexofmax1(b);
        if(mindex<vindex)
        {
            remainvalue=m;
            return 0;
        }
        int c=mindex-vindex;
        int tmp=b;
        tmp=tmp<<c;
        m=m^tmp;
        return powofvalue(c)|divide(m,b,remainvalue);
    }
    int Tx(int ax,int q,int bx)
    {
        //cout<<endl;
        //cout<<ax<<"\t"<<bx<<"\t";
        int tmp=1;
        int value=0;
        for(int i=0;i<sizeof(int)*8;++i)
        {
              if((q&tmp))
              {
                 value=value^((bx<<i));   
              }   
              tmp=tmp*2;
        }
        //cout<<ax<<"\t"<<value<<"\t";
        //cout<<endl;
        return ax^(value);
    }
    int extent_gcd(int m,int b,int &x,int &y)
    {
       int a1=1,a2=0,a3=m;
       int b1=0,b2=1,b3=b;
       int remainvalue=0;
       while(1)
       {
               polynomialtostring(a1);
               cout<<"    ";
               polynomialtostring(a2);
               cout<<"    ";
               polynomialtostring(a3);
               cout<<"    ";
               polynomialtostring(b1);
               cout<<"    ";
               polynomialtostring(b2);
               cout<<"    ";
               polynomialtostring(b3);
               cout<<"    ";
              if(b3==0)
                  return a3;
              if(b3==1)
                   return b3;
              int q=divide(a3,b3,remainvalue);
              int t1=Tx(a1,q,b1);
              int t2=Tx(a2,q,b2);
              int t3=remainvalue;
              cout<<t1<<endl;
              cout<<t2<<endl;
              a1=b1;a2=b2;a3=b3;
              b1=t1;b2=t2;b3=t3;
              x=b2;y=b3;
              polynomialtostring(q);
              cout<<endl;
       } 
    }
    int main(void)
    {
    int m=283,b=83,x=0,y=0;
    cout<<"中间结果如下:"<<endl; 
    cout<<"a1   a2    a3    b1    b2    b3    q"<<endl;
    int reault=extent_gcd(m,b,x,y);
    cout<<endl;
    cout<<"多项式(";polynomialtostring(b);cout<<")mod(";polynomialtostring(m);cout<<")的乘法逆元是(";polynomialtostring(x);cout<<")"<<endl;
    system("pause"); 
    return 0;
    }
    

    运行结果如下图


    展开全文
  • 有限域的本原多项式

    千次阅读 2020-01-02 09:07:27
    查询有限域上的本原多项式: http://wims.unice.fr/wims/cn_tool~algebra~primpoly.cn.html
  • 有限域的构造之常见本原多项式

    万次阅读 2018-11-28 14:26:14
  • 不可约多项式
  • 共回答了14个问题采纳率:...“有限域GF(2)上的多项式”,说明:(1)多项式的系数只能是0或1,(不能是2,3,.也不能是-1,-2.)(2)同类项合并时的运算按照上面的GF(2)上的加法运算例子:F(X)*G(X):(1)先做普通多项是乘法:...
  • 利用更为简明的方法证明了有关有限域上的以 为周期的次不可约多项式的个数的结论,另外,本文结合初等数论知识得到了前面这个结论的几 个推论,并对利用低次不可约多项式构造高次不可约多项式进行了研究
  • #在有限域GF(2^n)下求多项式乘法 Python代码实现 一、了解运算规则 二、例题展示 三、直接上代码(代码有详细备注,不做一一解释) def yxydxscf(a,b,c): # 不可约多项式系数模二运算 c = str(c)[1:] # 不够8...
  • 最近在复习现代密码理论中的AES,AES中的字节变换的核心操作就是求GF(28)GF(2^8)GF(28)上的多项式逆元,这个问题困扰了我一段时间,今天终于得到解决,其实计算方式和数论中求两个数的Bezout算法是一样的,这里感谢...
  • 多项式的逆元

    千次阅读 2020-05-12 19:04:54
    F=Zp[x]f(x)\mathbb{F}=\mathbb{Z}_p[x]_{f(x)}F=Zp​[x]f(x)​,其中f(x)f(x)f(x)为Zp\mathbb{Z}_pZp​上的不可约多项式多项式g(x)∈Fg(x)\in\mathbb{F}g(x)∈F,求g(x)g(x)g(x)在F\mathbb{F}F上的逆元。...
  • 文章目录一、前言二、数学基础1、GF(2⁸)有限域内的多项式2、不可约多项式3、多项式模运算3、乘法逆元三、算法步骤1、扩展欧几里得算法2、多项式除法3、多项式乘法四、代码实现1、多项式除法2、多项式乘法3、...
  • 有限域上环多项式的显式分解
  • GF⁡(2r)\operatorname {GF}(2^r)GF(2r)中元素的α\alphaα方幂与向量表示 例如:r=4r=4r=4,GF⁡(24)=F2[x]/(x4+x+1)\operatorname {GF}(2^4)=F_2[x]/(x^4+x+1)GF(24...根据本原多项式x4+x+1x^4+x+1x4+x+1,得α4+α+1
  • 证明了对任意给定正整数n和k,如果满足k
  • 对计算有限域上切比雪夫多项式的特征多项式算法(CPA)进行了改进,以提高算法的执行速度。首先用蒙哥马利模乘代替原有算法中的普通模乘运算,从而降低单次模乘运算的平均运行时间;其次对蒙哥马利模平方运算的算法流程...
  • 对于给定的置换多项式f(x)∈Fq[x],研究f(x)是否Fqr(r>1)上的置换多项式,是研究有限域上置换多项式的主要问题之一.本文改进了Carlitz和万大庆的方法,完全解决了形如x q-1/4+1+ax的多项式是否Fqr上置换多项式的...
  • 关于有限域多项式因式分解,和不错的,对有限域多项式概括很完整
  • 本文证明了定理设F是一个特征为P的含P~a个元的有限域.f(x)=f_1(x)~l1…f_k(x)~lk是f(x)在多项式环F[x]中的标准分解式,f_i(x)是最高系数为1、次数为n_i的不可约多项式.那么f(x)有原根的充分必要条件为当p≥3时:k=1...
  • 定义:线性码C\mathscr CC,是有限域上线性空间GF(q)nGF(q)^nGF(q)n的子空间。长度为nnn,子空间维度为kkk,最小距离为ddd,那么记做(n,k,d)(n,k,d)(n,k,d)码。 性质最好的编码,几乎都是线性码,因此编码理论只研究...
  • 有限域多项式的多发送者认证码的新构造
  • 大数据-算法-有限域多项式形式的约化RSA公钥密码算法.pdf
  • 当然也可以参考《密码编码学与网络安全》这书的有限域一章。形象地说,域有这样一个性质:在加法和乘法上具有封闭性。也就是说对域中的元素进行加法或乘法运算后的结果仍然是域中的元素。有一点要注意,域里面的乘法...
  • 关于素域Fp上多项式x^p-ax-b的阶的注记,武跟强,李研超,本文研究了有限域Fp上多项式x^p-ax-b的阶。首先给出并证明了x^p-x-b与x^p-x-b’两者的阶之间的一个关系;其次探讨了x^p-x-b的阶必须满足的�
  • 本文利用具有线性结构的多项式和线性化多项式得到了一种形式为L1(x)+ L-1(γ)h(f(x))的置换多项式,该结果推广了Kyureghyan在2011年得到的一个结果.本文还利用具有线性结构的多项式和核的维数为k+1线性化多项式构造...
  • 1) * ['0'])]))) # target: 两多项式相加 # params: 多项式0,多项式1 # return: 多项式之和 def add(self, poly0, poly1): poly0_lt = poly0.poly poly1_lt = poly1.poly align_lt = [0] * abs(poly0.length-poly1....
  • 该方案基于有限域上Chebyshev多项式良好的半群特性, 运用RSA算法巧妙地隐藏通信双方产生的有限域上的Chebyshev多项式值, 从而避免了以往的种种主动攻击, 保证了密钥协商的安全; 同时, 该密钥协商方案还实现了身份...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 6,470
精华内容 2,588
关键字:

有限域的多项式