题目链接
对于有向图，首先考虑\(bfs\)拓扑，然后就剩下一些环。由于题目中每个点只有一条出边，故剩下的是一些简单环，所以整个图就可以化为一些有向基环树（森林）。

然后考虑计算答案，对于环上的点，考虑该基环树内每个点都可以贡献环上的一个区间，于是直接差分后做前缀和即可。对于树上的点，答案即为子树内与自身距离不超过\(k\)的点的个数。比较直观的做法是线段树合并（本来以为可以长剖，然而要算前缀和），但感觉空间开不下，VP时糊了一个按照\(dfs\)序和深度两个偏序用主席树的做法，写得乱七八糟的之后才调过。。。（其实时空复杂度都要\(O(nlogn)\)，依然很不优美）

对于树上的点，我们不考虑直接计算其子树内的答案（限制相对比较复杂）；而是仍然像计算环上的点那样，考虑每个点的贡献。那么会发现每个点贡献到往上的长度为\(k\)的链上，这就比直接计算子树答案要简单得多。

这时候可以直接树剖，但利用树上差分会更简洁。对链顶的父亲\(-1\)，链底\(+1\)，然后查询每个点的子树和即为答案，用dfs序离线查询区间和即可，时空复杂度都是\(O(n)\)的。

所以其实按照环上的点的思路计算树上的点就很简单了，原本还想了各种数据结构，还是VP的时候降智了吧。

一些错误点
1.有向图区分树边和环边只需要\(bfs\)拓扑，不需要\(tarjan\)（即使写了\(tarjan\)也要记得访问完一个点要退栈！）
2.注意是有向图，环也是有向的！

带log做法

#include
using namespace std;
const int N=5e5+5;
int n,k,fa[N],du[N];
bool is[N];
int hd[N],to[N<<1],nx[N<<1],tt_edge;
void add(int u,int v){
	nx[++tt_edge]=hd[u];
	to[hd[u]=tt_edge]=v;
	du[u]++;
}
queueq;
void bfs(){
	for(int i=1;i<=n;i++) if(!du[i]) is[i]=1,q.push(i);
	while(!q.empty()){
		int u=q.front();
		q.pop();
		int v=fa[u];
		du[v]--;
		if(!du[v]) is[v]=1,q.push(v);
	}
}
int dn[N],sz[N],dp[N],rt[N],lc[N*20],rc[N*20],su[N*20],tt,ct,Rt;
//consider carefully about the apporopriate size__(RE,MLE)
#define mid ((l+r)>>1)
int cnt(int u,int l,int r,int x){
	if(x>=r) return su[u];
	if(xmid) sm+=cnt(rc[u],mid+1,r,x);
	return sm;
}
int bld(int v,int l,int r,int x){
	int u=++ct;
	lc[u]=lc[v];
	rc[u]=rc[v];
	su[u]=su[v]+1;
	if(l==r) return u;
	if(x<=mid) lc[u]=bld(lc[v],l,mid,x);
	else rc[u]=bld(rc[v],mid+1,r,x);
	return u;
}
int ans[N];
void dfs(int u){
	dn[u]=++tt; sz[u]=1;
	rt[tt]=bld(rt[tt-1],0,n,dp[u]);
	for(int e=hd[u];e;e=nx[e]){
		int v=to[e];
		if(!is[v]) continue;
		dp[v]=dp[u]+1;
		dfs(v);
		sz[u]+=sz[v];
	}
	if(u!=Rt) ans[u]+=cnt(rt[dn[u]+sz[u]-1],0,n,min(dp[u]+k,n))-cnt(rt[dn[u]-1],0,n,min(dp[u]+k,n));
}
int p[N],cu[N];
int main()
{
	cin>>n>>k;
	for(int i=1;i<=n;i++) scanf("%d",&fa[i]),add(fa[i],i);
	bfs();
	/*
	!!!You can find circles in a directed graph by just bfs!!!
	!!!Even if you use tarjan,the stack should just store the nodes from the root to u,
		whilch means you should pop u when you finish visiting it!!!
	*/
	for(int i=1;i<=n;i++) if(!is[i]){
		int m=0;
		p[++m]=i;
		for(int j=fa[i];j!=i;j=fa[j]) p[++m]=j;
		for(int j=1;j<=m;j++){
			int u=p[j];
			tt=0; Rt=u;
			dfs(u);
			for(int i2=0;i2<=min(k,sz[u]);i2++){
				int x;
				x=cnt(rt[sz[u]],0,n,i2)-cnt(rt[sz[u]],0,n,i2-1);
				int t=k-i2;
				//remember that the circle is directed!!!
				if(t>=m-1){
					cu[1]+=x;
					continue;
				}
				if(j+t<=m) cu[j]+=x,cu[j+t+1]-=x;
				else cu[j]+=x,cu[1]+=x,cu[j+t-m+1]-=x;
			}
		}
		for(int j=1;j<=m;j++) cu[j]+=cu[j-1],ans[p[j]]+=cu[j],is[p[j]]=1;
		for(int j=0;j<=m+1;j++) cu[j]=0;
	}
	for(int i=1;i<=n;i++) printf("%d\n",ans[i]);
	return 0;
}

