kk Blog —— 通用基础


date [-d @int|str] [+%s|"+%F %T"]
netstat -ltunp
sar -n DEV 1

生成树计数

http://blog.sina.com.cn/s/blog_a55522150102w6sp.html http://www.xuebuyuan.com/1622347.html

基尔霍夫定理

算法思想:

1
2
3
4
5
6
7
8
9
10
*(1)G的度数矩阵D[G]是一个n*n的矩阵,并且满足:当i≠j时,dij=0;当i=j时,dij等于vi的度数; 
*(2)G的邻接矩阵A[G]是一个n*n的矩阵,并且满足:如果vi,vj之间有边直接相连,则aij=1,否则为0; 
*定义图G的Kirchhoff矩阵C[G]为C[G]=D[G]-A[G]; 
*Matrix-Tree定理:G的所有不同的生成树的个数等于其Kirchhoff矩阵C[G]任何一个n-1阶主子式的行列式的绝对值; 
*所谓n-1阶主子式,就是对于r(1≤r≤n),将C[G]的第r行,第r列同时去掉后得到的新矩阵,用Cr[G]表示; 
* 
*Kirchhoff矩阵的特殊性质: 
*(1)对于任何一个图G,它的Kirchhoff矩阵C的行列式总是0,这是因为C每行每列所有元素的和均为0; 
*(2)如果G是不连通的,则它的Kirchhoff矩阵C的任一个主子式的行列式均为0; 
*(3)如果G是一颗树,那么它的Kirchhoff矩阵C的任一个n-1阶主子式的行列式均为1; 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <iostream>
#include <algorithm>
#define  zero(x) (((x)>0?(x):(-x))<1e-15)
using namespace std;
const int MAXN = 110;
double a[MAXN][MAXN], b[MAXN][MAXN];
int G[MAXN][MAXN];
int N, M;

/*
*生成树计数
*1、G的度数矩阵D[G]是一个n*n的矩阵,并且满足:当i≠j时,dij=0;当i=j时,dij等于vi的度数。
*2、G的邻接矩阵A[G]也是一个n*n的矩阵, 并且满足:如果vi、vj之间有边直接相连,则aij=1,否则为0。
*我们定义G的Kirchhoff矩阵(也称为拉普拉斯算子)C[G]为C[G]=D[G]-A[G],则Matrix-Tree定理可以描述为:
*G的所有不同的生成树的个数等于其Kirchhoff矩阵C[G]任何一个n-1阶主子式的行列式的绝对值。
*所谓n-1阶主子式,就是对于r(1≤r≤n),将C[G]的第r行、第r列同时去掉后得到的新矩阵,用Cr[G]表示。
*/

double Det(double a[MAXN][MAXN], int n)
{
	int i, j, k, sign = 0;
	double ret = 1, t;
	for(i = 0; i < n; ++i)
		for(j = 0; j < n; ++j)
			b[i][j] = a[i][j];
	for(i = 0; i < n; ++i)
	{
		if(zero(b[i][i]))
		{
			for(j = i + 1; j < n; ++j)
			{
				if(!zero(b[j][i]))
					break;
			}
			if(j == n)
				return 0;
			for(k = i; k < n; ++k)
				t = b[i][k], b[i][k] = b[j][k], b[j][k] = t;
			sign++;
		}
		ret *= b[i][i];
		for(k = i + 1; k < n; ++k)
			b[i][k] /= b[i][i];
		for(j = i + 1; j < n; ++j)
			for(k = i + 1; k < n; ++k)
				b[j][k] -= b[j][i] * b[i][k];
	}
	if(sign & 1)
		ret = - ret;
	return ret;
}

int main()
{
	int T, u, v;
	scanf("%d", &T);
	while (T--)
	{
		scanf("%d %d", &N, &M);
		memset(G, 0, sizeof(G));
		memset(a, 0, sizeof(a));
		while(M--)
		{
			scanf("%d %d", &u, &v);
			G[u-1][v-1] = G[v-1][u-1] = 1;
		}
		for(int i = 0; i < N; ++i)
		{
		   int d = 0;
		   for (int j = 0; j < N; ++j) if(G[i][j]) d++;
		   a[i][i] = d;
		}
		for(int i = 0; i < N; ++i)
		{
			for (int j = 0; j < N; ++j)
			{
				if (G[i][j]) a[i][j] = -1;
			}
		}
	   double ans = Det(a, N - 1);
	   printf("%0.0lf\n", ans);
	}
	return 0;
}

