北京Day 13

今日得分：85+0+0。

T1

题目大意：给定n个点n条有向边，保证每个点的入度和出度都为1，求随机选定k个点，从这k个点出发能遍历整个图的概率。答案对998244353取模。n,k<=152501。

题解：

显然题目给定的图结构是一个置换的若干循环，而每个循环只要分配 1 个点即可。我们可以预处理所有循环的大小，记为ai。

其实这题存在乱搞做法。

对于k比较小的数据可以背包，f(i,j)表示第i个环，已经分配j个点的合法方案数，f(i,j)=Σf(i-1,j-l)*C(a i,l)。时间复杂度O(cnt*k)。

对于k比较大的数据可以容斥，合法方案=C(n,k)-ΣC(n-ai,k)+ΣC(n-ai-aj,k)-…… 用背包预处理一下每个数对答案的贡献，最后一起统计答案即可，时间复杂度O(cnt*(n-k))。

真正正解NTT。

考虑每个循环的生成函数 A(x)，x i 的系数为该循环分配 i 个点的方案数贡献。我们将所有的 A(x) 相乘得到 G(x)，则 G(x) 的 x i 系数即为分配 i 个点满足任意一个循环均有至少 1 个点的方案数。

经过一些分析可以得到 A(i) = ΣC(ai,j) x^j。考虑不同的 a i 至多有 sqrt(n))种，我们可以对于不同的 a i 使 用 NTT 计算 A(x) 的点值再快速幂合并起来。

接着对上述做法进行优化。

我们在使用 NTT 合并两个大小为 Sa 与 Sb 的循环时，复 杂度其实是 O((Sa + Sb)log(Sa + Sb)) 的。

这与哈夫曼编码的过程其实很相似，我们使用哈夫曼编码的顺序 来合并，可以保证合并的复杂度最优。

另一种的合并方法是每轮将相邻的两个多项式合并，这样每轮的 复杂度均为 O(log n)，而至多进行 O(log n) 轮（因为每轮多项式 数量减半）

可以通过本题。

T1AC代码（乱搞）

​
#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
#include<algorithm>
#include<queue>
#include<vector>
#include<cmath>
#include<set>
#include<bitset>
using namespace std;
inline int re_ad()
{
	int x=0,f=1;char ch=getchar();
	while(ch>'9'||ch<'0'){if(ch=='-')f=-1;ch=getchar();}
	while(ch>='0'&&ch<='9')x=x*10+ch-48,ch=getchar();
	return x*f;
}
int T,n,k,P[160010],pre[160010],ni[160010],num[160010],cnt,sum[160010];
unsigned long long a[2][160010],ans,f1[3010],f2[3010];
bool vis[160010];
const int mo=998244353;
inline long long ksm(int a,int b){int ret=1;while(b){if(b&1)ret=1ll*ret*a%mo;b>>=1;a=1ll*a*a%mo;}return ret;}
inline long long c(int x,int y){if(y>x)return 0;return 1ll*pre[x]%mo*ni[y]%mo*ni[x-y]%mo;}
inline long long solve(int k)
{
	register int su=sum[cnt];
	long long ret=0;
	for(register int i=0;i<=n-k;i++)
	ret+=1ll*c(su-i,k)*(a[0][i]-a[1][i]),ret%=mo;
	return ret;
}
signed main()
{
	register int i,j,l,x;
	T=re_ad();
	while(T--)
	{
	n=re_ad();k=re_ad();
	memset(vis,0,sizeof(vis));
	memset(num,0,sizeof(num));
	memset(a,0,sizeof(a));
	memset(sum,0,sizeof(sum));
	cnt=0;ans=0;
	for(i=1;i<=n;++i)P[i]=re_ad();
	pre[0]=1;
	for(i=1;i<=n;++i)pre[i]=1ll*pre[i-1]*i%mo;
	ni[n]=ksm(pre[n],mo-2);
	for(i=n-1;i>=0;--i)ni[i]=1ll*ni[i+1]*(i+1)%mo;
	for(i=1;i<=n;++i)
	{
	if(!vis[i])
	{++cnt;++num[cnt];vis[i]=1;x=P[i];while(!vis[x])vis[x]=true,++num[cnt],x=P[x];sum[cnt]=sum[cnt-1]+num[cnt];}
	}
	if(k<=3000)
	{
	memset(f1,0,sizeof(f1));f1[0]=1;
	for(i=1;i<=cnt;++i)
	{
	swap(f1,f2);memset(f1,0,sizeof(f1));
	for(j=i;j<=k-cnt+i;++j)for(l=1;l<=min(j,num[i]);++l)f1[j]=(1ll*f1[j]+(f2[j-l])*1ll*c(num[i],l)%mo)%mo;
	}
	cout<<1ll*f1[k]*ksm(c(n,k),mo-2)%mo<<endl;
	}
	else {
	a[0][0]=1;
	for(i=1;i<=cnt;++i)
	{
	for(j=n-k;j>=num[i];--j)a[0][j]=(a[0][j]+a[1][j-num[i]])%mo,a[1][j]=(a[1][j]+a[0][j-num[i]])%mo;
	}
	ans=solve(k)%mo;
	ans=(ans+mo)%mo;
	cout<<ans*ksm(c(n,k),mo-2)%mo<<endl;}
	}
	
}