线性做法

#include
using namespace std;
const int N=5e5+5;
int n,k,fa[N],du[N];
bool is[N];
int hd[N],to[N<<1],nx[N<<1],edge;
void add(int u,int v){
	nx[++edge]=hd[u];
	to[hd[u]=edge]=v;
	du[u]++;
}
void bfs(){
	queueq;
	for(int i=1;i<=n;i++) if(!du[i]) is[i]=1,q.push(i);
	while(!q.empty()){
		int u=q.front();
		q.pop();
		int v=fa[u];
		du[v]--;
		if(!du[v]) is[v]=1,q.push(v);
	}
}
int dn[N],id[N],sz[N],dp[N],su[N],dfn,rt;
int qu[N],tp;
void dfs(int u){
	qu[++tp]=u;
	dn[u]=++dfn;
	id[dfn]=u;
	sz[u]=1;
	su[dn[u]]++;
	if(tp>k+1) su[dn[qu[tp-k-1]]]--;
	for(int e=hd[u];e;e=nx[e]){
		int v=to[e];
		if(!is[v]) continue;
		dp[v]=dp[u]+1;
		dfs(v);
		sz[u]+=sz[v];
	}
	tp--;
}
int ans[N],p[N],c[N];
int main()
{
	cin>>n>>k;
	for(int i=1;i<=n;i++) scanf("%d",&fa[i]),add(fa[i],i);
	bfs();
	for(int i=1;i<=n;i++) if(!is[i]){
		int m=0;
		p[++m]=i;
		for(int j=fa[i];j!=i;j=fa[j]) p[++m]=j;
		for(int j=1;j<=m;j++){
			int u=p[j];
			dfn=0; rt=u;
			dfs(u);
			for(int l=1;l<=sz[u];l++) su[l]+=su[l-1];
			for(int l=2;l<=sz[u];l++){
				int v=id[l];
				ans[v]+=su[l+sz[v]-1]-su[l-1];
			}
			for(int l=1;l<=sz[u];l++) su[l]=0;
			for(int l=1;l<=sz[u];l++){
				int v=id[l],t=k-dp[v];
				if(t<0) continue;
				if(t>=m-1){
					c[1]++;
					continue;
				}
				if(j+t<=m) c[j]++,c[j+t+1]--;
				else c[j]++,c[1]++,c[j+t-m+1]--;
			}
		}
		for(int j=1;j<=m;j++) c[j]+=c[j-1],ans[p[j]]+=c[j],is[p[j]]=1;
		for(int j=0;j<=m+1;j++) c[j]=0;
	}
	for(int i=1;i<=n;i++) printf("%d\n",ans[i]);
	return 0;
}