SRM733-DIV1-B

https://community.topcoder.com/stat?c=round_stats&rd=17140&dn=1

Problem Statement

     Consider an undirected, complete, simple graph G on n vertices. The vertices of G are labeled from 1 to n. Specifically, each pair of distinct vertices is connected by a single undirected edge, so there are precisely n*(n-1)/2 edges in this graph.

You are given a set E that contains some edges of the graph G. More precisely, you are given the vector s x and y. For each valid i, (x[i], y[i]) is one of the edges in E.

A spanning tree of G is a subset of exactly n-1 edges of G such that the edges connect all n vertices. You may note that the edges of a spanning tree do indeed form a tree that is a subgraph of G and spans all its vertices.

We are interested in spanning trees that contain all of the edges in the provided set E. Compute and return the number of such spanning trees modulo 987,654,323. Two spanning trees are different if there is an edge of G that is in one of them but not in the other. Definition     

Class:

BuildingSpanningTreesDiv1

Method:

getNumberOfSpanningTrees

Parameters:

int, vector , vector

Returns:

int

Method signature:

int getNumberOfSpanningTrees(int n, vector x, vector y)

(be sure your method is public)

Limits

Time limit (s): 2.000

Memory limit (MB): 256

Notes

987,654,323 is a prime number.

Constraints

n will be between 1 and 1,000, inclusive.
x will contain between 1 and 1,000 elements, inclusive.
y will contain between 1 and 1,000 elements, inclusive.
x and y will contain the same number of elements.
Each element of x will be between 1 and n-1, inclusive.
Each element of y will be between 2 and n, inclusive.
For each valid i, x[i] will be less than y[i].
All edges in E will be distinct.

Examples

0)
3
{1,2}
{2,3}
Returns: 1
The edges in the provided set E alredy form a spanning tree, so there is exactly one spanning tree that contains them.

1)
5
{1,3,4}
{2,4,5}
Returns: 6
There are six ways to add the one missing edge: one endpoint of the new edge must lie in the set {1,2} and the other in the set {3,4,5}.

2)
4
{1}
{2}
Returns: 8
There are eight spanning trees that contain the edge (1, 2):
{(1, 2), (1, 3), (1, 4)}
{(1, 2), (1, 3), (2, 4)}
{(1, 2), (1, 3), (3, 4)}
{(1, 2), (2, 3), (2, 4)}
{(1, 2), (1, 4), (2, 3)}
{(1, 2), (1, 4), (3, 4)}
{(1, 2), (2, 3), (3, 4)}
{(1, 2), (2, 4), (3, 4)}

3)
4
{1,2,1}
{2,3,3}
Returns: 0
The set E contains a cycle, and thus there is no spanning tree that contains all these edges.

4)
8
{1,4,2,3,5}
{4,7,6,5,8}
Returns: 144

5)
1000
{1}
{2}
Returns: 529013784

Don’t forget to take the modulo.

This problem statement is the exclusive and proprietary property of TopCoder, Inc. Any unauthorized use or reproduction of this information without the prior written consent of TopCoder, Inc. is strictly prohibited. ©2003, TopCoder, Inc. All rights reserved.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
#include <algorithm>
#include <cassert>
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <iostream>
#include <map>
#include <queue>
#include <set>
#include <sstream>
#include <stdint.h>
#include <string>
#include <utility>
#include <vector>
 
using namespace std;
 
int const MAX_N = 1003;
int const MOD = 987654323;
 
int p [MAX_N];
int s [MAX_N];
 
int root (int v)
{
	if (p[v] != v)
	{
		p[v] = root (p[v]);
	}
	return p[v];
}
 
bool unite (int u, int v)
{
	u = root (u);
	v = root (v);
	if (u == v)
	{
		return false;
	}
	p[v] = u;
	s[u] += s[v];
	return true;
}
 
int powmod (int a, int b)
{
	int res = 1;
	for ( ; b > 0; b >>= 1)
	{
		if (b & 1)
		{
			res = (res * 1LL * a) % MOD;
		}
		a = (a * 1LL * a) % MOD;
	}
	return res;
}
 
int a [MAX_N] [MAX_N];
int n;
 
