[2018.6.21集训]走路-分块高维前缀和-Pollard-Rho

题目大意

给一棵树，每个节点有一个权值$val$。
如果两个点$a$和$b$满足$a$为$b$的祖先且$val[b]$为$val[a]$的约数，那么可以从$a$一步跳到$b$。
求从$1$号节点走到各每个节点的路径数。

$n leq 10^5 , val[i] leq 10^{18},$保证对于任意节点$i$，$val[i]$为$val[1]$的约数。

题解

首先有定理:一个数$n$的约数个数大概是$n^{frac{1}{3}}$的级别
然后得到有理有据的一个部分分（$val[i] leq 10^9,40$分）算法:
对$val[1]$分解质因数，求出其所有约数，然后开等同于约数个数个桶，每到一个节点就枚举每个约数看是否为自己的倍数，如果是就取出对应桶中的值加入当前点的$dp$值，最后再把当前点也加入桶中。
（事实上对所有$val$进行一次$unique$就可以不用分解质因数了）

复杂度$O(n*10^6)$，显然不能过。

于是来一发神奇的分块高维前缀和。
具体而言，将每个数分解为若干个质因数的乘积，那么就可以将这个数看成一个向量，每一维的值对应其某一个质因数的数量。
将向量从中间划开，分成尽量相等的两组维度，然后对于其中一组维度做前缀和。

也就是说在此题中，假如是一个二维向量$(i,j)$，对应质因数$2^i*3j$。
此时，将向量划成$i$和$j$两组维度，那么将向量$i,j$的贡献$v$加入的过程如下:

for(int k=i;k<=max_i;k++)
    sum[k][j]+=v

查询向量$(i,j)$的前缀和则用:

for(int k=j;k<=max_j;k++)
    ans+=sum[i][k]

此时查询复杂度和修改复杂度均为$O(sqrt(n))$，$n$为所有维度的$max_i$之积，在这里为$10^6$。

于是复杂度降为$O(n*10^3)=O(108)$，玄学卡常一波即可通过~
由于要分解$10^{18}$级别的大数，还需要来一发Pollard-Rho......

代码:

#include<map>
#include<vector>
#include<cstdio>
#include<algorithm>
using namespace std;

typedef long long ll;
const int N=1e5+9;
const ll md=1e9+7;

inline ll read()
{
	ll x=0,f=1;char ch=getchar();
	while(ch<'0' || '9'<ch){if(ch=='-')f=-1;ch=getchar();}
	while('0'<=ch && ch<='9')x=x*10+(ch^48),ch=getchar();
	return x*f;
}

inline void write(ll x)
{
	if(x>=10)write(x/10);
	putchar('0'+x%10);
}

int n;
int to[N<<1],nxt[N<<1],beg[N],tot;
int fa[N],id[N],ed[N],seg[N],dfn;
ll val[N],f[N];

inline void add(int u,int v)
{
	to[++tot]=v;
	nxt[tot]=beg[u];
	beg[u]=tot;
}

namespace di
{
	ll stk[100],cnt[100],top,c;
	vector<ll> vec,vecs;
	vector<ll> vf;
	map<ll,int> ha;

	inline ll gcd(ll a,ll b){return !b?a:gcd(b,a%b);}

	inline ll mul(ll a,ll b,ll p)
	{
		return ((a*b-(ll)((long double)a*b/p)*p)%p+p)%p;
	}

	inline ll qpow(ll a,ll b,ll p)
	{
		ll ret=1;
		while(b)
		{
			if(b&1)ret=mul(ret,a,p);
			a=mul(a,a,p);b>>=1;
		}
		return ret;
	}

	inline bool check(ll a,ll n)
	{
		ll x=n-1,t=0;
		for(;!(x&1);x>>=1)t++;
		x=qpow(a,x,n);
		if(x==1 || x==n-1)return 1;
		x=mul(x,x,n);
		for(int i=1;i<=t;i++,x=mul(x,x,n))
		{
			if(x==n-1)return 1;
			if(x==1)return 0;
		}
		return 0;
	}

