A075464

mathe

添加A075463的代码，运行程序时只要输入一个参数n即可，不要管程序的提示

zhouguang

这个帖子钩起了我旧时的回忆。那时候，在智星看见了mathe，在东陆看见了hujunhua，呵呵。
为了纪念，关于这个问题写了个pdf，请大家看看。呵呵。

zhouguang

mathe在https://oeis.org/A159257中提到：
“Let f(k, 2x) be the Chebyshev Polynomial of the Second Kind in Limit Field F_2. So f(0,x)=1, f(1,x)=x, f(2,x)=(1+x)^2, f(n+1,x)=x*f(n,x)+f(n-1,x). The order of gcd(f(n,x), f(n, 1+x)) is this sequence. For example, f(5,x)=x^5+x=x(1+x)^4. f(5, 1+x)=x^4(1+x). So the GCD of them is x(1+x) and the order is 2. So the 5th element of the array is 2. - Du, Zhao Hui, Mar 17 2014”
而上面附件pdf中第3页提到的“（多项式d 和M 的特征多项式相关，是有名字的，呵呵）。”中d的名字就是“the Chebyshev Polynomial of the Second Kind in Limit Field F_2”。
下面代码中对639的操作如下图，而那第一行右面长长的空白部分，就是mathe提到的2^h的自由度。我想mathe的代码会比M快很多。呵呵。

1. 
   
2. cc[n_] := Module[{ds, sds, x, dn, sdn, ta, tb, a, b, ib, ms, m0, ans, js}, ds = {1, x}; sds = {1, 1 + x};
   
3.    Do[AppendTo[ds, PolynomialMod[x ds[[-1]] + ds[[-2]], 2]]; AppendTo[sds, PolynomialMod[Last[sds] + Last[ds], 2]];, {n - 1}];
   
4.    dn = ds /. x -> 2; sdn = sds /. x -> 2; ta = Reverse[IntegerDigits[dn[[n + 1]], 2]]; tb = Append[Reverse[IntegerDigits[sdn[[n]], 2]], 0];
   
5.    a = Table[0, {n}, {n}]; b = Table[0, {n}]; ib = Table[1, {n}]; ms = IdentityMatrix[n]; m0 = Table[If[Abs[i - j] <= 1, 1, 0], {i, n}, {j, n}];
   
6.    Do[If[ta[[i]] == 1, a = Mod[a + ms, 2]]; If[tb[[i]] == 1, b = Mod[b + ms.ib, 2]]; ms = Mod[ms.m0, 2], {i, n + 1}];
   
7.    ans = LinearSolve[a, b, Modulus -> 2];  js = {ans, 1 - Mod[m0.ans, 2]};
   
8.    Do[AppendTo[js, 1 - Mod[js[[-1]] + Append[Rest[js[[-1]]], 0] +  Prepend[Most[js[[-1]]], 0] + js[[-2]], 2]];, {n - 2}]; js];
   
9. cc[639] // MatrixPlot
   

复制代码

mathe

A075464第65项算出来了，是1867,计算机单CPU计算了一个星期

zhouguang

mathe 发表于 2014-3-29 16:58
A075464第65项算出来了，是1867,计算机单CPU计算了一个星期

厉害，呵呵。
打算计算79么？

mathe

这个艰巨的任务留给你如何？