int det (void)
{
	int res = 1;
	for (int i = 1; i < n; i++)
	{
		int j;
		for (j = i; j < n; j++)
		{
			if (a[j][i])
			{
				break;
			}
		}
		if (j == n)
		{
			return 0;
		}
		if (j != i)
		{
			for (int k = i; k < n; k++)
			{
				swap (a[i][k], a[j][k]);
			}
		}
		res = (res * 1LL * a[i][i]) % MOD;
		int inv = powmod (a[i][i], MOD - 2);
		for (int j = i + 1; j < n; j++)
		{
			int mult = (a[j][i] * 1LL * inv) % MOD;
			for (int k = i; k < n; k++)
			{
				a[j][k] = (a[j][k] -
						a[i][k] * 1LL * mult) % MOD;
			}
		}
	}
	return (res + MOD) % MOD;
}
 
class BuildingSpanningTreesDiv1
{
public:
	int getNumberOfSpanningTrees (int n, vector <int> x, vector <int> y)
	{
		for (int i = 0; i < n; i++)
		{
			p[i] = i;
			s[i] = 1;
		}
		int k = (int) (x.size ());
		for (int w = 0; w < k; w++)
		{
			int u = x[w] - 1;
			int v = y[w] - 1;
			if (!unite (u, v))
			{
				return 0;
			}
			a[u][v] += 1;
			a[v][u] += 1;
			a[u][u] -= 1;
			a[v][v] -= 1;
		}
 
		int t = 0;
		for (int i = 0; i < n; i++)
		{
			if (p[i] == i)
			{
				s[t] = s[i];
				t += 1;
			}
		}
 
		n -= k;
		::n = n;
 
		memset (a, 0, sizeof (a));
		for (int i = 0; i < n; i++)
		{
			for (int j = i + 1; j < n; j++)
			{
				int e = s[i] * s[j];
				a[i][j] -= e;
				a[j][i] -= e;
				a[i][i] = (a[i][i] + e) % MOD;
				a[j][j] = (a[j][j] + e) % MOD;
			}
		}
 
		int res = det ();
		return res;
	}
};

逆元

https://blog.csdn.net/baidu_35643793/article/details/75268911

1.什么是逆元

当求解公式:(a/b)%m 时,因b可能会过大,会出现爆精度的情况,所以需变除法为乘法:

设c是b的逆元,则有b*c≡1(mod m);

则(a/b)%m = (a/b)1%m = (a/b)bc%m = ac(mod m);

即a/b的模等于a*b的逆元的模;

逆元就是这样应用的;

2.求逆元的方法

(1).费马小定理

在是素数的情况下,对任意整数都有。 如果无法被整除,则有。 可以在为素数的情况下求出一个数的逆元,,即为逆元。

题目中的数据范围1<=x<=109,p=1000000007,p是素数;

所以x肯定就无法被p整除啊,所以最后就得出xp-2为x的逆元啦。

复杂度O(logn);

代码:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
const int mod = 1000000009;  
long long quickpow(long long a, long long b) {  
	if (b < 0) return 0;  
	long long ret = 1;  
	a %= mod;  
	while(b) {  
		if (b & 1) ret = (ret * a) % mod;  
		b >>= 1;  
		a = (a * a) % mod;  
	}  
	return ret;  
}  
long long inv(long long a) {  
	return quickpow(a, mod - 2);  
}  

(2)扩展欧几里得算法求逆元

可扩展欧几里得求逆元ax≡1(mod n)其中a,n互质;

复杂度:O(logn);

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
ll extend_gcd(ll a, ll b, ll &x, ll &y) {  
	if (b == 0) {  
		x = 1, y = 0;  
		return a;  
	}  
	else {  
		ll r = extend_gcd(b, a % b, y, x);  
		y -= x * (a / b);  
		return r;  
	}  
}  
ll inv(ll a, ll n) {  
	ll x, y;  
	extend_gcd(a, n, x, y);  
	x = (x % n + n) % n;  
	return x;  
}  

(3) 逆元线性筛 ( P为质数 )

求1,2,…,N关于P的逆元(P为质数)

复杂度:O(N)

代码:

1
2
3
4
5
6
const int mod = 1000000009;  
const int maxn = 10005;  
int inv[maxn];  
inv[1] = 1;  
for(int i = 2; i < 10000; i++)  
	inv[i] = inv[mod % i] * (mod - mod / i) % mod;  