​

T2

题目大意：一个行数为 M、列数为 N 的方格图，有 S1+S2 种不同的地砖，其中 S1 种地砖是横着的 1*2 型地砖，S2 种地砖是竖着的 2*1 型地砖。用这些地砖把她的地板铺满，要求地砖之间不能有重叠，将 M*N 的地板铺满的 方案数，记为 F(N) 。求F(L)+F(L+1)+F(L+2)+……+F(R)的结果，答案对 998244353 取模。1<=M<= 8，1 <= L<=R<=10^2501，1 <=S1,S2<=10^9。

题解：

对于列数较小的情况，这是一个经典的轮廓线 DP 问题。 

对于列数较大的情况，设 f(i, S) 表示前 i − 1 列已经填满，第 i 列的状态为 S 的方案数。 预处理 f(i, S) 到 f(i + 1, S2) 的转移矩阵，然后快速幂优化（需要 写十进制快速幂）。

实际上有用的状态数不足 2^m，把从初始状态不可能到达与不可能 从此状态出发到达最终状态的状态去除，当 m = 8 时仅有 71 个 有用状态。

接下来考虑优化。

矩阵乘法本身的复杂度很难得到优化（使用分治矩阵乘法优化并不明显）。

我们考虑一个 k × k 的矩阵 A 的特征多项式 f(x) = |λE − A|。

根据 Cayley-Hamilton 定理有 f(A) = 0，故有 A r ≡ (A r mod f(A))。

这样只需要求出 f(A) 后计算 A r mod f(x)就可以在 O(k^4 ) 的时间内求出来 A k。 考虑到直接求元素是多项式的矩阵的行列式很麻烦，我们可以先 随便代入 k + 1 个不同的求出的 k + 1 个点值，最后通过拉格朗日 插值还原出 f(x)。

T2代码咕了。

T3

题目大意：n个节点的树，每个节点有点权，每条边有边权，每次操作修改一个点的点权或者询问x到y路径上的最小代价。代价的定义是选定一个点，路径上所有的点的点权乘以到达这个选定点的边权的总和。n<=152501。

题解：对于每个询问，把路径上所有点的信息取出，问题转化为数轴上有若干个点，每个位置有 Ai 个点，你要在其中选择一个位置作为 集合点，让所有点到该位置集合的距离之和最小。

显然应该取这若干个点的中位数位置。

用数据结构优化这一过程，我们需要解决的问题有 3 个： 更改一个节点的权值。 查询一段路径 u → v 上带权的中位数节点。 查询一段路径 u → v 上所有点到某一点 x 带权距离之和。

对于问题 1、2，我们可以考虑在 ST 表上倍增，使用另一个数据 结构带修改维护路径权值和就可以解决。也可以使用 LCT 维护树结构，在 splay 树上查找答案。

对于问题 1、3，把路径拆成 3 段：u → x，x → lca，lca → v，每 一段都可以使用 DFS 序＋树状数组或树链剖分＋线段树解决。

T3AC代码

