/category/math/ zh-CN Thu, 25 May 2017 20:24:07 +0800 Thu, 25 May 2017 20:24:07 +0800 /archives/236/ /archives/236/ Thu, 25 May 2017 20:24:07 +0800 admin ###Description 给一棵树，每个点点权 $$a_i$$，保证 $$a_i$$各不相同，现在随机选两个点 $$u,v$$，求 $$f(u,v)=\phi(a_u\cdot a_v)\cdot dis(u,v)$$ 的期望，$$mod\ 10^9+7$$。 链接：[http://codeforces.com/contest/809/problem/E](http://codeforces.com/contest/809/problem/E "http://codeforces.com/contest/809/problem/E") ###Solution 首先，有一个结论：$$\phi(a\cdot b)=\frac{\phi(a)\cdot\phi(b)\cdot g}{\phi(g)},g=gcd(a,b)$$ 证明：考虑找一个 $$x$$，使得 $$x|b,gcd(\frac{b}{x},a)=1,gcd(x,\frac{b}{x})=1$$ 且 $$a$$ 包含 $$x$$ 中所有质因子。 那么 $$\phi(a\cdot b)=\phi(a)\cdot x\cdot\phi(\frac{b}{x})=\frac{\phi(a)\cdot\phi(b)\cdot x}{\phi(x)}$$ 显然，$$x$$ 是 $$g$$ 的倍数，所以 $$\frac{x}{g}=\frac{\phi(x)}{\phi(g)}$$，所以$$\phi(a\cdot b)=\frac{\phi(a)\cdot\phi(b)\cdot g}{\phi(g)}$$ 有了这个结论之后可以考虑点分，那么设 $$dis(x)$$ 表示 $$x$$ 到当前根的距离。 设 $$h(x)=\frac{x}{\phi(x)}$$，那么$$f(u,v)=\phi(a_u)\cdot\phi(a_v)\cdot h(g)\cdot(dis(u)+dis(v)),g=gcd(a_u,a_v)$$ 考虑当 $$u$$ 确定时，需要知道 $$\sum\phi(a_v)\cdot h(gcd(a_u,a_v))$$ 和 $$\sum\phi(a_v)\cdot h(gcd(a_u,a_v))\cdot dis(v)$$。 设 $$g(x)=h(x)-\sum[y|x,y\neq x]g(y)$$，那么 $$h(gcd(a_u,a_v))=\sum[x|a_u,x|a_v]g(x)$$ 这样就可以直接枚举约数算贡献了。 ###Code ```c++ #include<bits/stdc++.h> typedef long long ll; inline void repr(int&a,int b){if(a<b)a=b;} #define mset(a,b) memset(a,b,sizeof(a)) #define fo0(i,n) for(int i=0,i##end=n;i<i##end;i++) #define fo1(i,n) for(int i=1,i##end=n;i<=i##end;i++) #define fo(i,a,b) for(int i=a,i##end=b;i<=i##end;i++) #define foe(i,x)for(__typeof(x.end())i=x.begin();i!=x.end();++i) const int N=200007,P=1000000007; struct edge { int to;edge*ne; }_e[N*2],*e=_e,*p[N]; inline void add(int a,int b) { *e=(edge){b,p[a]};p[a]=e++; } std::vector<int>v[N]; int n,pr[N],pm,phi[N],inv[N],f[N],w[N],sz[N],rsz,nro,nsz,g[N],h[N],q[N],qe,ans; bool vis[N]; inline void dfs1(int x,int fa) { sz[x]=1; for(edge*i=p[x];i;i=i->ne) if(i->to!=fa&&!vis[i->to]) { dfs1(i->to,x); sz[x]+=sz[i->to]; } } inline void dfs2(int x,int fa) { int t=rsz-sz[x]; for(edge*i=p[x];i;i=i->ne) if(i->to!=fa&&!vis[i->to]) { repr(t,sz[i->to]); dfs2(i->to,x); } if(t<nsz)nsz=t,nro=x; } inline void set(int x,int y,int z) { if(!~g[x])q[qe++]=x,g[x]=h[x]=0; g[x]=(g[x]+(ll)y*f[x])%P; h[x]=(h[x]+(ll)z*f[x])%P; } inline void get(int x,int&y,int&z) { if(~g[x])y=g[x],z=h[x]; else y=z=0; } inline void dfs3(int x,int fa,int dis) { int a=phi[w[x]],b=(ll)a*dis%P; foe(i,v[w[x]])set(*i,a,b); for(edge*i=p[x];i;i=i->ne) if(i->to!=fa&&!vis[i->to]) dfs3(i->to,x,dis+1); } inline void dfs4(int x,int fa,int dis) { int a=phi[w[x]],b=(ll)a*dis%P,A,B; foe(i,v[w[x]]) { get(*i,A,B); ans=(ans+(ll)a*B+(ll)b*A)%P; } for(edge*i=p[x];i;i=i->ne) if(i->to!=fa&&!vis[i->to]) dfs4(i->to,x,dis+1); } inline void work(int x) { dfs1(x,0); rsz=nsz=sz[x],nro=x; dfs2(x,0); vis[x=nro]=1; foe(i,v[w[x]])set(*i,phi[w[x]],0); for(edge*i=p[x];i;i=i->ne) if(!vis[i->to]) { dfs4(i->to,x,1); dfs3(i->to,x,1); } fo0(i,qe)g[q[i]]=h[q[i]]=-1; qe=0; for(edge*i=p[x];i;i=i->ne) if(!vis[i->to])work(i->to); } int main() { fo1(i,N-1)for(int j=i;j<N;j+=i)v[j].push_back(i); fo(i,2,N-1) { if(!vis[i])pr[pm++]=i,phi[i]=i-1; for(int j=0;i*pr[j]<N;j++) { vis[i*pr[j]]=1; if(i%pr[j]==0) { phi[i*pr[j]]=phi[i]*pr[j]; break; } phi[i*pr[j]]=phi[i]*(pr[j]-1); } } phi[1]=1; inv[1]=1; fo(i,2,N-1)inv[i]=(ll)(P-P/i)*inv[P%i]%P; fo1(i,N-1) { f[i]=(ll)i*inv[phi[i]]%P; foe(j,v[i])if(i^*j)(f[i]+=P-f[*j])%=P; } scanf("%d",&n); fo1(i,n)scanf("%d",w+i); fo0(i,n-1) { int x,y; scanf("%d%d",&x,&y); add(x,y); add(y,x); } mset(vis,0); mset(g,0xff); mset(h,0xff); work(1); ans=(ll)ans*inv[n]%P*inv[n-1]%P*2%P; printf("%d",ans); } ``` ]]> 1 /archives/236/#comments /feed/archives/236/ /archives/181/ /archives/181/ Tue, 28 Feb 2017 00:16:00 +0800 admin ###Description ![fa(1).jpg][1] ###Solution $$\sum_i^n\sum_j^m[\mu(gcd(i,j))\neq 0]lcm(i,j)$$ $$=\sum_k^n|\mu(k)|\sum_i^{\frac{n}{k}}\sum_j^{\frac{n}{k}}[gcd(i,j)=1]i\cdot j\cdot k$$ $$=\sum_k^n|\mu(k)|\cdot k\sum_d^{\frac{n}{k}}\mu(d)\cdot d^2\sum_i^{\frac{n}{k\cdot d}}\sum_j^{\frac{m}{k\cdot d}}i\cdot j$$ $$=\sum_k^n\sum_d[d|k]\mu(d)\cdot d^2\cdot|\mu(\frac{k}{d})|\cdot\frac{p}{d}\cdot\frac{\frac{n}{k}\cdot(\frac{n}{k}+1)}{2}\cdot\frac{\frac{m}{k}\cdot(\frac{m}{k}+1)}{2}$$ 最后的式子，前面是积性函数，后面根号分块 ###Code ```c++ #define _GLIBCXX_IOSTREAM #include<bits/stdc++.h> #define N 4000001 int pri[N],pm,f[N],ps[N],pc[N]; bool np[N]; int main() { for(int i=2;i<N;i++) { if(!np[i])pri[pm++]=i,f[i]=1,ps[i]=i,pc[i]=1; int rf; if(pc[i]==1)rf=f[i]*(ps[i]-ps[i]*ps[i]); else if(pc[i]==2)rf=-f[i]*ps[i]*ps[i]*ps[i]; else rf=0; for(int j=0;i*pri[j]<N;j++) { np[i*pri[j]]=1; if(i%pri[j])f[i*pri[j]]=rf,ps[i*pri[j]]=pri[j],pc[i*pri[j]]=1; else { f[i*pri[j]]=f[i],ps[i*pri[j]]=pri[j],pc[i*pri[j]]=pc[i]+1; break; } } f[i]=rf; } f[1]=1; for(int i=1;i<N;i++) f[i]+=f[i-1]; int T; for(scanf("%d",&T);T--;) { int n,m,ans=0; scanf("%d%d",&n,&m); for(int i=1,t,p,q;i<=n&&i<=m;i=t+1) { p=n/i,q=m/i,t=std::min(n/p,m/q); ans+=(f[t]-f[i-1])*p*(p+1)*q*(q+1)>>2; } printf("%d\n",ans&0x3fffffff); } } ``` [1]: /usr/uploads/2017/02/3245901942.jpg [2]: /usr/uploads/2017/02/1279320617.png ]]> 0 /archives/181/#comments /feed/archives/181/ /archives/154/ /archives/154/ Fri, 16 Dec 2016 23:38:00 +0800 admin ###Description 给出正整数N（可能有前导0），请求出N!最右非零的数位的值。 ###Input 第一行一个数T表示数据组数 下接T行每行一个数N表示一组数据 ###Output 对于每组数据，输出一行一个数表示这组数据的答案 ###Sample Input 2 5 4 ###Sample Output 2 4 ###Solution 随便lucas一下就好了，然而懒得写高精，就用了python，于是实力rank last ###Code ```python T=int(raw_input()) for I in range(T): n=int(raw_input());t=n;a=0;b=0;c=1 while t>0:t/=2;a+=t while n>0:c=c*[1,1,2,1,4,4,4,3,4,1][n%10]%5;n/=5;b+=n a=a-b+1 if a>1:a=0 t=3 while b>0: if b&1:c=c*t%5 t=t*t%5 b/=2 print (a*5+c*6)%10 ``` ]]> 0 /archives/154/#comments /feed/archives/154/ /archives/87/ /archives/87/ Tue, 29 Nov 2016 23:24:00 +0800 admin ###Description 外面有一圈N个结点，中心有一个结点与N个结点都相连，总共就是2*N条边，删除N条边，使N+1个点连通，旋转相同视为等价，问有多少种情况。 ![1.jpg][1] ###Input 输入N，M 3<=N<=10^9, 2<=M<=10^9 ###Output 输出方案数 Mod M的结果 ###Sample Input 3 10000 4 10000 4 10 ###Sample Output 6 13 3 ###Solution 考虑burnside引理，假设分成d块，那么每块的方案数是和bzoj1002一样的，可以用公式 f(n+1) = 3\*f(n) - f(n-1) + 2 计算 所以得到这个公式：1/n\*Sum_{d divides n} phi(n/d)*f(d) 于是直接根号n枚举d，然后根号d求phi，矩阵快速幂求f ###Code ```c++ #include<bits/stdc++.h> typedef long long ll; typedef long double ldb; ll mod; inline ll fm(ll x,ll y) { ll tmp=(x*y-(ll)((ldb)x/mod*y+1e-8)*mod); return tmp<0?tmp+mod:tmp; } inline void mul(ll (&a)[3][3],ll (&b)[3][3],ll (&c)[3][3]) { ll t[3][3]={{0,0,0},{0,0,0},{0,0,0}}; for(int i=0;i<3;i++) for(int j=0;j<3;j++) for(int k=0;k<3;k++) t[i][k]+=fm(a[i][j],b[j][k]); for(int i=0;i<3;i++) for(int j=0;j<3;j++) c[i][j]=t[i][j]%mod; } inline ll f(int x) { if(x<=1)return x; x--; ll a[3][3]={{1,0,2},{0,0,mod-1},{0,1,3}},b[3][3]={{1,0,0},{0,1,0},{0,0,1}}; for(int i=1;i<=x;i<<=1) { if(i&x)mul(b,a,b); mul(a,a,a); } ll ret=b[0][2]+b[2][2]; if(ret>=mod)return ret-mod;return ret; } inline int phi(int x) { int t=x; for(int i=2;i*i<=x;i++) if(x%i==0) { t/=i,t*=i-1; x/=i; while(x%i==0)x/=i; } if(x>1)t/=x,t*=x-1; return t; } int main() { int n; while(~scanf("%d%lld",&n,&mod)) { mod*=n; ll a2=0; for(int i=1;i*i<=n;i++) if(n%i==0) { if(i*i==n) { a2+=fm(phi(i),f(i)); } else { a2+=fm(phi(i),f(n/i))+fm(phi(n/i),f(i)); } a2%=mod; } printf("%lld\n",a2/n); } } ``` [1]: /usr/uploads/2016/11/709092212.jpg ]]> 0 /archives/87/#comments /feed/archives/87/ /archives/54/ /archives/54/ Sun, 16 Oct 2016 00:59:00 +0800 admin ###Description 根据一些书上的记载，上帝的一次失败的创世经历是这样的： 第一天， 上帝创造了一个世界的基本元素，称做“元”。 第二天， 上帝创造了一个新的元素，称作“α”。“α”被定义为“元”构成的集合。容易发现，一共有两种不同的“α”。 第三天， 上帝又创造了一个新的元素，称作“β”。“β”被定义为“α”构成的集合。容易发现，一共有四种不同的“β”。 第四天， 上帝创造了新的元素“γ”，“γ”被定义为“β”的集合。显然，一共会有16种不同的“γ”。 如果按照这样下去，上帝创造的第四种元素将会有65536种，第五种元素将会有2^65536种。这将会是一个天文数字。 然而，上帝并没有预料到元素种类数的增长是如此的迅速。他想要让世界的元素丰富起来，因此，日复一日，年复一年，他重复地创造着新的元素…… 然而不久，当上帝创造出最后一种元素“θ”时，他发现这世界的元素实在是太多了，以致于世界的容量不足，无法承受。因此在这一天，上帝毁灭了世界。 至今，上帝仍记得那次失败的创世经历，现在他想问问你，他最后一次创造的元素“θ”一共有多少种？ 上帝觉得这个数字可能过于巨大而无法表示出来，因此你只需要回答这个数对p取模后的值即可。 你可以认为上帝从“α”到“θ”一共创造了10^9次元素，或10^18次，或者干脆∞次。 一句话题意： ![1.png][1] ###Input 接下来T行，每行一个正整数p，代表你需要取模的值 ###Output T行，每行一个正整数，为答案对p取模后的值 ###Sample Input 3 2 3 6 ###Sample Output 0 1 4 ###HINT 对于100%的数据，T<=1000,p<=10^7 ###Solution 首先考虑欧拉定理，如果把这个数设为S，可以发现$$2^S=2^{S\ mod\ \phi(p)}$$ （当p是2的倍数时可以先除去所有2处理） 那么问题就转化为了求$$S\ mod\ \phi(p)$$，同样考虑，最终$$\phi(p)$$会等于1，直接处理即可 由于除了第一次，每次都会除去至少一个2，所以单次时间复杂度为$$O(logp)$$ ###Code ```c++ #include<cstdio> #define pmax 10000000 int prime[1000000],pm,fi[pmax+1],ans[pmax+1]; bool np[pmax+1],ga[pmax+1]; inline int pow(int a,long long t,int mod) { long long p=1,q=a,f=1; for(;f<=t;f<<=1) { if(f&t)p=p*q%mod; q=q*q%mod; } return p; } void exgcd(int a,int b,int &x,int &y) { if(!b) { x=1,y=0; return; } exgcd(b,a%b,y,x); y-=(a/b)*x; } inline int gen(int x,int t) { if(x==1)return 0; int a,b,p; if(x>t)exgcd(x,t,a,b);else exgcd(t,x,b,a); p=(long long)ans[x]*t%(x*t)*b%(x*t); if(p<0)p+=x*t; return p; } int getans(int x) { int t=1; while(~x&1)x>>=1,t<<=1; if(!ga[x])ga[x]=1,ans[x]=pow(2,1000ll*fi[x]+getans(fi[x]),x); return gen(x,t); } int main() { for(int i=2;i<=pmax;i++) { if(!np[i])prime[pm++]=i,fi[i]=i-1; for(int j=0;j<pm&&i*prime[j]<=pmax;j++) { np[i*prime[j]]=1; if(i%prime[j]) { fi[i*prime[j]]=fi[i]*(prime[j]-1); } else { fi[i*prime[j]]=fi[i]*prime[j]; break; } } } ga[1]=1,ans[1]=0; int T,n; scanf("%d",&T); while(T--) { scanf("%d",&n); printf("%d\n",getans(n)); } } ``` [1]: /usr/uploads/2016/10/774614809.png ]]> 0 /archives/54/#comments /feed/archives/54/ /archives/33/ /archives/33/ Sun, 18 Sep 2016 23:41:00 +0800 admin ###Description 求一个给定的圆(x^2+y^2=r^2)，在圆周上有多少个点的坐标是整数。 ###Input 只有一个正整数n,n<=2000 000 000 ###Output 整点个数 ###Sample Input 4 ###Sample Output 4 ###Solution 只考虑x,y均大于0的情况 首先，暴力枚举肯定是O(n)的，T 对于勾股数有一个结论，就是如果$$x^2+y^2=n^2$$，那么$$x=t(a^2-b^2),y=2tab,n=t(a^2+b^2)$$ 那么我们可以枚举t，显然，只需要枚举n的每个质因数是否在t中就可以了（t中不需要有重复质因子） 枚举t后，需要将n/t写成a^2+b^2的形式，枚举a从1到sqrt(n/t)，时间复杂度O(sqrt(n/t)) 注意去掉重复的a,b，用map实现即可 总时间复杂度应该是根号乘log吧。。 <del>这显然不是正解，不过跑的还是挺快的</del> ###Code ```c++ #include<cstdio> #include<algorithm> #include<cmath> #include<map> std::map<int,bool> f; #define pmax 65536 int prime[10000],pm,x[10],y[10],cc,ans=0,n; bool np[pmax+1]; inline int solve(int nn) { int ans=0; for(int i=sqrt(nn);i>0;i--) { int a=nn-i*i,b=sqrt(nn-i*i); if(b*b==a) { int c=2*b*i,d=std::abs(b*b-i*i); if(!(c&&d))continue; if(!(f[c*(n/nn)]||f[d*(n/nn)])) { f[c*(n/nn)]=f[d*(n/nn)]=1; ans+=8; } } } return ans; } void dfs(int id,int nn) { if(id==cc) { ans+=solve(n/nn); return; } dfs(id+1,nn); dfs(id+1,nn*x[id]); } int main() { for(int i=2;i<=pmax;i++) { if(!np[i])prime[pm++]=i; for(int j=0;j<pm&&i*prime[j]<=pmax;j++) { np[i*prime[j]]=1; if(i%prime[j]==0)break; } } int m; scanf("%d",&n); if(n==0) { printf("1\n"); return 0; } m=n; for(int i=0;i<pm;i++) if(n%prime[i]==0) { x[cc]=prime[i]; while(n%prime[i]==0)n/=prime[i],y[cc]++; cc++; } if(n>1)x[cc++]=n; n=m; dfs(0,1); printf("%d\n",ans+4); return 0; } ``` ]]> 0 /archives/33/#comments /feed/archives/33/ /archives/27/ /archives/27/ Wed, 14 Sep 2016 19:14:00 +0800 admin ###Description ![1407.jpg][1] ###Input 第1行为一个整数N(1<=N<=15)，即野人的数目。 第2行到第N+1每行为三个整数Ci, Pi, Li表示每个野人所住的初始洞穴编号，每年走过的洞穴数及寿命值。 (1<=Ci,Pi<=100, 0<=Li<=10^6) ###Output 仅包含一个数M，即最少可能的山洞数。输入数据保证有解，且M不大于10^6。 ###Sample Input 3 1 3 4 2 7 3 3 2 1 ###Sample Output 6 ###Solution 枚举m，每次枚举两个野人，将他们的ci,pi相减，记为c,p，那么如果有一个k在他们的年龄范围内，且c+pk模m等于0，则这个m不能取，k的最小值可以用扩展欧几里得求出 预处理每两个野人的ci,pi之差，li的最小值，可以优化时间 ###Code ```c++ #include<cstdio> int n,c[15],p[15],l[15],c2[105],p2[105],l2[105],n2=0; inline int min(int a,int b) { if(a<b)return a;else return b; } inline int max(int a,int b) { if(a>b)return a;else return b; } int sol,gcd; void exgcd(int a,int b,int &x,int &y) { if(!b) { x=sol/a;y=0; gcd=a; return; } exgcd(b,a%b,y,x); y-=x*(a/b); } inline int safeexgcd(int a,int b,int &x,int &y) { if(a<b)exgcd(b,a,y,x);else exgcd(a,b,x,y); } int main() { int maxc=0; scanf("%d",&n); for(int i=0;i<n;i++)scanf("%d%d%d",c+i,p+i,l+i),maxc=max(c[i],maxc); for(int i=0;i<n;i++) for(int j=0;j<i;j++) { c2[n2]=c[i]-c[j]; p2[n2]=p[i]-p[j]; l2[n2]=min(l[i],l[j]); n2++; } for(int i=maxc;1;i++) { bool flag=true; int x,y; for(int j=0;j<n2;j++) { sol=-c2[j]; safeexgcd(p2[j],i,x,y); if(p2[j]*x+i*y!=sol)continue; int tmp=i/gcd; x%=tmp; if(x<0)x+=tmp; if(x<=l2[j]) { flag=false; break; } } if(flag) { printf("%d\n",i); return 0; } } } ``` 另外随便写的代码就第4了 ![puh.png][2] [1]: /usr/uploads/2016/09/2556788862.jpg [2]: /usr/uploads/2016/09/909201119.png ]]> 0 /archives/27/#comments /feed/archives/27/