伽罗华域(Galois Field)上的四则运算

https://blog.csdn.net/shelldon/article/details/54729687

伽罗华域(Galois Field)上的四则运算

Évariste Galois ,伽罗华(也译作伽瓦罗),法国数学家,群论的创立者。用群论彻底解决了根式求解代数方程的问题,而且由此发展了一整套关于群和域的理论。 本文介绍伽罗华域,以及在伽罗华域上的四则运算方式。伽罗华域上的四则运算实际上是多项式计算,后文中详细介绍。

一、相关数学概念

1、域

一组元素的集合,以及在集合上的四则运算,构成一个域。其中加法和乘法必须满足交换、结合和分配的规律。加法和乘法具有封闭性,即加法和乘法结果仍然是域中的元素。
域中必须有加法单位元和乘法单位元,且每一个元素都有对应的加法逆元和乘法逆元。但不要求域中的 0有乘法逆元。

2、有限域

仅含有限多个元素的域。因为它由伽罗华所发现,因而又称为伽罗华域。

所以当我们说伽罗华域的时候,就是指有限域。 GF(2w)表示含有2w个元素的有限域。

3、单位元

Identity Element,也叫幺元(么元),通常使用e来表示单位元。单位元和其他元素结合时,并不会改变那些元素。 对于二元运算,若ae=a,e称为右单位元;若ea=a,e称为左单位元,若ae=e*a=a,则e称为单位元。

4、逆元

对于二元运算,若ab=e,则a称为b的左逆元素,b称为a的右逆元素。若ab=ba=e,则称a为b的逆元,b为a的逆元。

5、本原多项式

域中不可约多项式是不能够进行因子分解的多项式, 本原多项式 (primitive polynomial)是一种特殊的不可约多项式。当一个域上的本原多项式确定了,这个域上的运算也就确定了。本原多项式一般通过查表可得,同一个域往往有多个本原多项式。

通过将域中的元素化为多项式形式,可以将域上的乘法运算转化为普通的多项式乘法再模本原多项式。

二、多项式运算

由于GF(2w)上的四则运算是基于多项式运算的,这里先介绍多项式运算。 多项式一般长这个样子:f(x) = x6 + x^ 4 + x2 + x + 1。

1、多项式加减法

将两个多项式中相同阶数的项系数相加或相减。 例如 (x2 + x ) + (x + 1) = x2 + 2x +1

2、多项式乘法

将其中一个多项式的各项分别与另一个多项式的各项相乘,然后把相同指数的项的系数相加。 例如 (x2 + x) * (x + 1) = x2 * (x + 1) + x * (x + 1) = x3 + x2 + x2 + x

3、多项式除法

使用长除法。例如计算x3-12x2-42,除以x-3。使用长除法计算,商x2-9x-27,余数-123。

4、GF(2w)上的多项式运算

对于GF(2w)上的多项式计算,多项式系数只能取 0或1。(如果是GF(3w),那么系数可以取 0、 1、 2) GF(2w)的多项式加法中,合并阶数相同的同类项时,由于0+0=0,1+1=0,0+1=1+0=1,因此系数不是进行加法操作,而是进行异或操作。

GF(2w)的多项式减法等于加法,例如x ^4 – x4就等于x4 + x4

三、伽罗华域

1、有限域GF(p):

有限域GF(p),其中p为素数。GF(p)里面的加法和乘法与一般的加法和乘法差不多,区别是结果需要mod p,以保证结果都是域中的元素。GF(p)的加法和乘法单位元分别是 0和1。 GF(p)加法是(a+b) mod p,乘法是(a*b)mod p。

对于域中的乘法,当p为素数时,才能保证集合中的所有的元素都有乘法逆元(0除外)。即对于域中的任一个元素a,总能在域中找到另外一个元素b,使得a*b mod p 等于1。

说明:假如p等于10,其乘法单位元为1。对于元素2,找不到一个数a,使得2*a mod 10 等于1,即2没有乘法逆元。这时,在域上就不能进行除2运算。

2、有限域GF(2w):

GF(p)中p必须是一个素数,才能保证集合中的所有元素都有加法和乘法逆元(0除外)。但实际应用中,很多场合需要 0到255这256个数字能组成一个域。但256不是素数,小于256的最大素数为251,如果直接把大于等于251的数截断为250,则会丢失一部分数据。

