自然数前段的均衡样本

mathe

对于偶数的情况，偶数项和奇数项的特征方程是一样的，但表达式不是(x-1)SK(x)
比如：K=4

偶数项的特征方程是：
$(-1+x)^4 (1+x+x^2)$
wayne 发表于 2010-4-30 23:54

$(-1+x)^4 (1+x+x^2)$
是
$(1-x)(1-x^2)^2(1-x^3)$的因子，是不是K是偶数时，特征方程都是
$(1-x)(1-x^2)^2(1-x^3)...(1-x^{K-1})$的因子呢（或者对于K>4相等？）

wayne

算了一下14元的

$(x^2-1)Y_14(x)$

补充一下代码：

1. data = Block[{k = 14}, Table[SeriesCoefficient[Series[Product[(x^s - x^(n + 1))/(1-x^s), {s, 1,k}], {x,0,n k}], (n+1) k/2], {n, k,200}]]

复制代码

1. d = FindLinearRecurrence[data]; dlen = Length[d]; (x^dlen - d.x^Range[dlen - 1, 0, -1])/Product[x^ii - 1, {ii, 1, 11}] // Simplify

复制代码

mathe

所以看来都是$(x^2-1)Y_K(x)$的因子

mathe

于是现在我们可以来如下计算母函数，假设数列为
$a_0,a_1,...,a_n,...$
我们计算
$(a_0+a_1*x+...+a_n*x^n+..)*(x^2-1)Y_K(x)$
对于得到的乘法结果，我们只取所有次数低于$(x^2-1)Y_K(x)$的最高次数的部分，得到$T_K(x)$
那么${T_K(x)}/{(x^2-1)Y_K(x)}$就是结果了

wayne

81# mathe

是的，这个可以根据80楼给出的表达式，把系数方程倒过来乘以x^2-1一下即可

mathe

不过我计算过K=6时所有项（同时包含奇数和偶数项）时的母函数为
$ {x^10 + x^8 + 3*x^7 + 4*x^6 + 4*x^5 + 4*x^4 + 3*x^3 + x^2 + 1}/{x^17 - x^16 - 2*x^15 + x^14 + x^13 + x^12 + x^11 - 2*x^9 - 2*x^8 + x^6 + x^5 + x^4 + x^3 - 2*x^2 - x + 1}$
好像不是很规律

wayne

我算过通项公式，很复杂，所以，我感觉母函数也不会很简洁的

wayne

算了16元的特征方程,把前面的总结一下，偶数的规律应该可以确定了：

4元： $ (-1+x)Y_4(x)$
6元： $ (x^2-1)Y_6(x)$
8元：$\frac{-1+x}{1-x+x^2}Y_8(x)$
10元：$ (x^2-1)Y_10(x)$
12元：$\frac{-1+x}{1-x+x^2-x^3+x^4}Y_12(x)$
14元：$(x^2-1)Y_14(x)$
16元：$\frac{-1+x}{1-x+x^2-x^3+x^4-x^5+x^6}Y_16(x)$

mathe

我也计算了以下，发现对于偶数情况，应该还能够化简，
比如
4元:${x-1}/{x+1}Y_4(x)$
6元:${x-1}/{(x+1)(x^2+1)}Y_6(x)$
8元:${x-1}/{(x+1)(x^2+1)(x^3+1)}Y_8(x)$
10元:${x-1}/{(x+1)(x^2+1)(x^3+1)(x^4+1)}Y_10(x)$
12元:${x-1}/{(x+1)(x^2+1)(x^3+1)(x^4+1)(x^5+1)}Y_12(x)$
14元:${x-1}/{(x+1)(x^2+1)(x^3+1)(x^4+1)(x^5+1)(x^6+1)}Y_14(x)$
16元:${x-1}/{(x+1)(x^2+1)(x^3+1)(x^4+1)(x^5+1)(x^6+1)(x^7+1)}Y_16(x)$
规律非常之强

mathe

1. 
   
2. S(K,n)=
   
3. {
   
4.     local(V,L,s,R);
   
5.     L=K*n;
   
6.     V=matrix(K,L);
   
7.     R=vector(n);
   
8.     V[1,1]=1;
   
9.     R[1]=0;
   
10.     for(u=2,2*n-1,
    
11.        for(h=1,K-1,
    
12.           for(v=1,L,
    
13.             s=V[K+1-h,v];
    
14.             if(v>u,s=s+V[K-h,v-u]);
    
15.             V[K+1-h,v]=s
    
16.           )
    
17.        );
    
18.        V[1,u]=1;
    
19.        if(bitand(u,1)==1,
    
20.             s=(u+1)/2;
    
21.             R[s]=V[K,K*s]
    
22.        )
    
23.     );
    
24.     R
    
25. }
    
26. 
    
27. EstPol(V,K)=
    
28. {
    
29.     local(len,p,r);
    
30.     len = length(V);
    
31.     p=V[1];
    
32.     for(u=2,len,
    
33.        p=p+V[u]*x^(u-1)
    
34.     );
    
35.     p=p*(1-x^2);
    
36.     for(u=1,K-1,
    
37.        p=p*(1-x^u)
    
38.     );
    
39.     r=polcoeff(p,0);
    
40.     for(u=1,len-1,
    
41.       r=r+polcoeff(p,u)*x^u
    
42.     );
    
43.     r
    
44. }
    
45. MS(k,n)=EstPol(S(k,n),k)
    

复制代码

这个代码可以计算出84#定义的$T_K(x)$
然后我们将$T_K(x)$和$(1-x^2)Y_K(x)$的公因子去掉就可以了。
其中计算时，n需要充分大，大于${K*(K+1)}/2$就可以了