	inline bool miller_rabin(ll x)
	{
		if(x==1 || x==2)return 1;
		for(int i=1;i<=15;i++)
			if(!check(rand()%(x-1)+1,x))
				return 0;
		return 1;
	}

	inline void find(ll x);

	inline bool pollard_rho(ll x)
	{
		int cnt=0,k=0;
		ll x1=rand()%x,x2=x1,d;
		while(1)
		{
			x1=(mul(x1,x1,x)+c)%x;
			d=gcd(abs(x1-x2),x);
			if(d!=1 && d!=x)
				return find(x/d),find(d),1;
			if(x1==x2)return 0;
			if((++cnt)==(1<<k))
				k++,x2=x1;
		}
	}

	inline void find(ll x)
	{
		if(miller_rabin(x))
		{
			stk[++top]=x;
			cnt[top]=0;
			return;
		}
		while(!pollard_rho(x))
			c=rand();
	}

	inline void work(ll x)
	{
		top=0;find(x);int e=top;
		sort(stk+1,stk+top+1);
		top=0;stk[0]=-1;
		for(int i=1;i<=e;i++)
		{
			if(stk[i]!=stk[top])
			{
				stk[++top]=stk[i];
				cnt[top]=0;
			}
			cnt[top]++;
		}
	}
	
	inline void dfss(int x,ll v)
	{
		if(x==top+1)
		{
			vec.push_back(v);
			vf.push_back(0);
			ha[v]=vec.size()-1;
			return;
		}

		ll mul=1;
		for(int i=0;i<=cnt[x];i++,mul*=stk[x])
			dfss(x+1,v*mul);
	}
}

using namespace di;

inline void chk(ll &x){if(x>=md)x-=md;}

namespace sbds
{
	const int M=5009;
	int mid;
	ll ds[M][M];

	inline void init(ll x)
	{
		work(x);
		mid=top>>1;
	}

	inline void modifys(ll x,int nid,int id,ll v,int p)
	{
		if(p==mid+1)
		{
			chk(ds[nid][id]+=v);
			return;
		}
		int cnts=0;
		while(x/stk[p]*stk[p]==x)
			x/=stk[p],cnts++;
		for(int i=cnts;i>=0;i--)
			modifys(x,nid*(cnt[p]+1)+i,id,v,p+1);
	}
		
	inline void modify(ll x,ll v)
	{
		int id=0;
		for(int i=mid+1;i<=top;i++)
		{
			int cnts=0;
			while(x/stk[i]*stk[i]==x)
				x/=stk[i],cnts++;
			id=id*(cnt[i]+1)+cnts;
		}
		modifys(x,0,id,v,1);
	}

	inline ll querys(ll x,int nid,int id,int p)
	{
		if(p==top+1)
			return ds[id][nid];
		int cnts=0;ll ret=0;
		while(x/stk[p]*stk[p]==x)
			x/=stk[p],cnts++;
		for(int i=cnts;i<=cnt[p];i++)
			chk(ret+=querys(x,nid*(cnt[p]+1)+i,id,p+1));
		return ret;
	}

	inline ll query(ll x)
	{
		int id=0;
		for(int i=1;i<=mid;i++)
		{
			int cnts=0;
			while(x/stk[i]*stk[i]==x)
				x/=stk[i],cnts++;
			id=id*(cnt[i]+1)+cnts;
		}
		return querys(x,0,id,mid+1);
	}
}

inline void dfs(int u)
{
	if(u!=1)
		f[u]=sbds::query(val[u]);
	sbds::modify(val[u],f[u]);
	for(int i=beg[u];i;i=nxt[i])
		if(to[i]!=fa[u])
			fa[to[i]]=u,dfs(to[i]);
	sbds::modify(val[u],md-f[u]);
}

int main()
{
	srand(514);n=read();
	for(int i=2,u,v;i<=n;i++)
	{
		u=read();v=read();
		add(u,v);add(v,u);
	}
	for(int i=1;i<=n;i++)
		val[i]=read();

	sbds::init(val[1]);
	dfs(f[1]=1);
	for(int i=1;i<=n;i++)
		write(f[i]),putchar('
');
	return 0;
}