因此引入了GF(pw),其中p为素数,通常取p=2。计算机领域中经常使用的是GF(28),8刚好是一个字节的比特数。为了保证单位元性质,GF(2w)上的加法运算和乘法运算,不再使用一般的加法和乘法,而是使用多项式运算。

四、本原多项式

伽罗华域的元素可以通过该域上的本原多项式生成。通过本原多项式得到的域,其加法单位元都是 0,乘法单位元是1。

以GF(23)为例,指数小于3的多项式共8个: 0, 1, x, x+1, x2, x2+1, x2 + x, x2+x+1。其系数刚好就是000,001, 010, 011, 100, 101, 110, 111,是0 到7这8个数的二进制形式。

GF(23)上有不只一个本原多项式,选一个本原多项式x3+x+1,这8个多项式进行四则运算后 mod (x3+x+1)的结果都是8个之中的某一个。因此这8个多项式构成一个有限域。

对于GF(23),取素多项式为x3 + x+1,那么多项式x2+x的乘法逆元就是x+1。系数对应的二进制分别为110和011。此时,我们就认为对应的十进制数6和3互为逆元。

部分 GF(2w)域经常使用的本原多项式如下:

通过本原多项式生成元素

设P(x)是GF(2w)上的某一个本原多项式,GF(2w)的元素产生步骤是:
1、给定一个初始集合,包含0,1和元素x,即 {0,1,x};
2、将这个集合中的最后一个元素,即x,乘以x,如果结果的度大于等于w,则将结果mod P(x)后加入集合;
3、直到集合有2w个元素,此时最后一个元素乘以x再mod P(x)的值等于1。

例如,GF(24)含有16个元素,本原多项式为P(x)=x4+x+1,除了 0、1外,另外14个符号均由本原多项式生成。 注意到x14=x3+1,此时计算x15=x14x=(x3+1)x=x4+x=1,生成结束。

生成元素 多项式表示 二进制表示 数值表示 推导过程
0 0 0000 0
x^0 x^0 0001 1
x^1 x^1 0010 2
x^2 x^2 0100 4
x^3 x^3 1000 8
x^4 x+1 0011 3 x^3*x = x^4 mod P(x) = x+1
x^5 x^2+x 0110 6 x^4*x = (x+1)*x = x^2+x
x^6 x^3+x^2 1100 12
x^7 x^3+x+1 1011 11 x^6*x = (x^3+x^2)*x = x^4 +x^3 mod P(x) = x^3+x+1
x^8 x^2+1 0101 5
x^9 x^3+x 1010 10
x^10 x^2+x+1 0111 7 x^9*x=(x^3+x)*x = x^4+x^2 mod P(x) = x^2+x+1
x^11 x^3+x^2+x 1110 14
x^12 x^3+x^2+x+1 1111 15 x^11*x=(x^3+x^2+x)*x = x^4+x^3+x^2 mod P(x) = x^3+x^2+x+1
x^13 x^3+x^2+1 1101 13 x^12*x=(x^3+x^2+1 )*x = x^4+x^3+x mod P(x)= x^3+1
x^14 x^3+1 1001 9 x^13*x=(x^3+x^2+1)*x = x^4+x^3+x mod P(x) = x^3+1
x^15 1 0001 1 x^14*x = (x^3+1)*x = x^4+x mod P(x) = 1

五、伽罗华域上的运算

加法和减法:

加法和减法就是多项式计算中说的异或运算。

乘法和除法:

伽罗华域上的多项式乘法,其结果需要mod P(x),可以通过以下方式简化计算。

首先,考虑x8,x8 mod P(x) = P(x) – x8 = x4 + x3 +x2 +1 。

对于一般形式的多项式f(x)=a7x7 + a6x6 + a5x5 + a4x4 + a3x3 + a2x2 + a1x + a0,乘以x得到 xf(x) = (a7x8 + a6x7 + a5x6 + a4x5 + a3x4 + a2x3 + a1x1 + a0x) mod P(x)

这时有两种情况:

1)a7 == 0,此时结果是一个小于指数小于8的多项式,不需要取模计算。

2)a7 == 1,则将x8替换为x4 + x3 + x2 +1,而不用进行除法取模计算,结果为:

