• 【PE806】Nim on Towers of Hanoi(DP)(生成函数)


    Nim on Towers of Hanoi

    题目链接:PE806

    题目大意

    一个有 n 个盘子的汉诺塔,在第 i 个状态的时候如果三个柱子的盘子个数的异或和是 0,就会给 i 的贡献。
    求 n=100000 时候的贡献和。

    思路

    发现这个给贡献很难搞,因为你几乎很难确定它是在第几步。
    但是发现这个很诺塔它有一个对称性,也就是把第 i i i 个状态的两边的柱子换一下位置,就是第 2 n − 1 − i 2^n-1-i 2n1i 个状态了。

    这个有什么用呢?因为是异或,所以如果 i i i 满足, 2 n − 1 − i 2^n-1-i 2n1i 也满足。
    所以我们只需要求出满足的个数 k k k,答案就是 k ( 2 n − 1 ) 2 \dfrac{k(2^n-1)}{2} 2k(2n1)

    于是问题变成求 k k k,条件是 a + b + c = 0 , a ⊕ b ⊕ c = 0 a+b+c=0,a\oplus b\oplus c=0 a+b+c=0,abc=0
    会发现如果 n n n 的某一位是 1 1 1,它一定要拆成比它第一位的两个数(所以奇数没有解)。
    那分配给那两个呢?有 3 3 3 种选择的。
    所以满足条件的 a , b , c a,b,c a,b,c 的情况个数就是 3 p o p c o u n t ( n ) 3^{{\rm popcount}(n)} 3popcount(n)

    看起来挺多,但是 n n n p o p c o u n t ( n ) = 6 {\rm popcount}(n)=6 popcount(n)=6,所以其实很少。
    但是问题是,一个 a , b , c a,b,c a,b,c 不一定只对应一种情况,所以你还要看一个三元组 ( a , b , c ) (a,b,c) (a,b,c) 对应多少方案。
    设为 f i , j , k f_{i,j,k} fi,j,k,生成函数的话就是 F ( x , y , z ) F(x,y,z) F(x,y,z)
    (由于对称你有 f i , j , k = f k , j , i f_{i,j,k}=f_{k,j,i} fi,j,k=fk,j,i

    考虑怎么 DP,然后你每次考虑多一个最大的圆盘。
    那应该是把上面的一坨丢到第二个碟子,然后把最小面的丢过去,然后再把第二个柱子那一坨再扔到第三个柱子。
    那中间的两个过程都会分别给贡献,所以你最大的柱子可以放在第一个柱子的最下面,也可以放在第三个柱子的最小面,分别会给贡献。
    那就要看上面那坨是怎样的,所以有这样的转移:
    ( x − 1 , z , y ) → ( x , y , z ) (x-1,z,y)\rightarrow(x,y,z) (x1,z,y)(x,y,z)(上面那坨先到第三个,再到第二个)
    ( y , z , x − 1 ) → ( x , y , z ) (y,z,x-1)\rightarrow(x,y,z) (y,z,x1)(x,y,z)(上面的步骤反过来)
    那在生成函数上:
    F ( x , y , z ) = 1 + x F ( x , z , y ) + x F ( y , z , x ) F(x,y,z)=1+xF(x,z,y)+xF(y,z,x) F(x,y,z)=1+xF(x,z,y)+xF(y,z,x) 1 1 1 是修正项)
    然后由于 F ( x , y , z ) = F ( z , y , x ) F(x,y,z)=F(z,y,x) F(x,y,z)=F(z,y,x),我们用中间的做标志,设为 F y F_y Fy

    然后可以得到三个式子:
    F x = 1 + z F y + y F z F_x=1+zF_y+yF_z Fx=1+zFy+yFz
    F y = 1 + z F x + x F z F_y=1+zF_x+xF_z Fy=1+zFx+xFz
    F z = 1 + y F x + x F y F_z=1+yF_x+xF_y Fz=1+yFx+xFy
    解方程可以得到三个的表示,我们只需要用到 F y F_y Fy
    F y = ( 1 + y ) ( 1 + x + z − y ) 1 − x 2 − y 2 − z 2 − 2 x y z F_y=\dfrac{(1+y)(1+x+z-y)}{1-x^2-y^2-z^2-2xyz} Fy=1x2y2z22xyz(1+y)(1+x+zy)

    那我们只需要取它 [ x a y b z c ] [x^ay^bz^c] [xaybzc] 项的系数即可!

    考虑上面的拆开每个系数分开算,下面的话我们设 p = x 2 + y 2 + z 2 + 2 x y z p=x^2+y^2+z^2+2xyz p=x2+y2+z2+2xyz
    那下面的 1 1 − p = 1 + p + p 2 + . . . \dfrac{1}{1-p}=1+p+p^2+... 1p1=1+p+p2+...
    然后我们看到只有一个有多元,考虑枚举它的个数,那剩下三个用多少次也就固定了,组合数安排一下顺序就能求。

    那我们看看复杂度,上面枚举 ( a , b , c ) (a,b,c) (a,b,c) 3 p o p c o u n t ( n ) 3^{{\rm popcount(n)}} 3popcount(n),这里算有一个 2 × 4 2\times 4 2×4 的常数(上面的系数每个都要),然后里面再一个 n n n,所以复杂度是 3 p o p c o u n t ( n ) n 3^{{\rm popcount(n)}}n 3popcount(n)n,在 n = 100000 n=100000 n=100000 的时候 p o p c o u n t ( n ) = 6 {\rm popcount(n)}=6 popcount(n)=6,跑的很快。

    代码

    #include
    #include
    #define mo 1000000007
    
    using namespace std;
    
    const int N = 1e5 + 100;
    int n = 100000, ans;
    int xl[10], jc[N], inv[N], invs[N];
    
    int add(int x, int y) {return x + y >= mo ? x + y - mo : x + y;}
    int dec(int x, int y) {return x < y ? x - y + mo : x - y;}
    int mul(int x, int y) {return 1ll * x * y % mo;}
    int ksm(int x, int y) {int re = 1; while (y) {if (y & 1) re = mul(re, x); x = mul(x, x); y >>= 1;} return re;}
    
    int slove_(int a, int b, int c) {
    	if (a < 0 || b < 0 || c < 0) return 0;
    	if (((a & 1) ^ (b & 1)) || ((a & 1) ^ (c & 1))) return 0;
    	int ans = 0;
    	for (int i = 0, di = 1; i <= min(min(a, b), c); i++, di = mul(di, 2)) {
    		if ((a - i) & 1) continue;
    		int n = (a - i + b - i + c - i) / 2 + i;
    		ans = add(ans, mul(di, mul(jc[n], mul(invs[i], mul(invs[(a - i) / 2], mul(invs[(b - i) / 2], invs[(c - i) / 2]))))));
    	}
    	return ans;
    }
    
    int slove(int a, int b, int c) {
    	int ans = 0;
    	ans = add(ans, slove_(a, b, c));
    	ans = add(ans, slove_(a - 1, b, c));
    	ans = add(ans, slove_(a, b, c - 1));
    	ans = dec(ans, slove_(a, b - 1, c));
    	ans = add(ans, slove_(a, b - 1, c));
    	ans = add(ans, slove_(a - 1, b - 1, c));
    	ans = add(ans, slove_(a, b - 1, c - 1));
    	ans = dec(ans, slove_(a, b - 2, c));
    	return ans;
    } 
    
    void dfs(int now, int a, int b, int c) {
    	if (now > xl[0]) {
    		ans = add(ans, slove(a, b, c));
    		return ;
    	}
    	dfs(now + 1, a + (1 << (xl[now] - 1)), b + (1 << (xl[now] - 1)), c);
    	dfs(now + 1, a + (1 << (xl[now] - 1)), b, c + (1 << (xl[now] - 1)));
    	dfs(now + 1, a, b + (1 << (xl[now] - 1)), c + (1 << (xl[now] - 1)));
    }
    
    int main() {
    	jc[0] = 1; for (int i = 1; i < N; i++) jc[i] = mul(jc[i - 1], i);
    	inv[0] = inv[1] = 1; for (int i = 2; i < N; i++) inv[i] = mul(inv[mo % i], mo - mo / i);
    	invs[0] = 1; for (int i = 1; i < N; i++) invs[i] = mul(invs[i - 1], inv[i]);
    	
    	if (n & 1) {printf("0"); return 0;}
    	for (int i = 0; i <= 20; i++)
    		if ((n >> i) & 1) xl[++xl[0]] = i;
    	dfs(1, 0, 0, 0);
    	
    	int val = dec(ksm(2, n), 1);
    	printf("%d", mul(mul(ans, val), (mo + 1) / 2));
    	
    	return 0;
    }
    

    答案自己编译一下就有。

  • 相关阅读:
    【算法题解】拓扑序计数+树形DP
    面试常问CAS和ABA问题,你懂了吗?
    LabVIEW中NIPackageManager功能介绍
    【ACWing】1401. 围住奶牛
    多因素方差分析(Multi-way Analysis of Variance) R实现
    全真互联不同于元宇宙,五大技术大公开
    R语言Sys.Date函数获取当前日期、获取指定日期自1970年1月1日以来经过的天数(Julian Date)
    基于BM1684X 架构实现 Faiss 的两个查询接口
    城市内涝解决方案:实时监测,提前预警,让城市更安全
    ubuntu安装和启动redis命令步骤及其配置文件redis.conf
  • 原文地址:https://blog.csdn.net/weixin_43346722/article/details/127086733