LOJ2462 「2018 集训队互测 Day 1」完美的集合

Link

\(\text{Part.1}\)

首先考虑dp求出完美集合的价值和。
设\(f_{i,j}\)表示考虑dfs序前\(i\)个点，集合重量限制为\(j\)的最大权值和。
假设当前节点为\(u\)，如果选择当前节点，那么接下来我们有可能选择\(u\)子树中的其它节点，即\(f_{i,j}+v_u\rightarrow f_{i+1,j+w_u}\)。
如果不选当前节点，那么我们应该跳过整个\(u\)子树，即\(f_{i,j}\rightarrow f_{i+size_u,j}\)。
最后的答案就在\(f_{n+1,m}\)。

\(\text{Part.2}\)

假如我们选定了一些完美集合，那么合法的放置测试装置的点也构成一个连通块。
为了保证计算不重不漏，我们可以计算所有点是合法放置点的方案数减去所有边的两端的点都是合法放置点的方案数。因为每个连通块都满足点数等于边数加一，因此恰好每种情况都只会被计算一遍。

\(\text{Part.3}\)

计算\({n\choose k}\pmod{5^{23}}\)。
类似于Ex-Lucas的过程，我们将要计算的转化为\(x=\alpha(5,n)\)和\(\frac{n!}{5^x}\bmod5^{23}\)，\(x\)的计算是很简单的。
\(n!\equiv 5^{\lfloor\frac n5\rfloor}\lfloor\frac n5\rfloor!\prod\limits_{j=1}^n[j\perp 5]j\pmod{5^{23}}\)
其中\(\lfloor\frac n5\rfloor!\)可以递归计算，记\(F_n(x)=\prod\limits_{i=1}^n[i\perp5](x+i)\)，那么我们现在需要求的就是\(F_n(0)\)。
我们知道\(F_{10n}(0)=F_{5n}(0)F_{5n}(5n)\)，而\((5n)^{23}\equiv0\pmod{5^{23}}\)。
因此多项式只有前\(23\)项系数是有用的，先递归计算\(F_{5n}(x)\)，然后利用二项式定理展开计算\(F_{5n}(x+5n)\)。除此之外还有\(O(1)\)个单项式，暴力乘即可。
这样计算一次\({n\choose k}\pmod{5^{23}}\)的时间复杂度为\(O(\log^2n)\)，大概会带\(10^3\)左右的常数。
一个常数优化是一开始可以\(\forall i\in[1,T]\)预处理\(F_i(x)\)。

总的时间复杂度为\(O(1000n^3+n^2m)\)，第一部分的常数不满。

#include<cstdio>
#include<vector>
#include<cstring>
#include<algorithm>
using i64=long long;
using pi=std::pair<int,int>;
using pii=std::pair<i64,i64>;
const int N=67,M=10007;const i64 P=11920928955078125;
int n,m,k,tim,W[N],V[N],dep[N],id[N],size[N];i64 lim,C[23][23];std::vector<pi>e[N];
int read(){int x;scanf("%d",&x);return x;}
void inc(i64&a,i64 b){a+=b-P,a+=a>>63&P;}
i64 mul(i64 a,i64 b){i64 d=a*b-(i64)((long double)a/P*b+0.5)*P;return d<0? d+P:d;}
i64 pow(i64 a,i64 b){i64 r=1;for(;b;b>>=1,a=mul(a,a))if(b&1)r=mul(a,r);return r;}
struct poly
{
    i64 a[23];
    poly(){memset(a,0,184);}
    poly(i64 x){a[0]=x,a[1]=1,memset(a+2,0,168);}
    i64&operator[](int x){return a[x];}
}pre[M];
poly operator*(poly a,poly b)
{
    poly c;
    for(int i=0;i<23;++i) for(int j=0;i+j<23;++j) inc(c[i+j],mul(a[i],b[j]));
    return c;
}
poly move(poly f,i64 k)
{
    static i64 pw[23];poly g;pw[0]=1,pw[1]=k;
    for(int i=2;i<23;++i) pw[i]=mul(pw[i-1],k);
    for(int i=0;i<23;++i) for(int j=0;j<=i;++j) inc(g[j],mul(mul(pw[i-j],C[i][j]),f[i]));
    return g;
}
poly get(i64 n)
{
    if(n<=10000) return pre[n];
    i64 k=n/10*10;poly f=get(k/2);f=f*move(f,k/2);
    for(i64 i=k+1;i<=n;++i) if(i%5) f=f*poly(i);
    return f;
}
pii fac(i64 n){i64 c=0,r=1;for(;n;n/=5)r=mul(r,get(n).a[0]),c+=n/5;return pii{r,c};}
i64 binom(i64 n){pii a=fac(n),b=fac(n-k),c=fac(k);return mul(mul(a.first,pow(5,a.second-b.second-c.second)),pow(mul(b.first,c.first),P/5*4-1));}
struct node
{
    i64 x,y;
    node(){x=-1e18,y=0;}
    node(i64 X,i64 Y){x=X,y=Y;}
    node cal(){return node(x,y<k? 0:binom(y));}
    node conj(){return node(x,P-y);}
};
node operator+(node&a,i64 b){return node(a.x+b,a.y);}
void update(node&a,node b){if(b.x>a.x)a.x=b.x,a.y=b.y;else if(a.x==b.x)inc(a.y,b.y);}
node solve(int w,int v,int flg)
{
    static node f[N][M];
    if(w>m) return node();
    for(int i=1;i<=tim+1;++i) std::fill(f[i]+1,f[i]+m+1,node());
    for(int i=w;i<=m;++i) f[1][i]=node(v,flg);
    for(int i=1,u;i<=tim;++i)
    {
	u=id[i];
	if(!flg||1ll*dep[u]*V[u]<=lim) for(int j=1;j+W[u]<=m;++j) update(f[i+1][j+W[u]],f[i][j]+V[u]);
	for(int j=1;j<=m;++j) update(f[i+size[u]][j],f[i][j]);
    }
    return f[tim+1][m].cal();
}
void dfs(int u,int fa){size[id[++tim]=u]=1;for(auto[v,w]:e[u]) if(v^fa) dep[v]=dep[u]+w,dfs(v,u),size[u]+=size[v];}
void bfs(int u,int fa){for(auto[v,w]:e[u]) if(v^fa) dep[v]=dep[u]+w,dfs(v,u);}
int main()
{
    n=read(),m=read(),k=read(),scanf("%lld",&lim),pre[0][0]=1;node ans;
    for(int i=0;i<23;++i) for(int j=C[i][0]=1;j<=i;++j) inc(C[i][j]=C[i-1][j],C[i-1][j-1]);
    for(int i=1;i<=10000;++i) pre[i]=i%5? pre[i-1]*poly(i):pre[i-1];
    for(int i=1;i<=n;++i) W[i]=read();
    for(int i=1;i<=n;++i) V[i]=read();
    for(int i=1,u,v,w;i<n;++i) u=read(),v=read(),w=read(),e[u].push_back({v,w}),e[v].push_back({u,w});
    for(int u=1;u<=n;++u)
    {
	bfs(u,dep[u]=tim=0),update(ans,solve(W[u],V[u],0)),update(ans,solve(W[u],V[u],1));
	for(auto[v,w]:e[u]) if(u<v&&1ll*w*std::max(V[u],V[v])<=lim) dep[u]=dep[v]=w,tim=0,bfs(u,v),bfs(v,u),update(ans,solve(W[u]+W[v],V[u]+V[v],1).conj());
    }
    printf("%lld",ans.y);
}