将过程逆序，问题即转换为以下形式——
从$(a,b,c)$出发，每步移动到周围六个格子之一，求$d$步内不离开长方体的方案数
显然每一维可以通过生成函数合并，不妨仅考虑其中一维，问题也即
从$(0,a)$出发，每步移动到右上/右下的格子，$\forall 0\le i\le d$求$i$步内与$y=0,n+1$不交的方案数
可以参考[hdu7069]，容斥得答案即
$$
f(a)-\sum_{l\ge 1,l\equiv 1(mod\ 2)}\Big(f(n'-ln'-a)+f(ln'+n'-a)\Big)+\sum_{l\ge 1,l\equiv 0(mod\ 2)}\Big(f(ln'+a)+f(a-ln')\Big)
$$
（其中$n'=n+1,f(a)$表示$i$步从$(0,a)$到$\forall j\in [1,n],(i,j)$的方案数）
关于$f(a)$，考虑枚举向上的步数$x$，不难得到
$$
f(a)=\sum_{0\le x\le i,1\le a+2x-i\le n}{i\choose x}=\sum_{\max(\lceil\frac{i-a+1}{2}\rceil,0)\le x\le \min(\lfloor\frac{i+n-a}{2}\rfloor,i)}{i\choose x}
$$
要求其范围非空，即解得$1-i\le a\le i+n$
换言之，$l$的范围是$o(\frac{i}{n})$的，暴力枚举的复杂度为$o(\frac{d^{2}}{n})$
另外，关于$f(a)$的计算，可以调换$i$和$l$的枚举顺序，并在$i$增加时维护区间
另一方面，当$d$较小时，可以直接暴力$o(nd)$递推
按$n\le \sqrt{d}$和$n>\sqrt{d}$分别处理，时间复杂度为$o(d\sqrt{d})$
结合合并三维时的生成函数，总复杂度为$o(d\sqrt{d}+d\log d)$，可以通过
1 #include<bits/stdc++.h> 2 using namespace std; 3 #define N 120005 4 #define M (1<<19) 5 #define K 400 6 #define mod 998244353 7 #define ll long long 8 int n,a,b,c,a0,b0,c0,L,R,s,ans,fac[N],inv[N],rev[M],A[M],B[M],C[M],mi[M],f[N][K]; 9 int qpow(int n,int m){ 10 int s=n,ans=1; 11 while (m){ 12 if (m&1)ans=(ll)ans*s%mod; 13 s=(ll)s*s%mod,m>>=1; 14 } 15 return ans; 16 } 17 void ntt(int *A,int n,int p=0){ 18 for(int i=0;i<n;i++) 19 if (i<rev[i])swap(A[i],A[rev[i]]); 20 for(int i=2;i<=n;i<<=1){ 21 int s=qpow(3,(mod-1)/i); 22 if (p)s=qpow(s,mod-2); 23 mi[0]=1; 24 for(int k=1;k<(i>>1);k++)mi[k]=(ll)mi[k-1]*s%mod; 25 for(int j=0;j<n;j+=i) 26 for(int k=0;k<(i>>1);k++){ 27 int x=A[j+k],y=(ll)A[j+k+(i>>1)]*mi[k]%mod; 28 A[j+k]=(x+y)%mod,A[j+k+(i>>1)]=(x-y+mod)%mod; 29 } 30 } 31 if (p){ 32 int s=qpow(n,mod-2); 33 for(int i=0;i<n;i++)A[i]=(ll)A[i]*s%mod; 34 } 35 } 36 int Com(int n,int m){ 37 return (ll)fac[n]*inv[m]%mod*inv[n-m]%mod; 38 } 39 void Calc(int *A,int a,int a0,int p=0){ 40 int i0=max(max(1-a0,a0-a),0); 41 L=max((i0-a0>>1)+1,0),R=min((i0+a-a0>>1),i0),s=0; 42 for(int x=L;x<=R;x++)s=(s+Com(i0,x))%mod; 43 if (!p)A[i0]=(A[i0]+s)%mod; 44 else A[i0]=(A[i0]-s+mod)%mod; 45 for(int i=i0+1;i<=n;i++){ 46 int L0=max((i-a0>>1)+1,0),R0=min((i+a-a0>>1),i); 47 s=(s<<1)%mod,s=(s-Com(i-1,R)+mod)%mod; 48 if (L)s=(s+Com(i-1,L-1))%mod; 49 while (L<L0)s=(s-Com(i,L++)+mod)%mod; 50 while (R<R0)s=(s+Com(i,++R))%mod; 51 if (!p)A[i]=(A[i]+s)%mod; 52 else A[i]=(A[i]-s+mod)%mod; 53 } 54 } 55 void calc(int *A,int a,int a0){ 56 if ((ll)a*a<=n){ 57 memset(f,0,sizeof(f)); 58 A[0]=f[0][a0]=1; 59 for(int i=1;i<=n;i++) 60 for(int j=1;j<=a;j++){ 61 if (j>1)f[i][j]=(f[i][j]+f[i-1][j-1])%mod; 62 if (j<n)f[i][j]=(f[i][j]+f[i-1][j+1])%mod; 63 A[i]=(A[i]+f[i][j])%mod; 64 } 65 return; 66 } 67 int n0=a+1; 68 Calc(A,a,a0); 69 for(int l=1;;l+=2){ 70 int k=n0-l*n0-a0; 71 if (k<1-n)break; 72 Calc(A,a,k,1); 73 } 74 for(int l=1;;l+=2){ 75 int k=l*n0+n0-a0; 76 if (k>n+a)break; 77 Calc(A,a,k,1); 78 } 79 for(int l=2;;l+=2){ 80 int k=l*n0+a0; 81 if (k>n+a)break; 82 Calc(A,a,k); 83 } 84 for(int l=2;;l+=2){ 85 int k=a0-l*n0; 86 if (k<1-n)break; 87 Calc(A,a,k); 88 } 89 } 90 int main(){ 91 fac[0]=inv[0]=inv[1]=1; 92 for(int i=1;i<N;i++)fac[i]=(ll)fac[i-1]*i%mod; 93 for(int i=2;i<N;i++)inv[i]=(ll)(mod-mod/i)*inv[mod%i]%mod; 94 for(int i=1;i<N;i++)inv[i]=(ll)inv[i-1]*inv[i]%mod; 95 for(int i=0;i<M;i++)rev[i]=(rev[i>>1]>>1)+((i&1)*(M>>1)); 96 scanf("%d%d%d%d%d%d%d",&n,&a,&b,&c,&a0,&b0,&c0); 97 calc(A,a,a0),calc(B,b,b0),calc(C,c,c0); 98 for(int i=0;i<=n;i++){ 99 A[i]=(ll)A[i]*inv[i]%mod; 100 B[i]=(ll)B[i]*inv[i]%mod; 101 C[i]=(ll)C[i]*inv[i]%mod; 102 } 103 ntt(A,M),ntt(B,M),ntt(C,M); 104 for(int i=0;i<M;i++)A[i]=(ll)A[i]*B[i]%mod*C[i]%mod; 105 ntt(A,M,1),ans=(ll)A[n]*fac[n]%mod; 106 printf("%d\n",ans); 107 return 0; 108 }