xf(x) = (x4 + x3 + x2 +1) + a6x7 + a5x6 + a4x5 + a3x4 + a2x3 + a1x1 + a0x

虽然可以通过替换减少除法计算,但还是过于复杂。尤其是在需要大量运算的场合,比如图像处理。于是牛人提出通过查表来减少计算。

六、查表的原理

首先介绍一个概念,生成元。

生成元是域上的一类特殊元素,生成元的幂可以遍历域上的所有元素。假设g是域GF(2w)上生成元,那么集合 {g0 ,g1 , ……,g(2w-1) } 包含了域GF(2w)上所有非零元素。在域GF(2w)中2总是生成元。

将生成元应用到多项式中, GF(2w)中的所有多项式都是可以通过多项式生成元g通过幂求得。即域中的任意元素a,都可以表示为a = gk

GF(2w)是一个有限域,就是元素个数是有限的,但指数k是可以无穷的。所以必然存在循环。这个循环的周期是2w-1(g不能生成多项式 0)。所以当k大于等于2w-1时,gk =gk%(2w-1)

对于gk = a,有正过程和逆过程。知道指数k求a是正过程,知道值a求指数k是逆过程。

对于乘法,假设a=gn,b=gm。那么ab = gn gm = gn+m。查表的方法就是根据a和b,分别查表得到n和m,然后查表gn+m即可。

因此需要构造正表和反表,在GF(2w)域上分别记为gflog和gfilog。gflog是将二进制形式映射为多项式形式,gfilog是将多项式形式映射为二进制形式。

注意:多项式0 ,是无法用生成元生成的。g0等于多项式1,而不是 0。

根据上文的GF(24)的元素表示,生成gflog表和gfilog表如下:

查表进行乘法和除法运算的例子

在GF(24)域上的乘法和除法,已知2w-1 = 24 -1 = 15:

乘法:

7 * 9 = gfilog[gflog[7] + gflog[9]] = gfilog[10 + 14] = gfilog[24 mod 15] = gfilog[9] = 10

除法:

13 / 11 = gfilog[gflog[13] - gflog[11]] = gfilog[13 - 7] = gfilog[6] = 12

五、生成GF(2w)gflog表和gfilog表的python代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
# coding=UTF-8  

# key : value => w : primitive_polynomial  
primitive_polynomial_dict = {4: 0b10011,                            # x**4  + x     + 1  
			    8: (1 << 8) + 0b11101,                 # x**8  + x**4  + x**3 + x**2 +1  
			    16: (1 << 16) + (1 << 12) + 0b1011,    # x**16 + x**12 + x**3 + x    + 1  
			    32: (1 << 32) + (1 << 22) + 0b111,     # x**32 + x**22 + x**2 + x    + 1  
			    64: (1 << 64) + 0b11011                # x**64 + x**4  + x**3 + x    + 1  
			    }  


def make_gf_dict(w):  
	gf_element_total_number = 1 << w  
	primitive_polynomial = primitive_polynomial_dict[w]  

	gfilog = [1]  # g(0) = 1  
	for i in xrange(1, gf_element_total_number - 1):  
		temp = gfilog[i - 1] << 1  # g(i) = g(i-1) * g  
		if temp & gf_element_total_number:  # if overflow, then mod primitive polynomial  
			temp ^= primitive_polynomial  # mod primitive_polynomial in GF(2**w) == XOR  
		gfilog.append(temp)  

	assert (gfilog[gf_element_total_number - 2] << 1) ^ primitive_polynomial  
	gfilog.append(None)  

	gflog = [None] * gf_element_total_number  
	for i in xrange(0, gf_element_total_number - 1):  
		gflog[gfilog[i]] = i  

	print "{:>8}\t{:>8}\t{:>8}".format("i", "gfilog[i]", "gflog[i]")  
	for i in xrange(0, gf_element_total_number):  
		print "{:>8}\t{:>8}\t{:>8}".format(i, gfilog[i], gflog[i])  


if __name__ == "__main__":  
	make_gf_dict(4)  

参考

http://blog.csdn.net/luotuo44/article/details/41645597 http://blog.csdn.net/mengboy/article/details/1514445 http://www.tuicool.com/articles/RZjAB3 http://ouyangmy.is-programmer.com/posts/41256.html