#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
#include<algorithm>
#include<queue>
#include<vector>
#include<cmath>
#include<set>
#include<bitset>
using namespace std;
inline int re_ad()
{
	int x=0,f=1;char ch=getchar();
	while(ch>'9'||ch<'0'){if(ch=='-')f=-1;ch=getchar();}
	while(ch>='0'&&ch<='9')x=x*10+ch-48,ch=getchar();
	return x*f;
}
inline long long mi(long long x,long long y){return x<y?x:y;}
void pu_t(long long x){if(x>9)pu_t(x/10);putchar(x%10+48);}
const long long inf=100000000000000000;
int n,q,a[160000],fa[20][160010],d[160010],tot=0,dfn[160010],en[160010];
long long ans,dep[160010];
struct node{int to,cost;};
vector<node>g[160010];
long long t1[160010],t2[160010];
inline void chan(long long t[],int y,long long x){while(y<=n)t[y]+=x,y+=y&-y;}
inline long long que(long long t[],int x){long long ret=0;while(x)ret+=t[x],x-=x&-x;return ret;}
void dfs(int x)
{
	for(register int i=1;i<=19;++i)
	fa[i][x]=fa[i-1][fa[i-1][x]];
	dfn[x]=++tot;
	for(register int i=0;i<g[x].size();++i)
	{
	int v=g[x][i].to;
	if(fa[0][x]==v)continue;
	fa[0][v]=x;d[v]=d[x]+1;dep[v]=dep[x]+g[x][i].cost;
	dfs(v);
	}en[x]=tot;
}
inline int lca(int x,int y)
{
	if(d[x]<d[y])swap(x,y);
	register int i=0,t=d[x]-d[y];
	for(;t;++i)
	{
	if(t&1){x=fa[i][x];}
	t>>=1;
	}
	if(x==y)return x;
	for(i=19;i>=0;--i)
	{
	if(fa[i][x]!=fa[i][y])
	{
	x=fa[i][x];y=fa[i][y];
	}
	}
	return fa[0][x];
}
inline long long ask1(int x,int y){int z=lca(x,y);return que(t1,dfn[x])+que(t1,dfn[y])-2ll*que(t1,dfn[z])+a[z];}
inline long long ask2(int x,int y){int z=lca(x,y);return que(t2,dfn[x])+que(t2,dfn[y])-2ll*que(t2,dfn[z])+a[z]*dep[z];}
inline long long ask(int x,int y)
{
	register int i,lc=lca(x,y),tp,t,an;
	long long su=ask1(x,y),zj=(su+1)>>1,ret=0;
	if(a[x]>=zj){tp=1;an=x;}
	else if(su-a[y]<zj){tp=2;an=y;}
	else
	{
	if(ask1(x,lc)>=zj)
	{t=x;
	for(i=19;i>=0;--i){if(d[fa[i][t]]>=d[lc]&&(ask1(fa[i][t],x)<zj))t=fa[i][t];an=fa[0][t];tp=1;}}
	else
	{t=y;
	for(i=19;i>=0;--i){if(d[fa[i][t]]>=d[lc]&&(ask1(fa[i][t],y)<zj))t=fa[i][t];an=fa[0][t];tp=2;}
	}
	}
	if(tp==2)swap(x,y);
	ret=ask2(x,an)-ask1(x,an)*dep[an];
	ret+=ask1(an,lc)*(dep[an])-ask2(an,lc);
	ret+=ask2(lc,y)+ask1(lc,y)*(dep[an]-2*dep[lc]);
	ret-=(dep[an]-dep[lc])*a[lc];
	return ret;
}
void change(int x,int nu)
{
	chan(t1,dfn[x],nu-a[x]);chan(t1,en[x]+1,a[x]-nu);
	chan(t2,dfn[x],1ll*(nu-a[x])*dep[x]);chan(t2,en[x]+1,1ll*(a[x]-nu)*dep[x]);
	a[x]=nu;
}
int main()
{
	register int i,x,y,z,op;
	n=re_ad();
	for(i=1;i<=n;++i)a[i]=re_ad();
	for(i=1;i<n;++i)
	{
	x=re_ad();y=re_ad();z=re_ad();
	g[x].push_back((node){y,z});g[y].push_back((node){x,z});
	}
	dfs(1);
	for(i=1;i<=n;i++){x=a[i];a[i]=0;change(i,x);}
	q=re_ad();
	while(q--)
	{
	op=re_ad();
	if(op==1){x=re_ad();y=re_ad();pu_t(ask(x,y));
	putchar('\n');}
	else{x=re_ad();y=re_ad();change(x,y);}
	}
}