数论

  1. 1 【bzoj 2186】[Sdoi2008]沙拉公主的困惑
  2. 2 【COGS 405】[NOIP2009] Hankson的趣味题
  3. 3 【COGS 1265】[NOIP2012] 同余方程
  4. 4 【COGS 630】【bzoj 2749】[NOIP2011] 计算系数
  5. 5 【COGS 792】[HAOI2012] 外星人
  6. 6 【COGS 662】【bzoj 1041】[HAOI2008] 圆上的整点
  7. 7 [【BZOJ 2705】【COGS 931】 [SDOI2012]Longge的问题]
  8. 8 【COGS 2057】[ZLXOI2015]殉国

数论的内容好多啊
看了一天vfk的PPT

1 【bzoj 2186】[Sdoi2008]沙拉公主的困惑

1~N!中与M!互质的数

\(ans=\varphi(m!)\times \frac {n!}{m!} \mod R\)
因为所有小于m!且与m!互质的数加上m!的整数倍都与m!互质,而其他数都不与m!互质。
\(f(m)=\frac {\varphi(m!)}{m!}=\sum_{gcd(pi,m)=1}\frac {pi-1}{pi}\)
\(\therefore ans=f(m)\times n! \mod R\)
欧拉筛法预处理f()函数
tips:
1 . 如果m<R<=n的话ans=0,因为ans为R的倍数
2 . 求f()时不乘R,不然没法求R的逆元,因为gcd(R,R)!=-1
3 . 求i在mod R意义下的逆元:\(ni[i]=i^{R-2}\),或者用exgcd

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
#include<iostream>
#include<cstdio>
using namespace std;
const int N=10000001;
int R,T,n[10001],m[10001],maxn,tot,prime[N],ni[N],fac[N],ans[N];//factorial 阶乘
bool not_p[N];
int power(int a,int b)
{
int re=1;
for(;b;b>>=1,a=(long long)a*a%R)
if(b&1) re=(long long)re*a%R;
return re;
}
void Pre()
{
not_p[1]=1;fac[1]=1;ans[1]=1;
for(int i=2;i<=maxn;i++)
{
fac[i]=(i!=R)?((long long)fac[i-1]*i%R):fac[i-1];
if(!not_p[i])
{
prime[++tot]=i;ni[i]=power(i,R-2);
ans[i]=(long long)ans[i-1]*(i-1)%R;
if(i!=R) ans[i]=(long long)ans[i]*ni[i]%R;
}
else ans[i]=ans[i-1];
for(int j=1,k;j<=tot&&prime[j]*i<=maxn;j++)
{
not_p[prime[j]*i]=1;
if(!(i%prime[j])) break;
}
}
}
int main()
{
scanf("%d%d",&T,&R);
for(int i=1;i<=T;i++) scanf("%d%d",&n[i],&m[i]),maxn=max(maxn,n[i]);
Pre();
for(int i=1;i<=T;i++)
{
if(R>m[i]&&R<=n[i]) {printf("0\n");continue;}
printf("%lld\n",(long long)fac[n[i]]*ans[m[i]]%R);
}
}

2 【COGS 405】[NOIP2009] Hankson的趣味题

已知正整数a0,a1,b0,b1,求正整数x的个数,满足:
gcd(x,a0)=a1;lim(x,b0)=b1;

首先表示暴力枚举也能过。。。(COGS开O2就是快)

正解:分解质因数
设质因数p
x=p^a×...
a0=p^c1×...
a1=p^c2×...
b0=p^c3×...
b1=p^c4×...
则min(a,c1)=c2,max(a,c3)=c4;

若c1==c2,则a>=c1; --①
若c1!=c2,则a=c2; --②
若c3==c4,则a<=c3; --③
若c3!=c4,则a=c4; --④

当①&&③时,需满足c1<=a<=c3,此时ans×=c3-c1+1;(乘法原理)
当①&&④时,需满足a=c4>=c1;
当②&&③时,需满足a=c2<=c3;
当②&&④时,需满足a=c2=c4;
如果不满足,则ans=0。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
#include<cstdio>
using namespace std;
int T,ans,a0,a1,b0,b1;
int gcd(int x,int y)
{
return y?gcd(y,x%y):x;
}
int main()
{
freopen("son.in","r",stdin);
freopen("son.out","w",stdout);
scanf("%d",&T);
while(T--)
{
ans=0;
scanf("%d%d%d%d",&a0,&a1,&b0,&b1);
for(int x=1;x*x<=b1;x++)
if(!(b1%x))
{
if(gcd(x,a0)==a1&&(long long)gcd(x,b0)*b1==(long long)x*b0) ans++;
if(x*x!=b1&&gcd(b1/x,a0)==a1&&gcd(b1/x,b0)==b0/x) ans++;
}
printf("%d\n",ans);
}
}

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
#include<cstdio>
using namespace std;
const int N=50001;
int T,ans,a0,a1,b0,b1,tot,prime[N];
bool not_p[N];
void pre()
{
not_p[1]=1;
for(int i=2;i<N;i++)
{
if(!not_p[i]) prime[++tot]=i;
for(int j=1;j<=tot&&prime[j]*i<N;j++)
{
not_p[prime[j]*i]=1;
if(!(i%prime[j])) break;
}
}
}
bool calc(int x)
{
int c1=0,c2=0,c3=0,c4=0;
while(!(a0%x)) a0/=x,c1++;
while(!(a1%x)) a1/=x,c2++;
while(!(b0%x)) b0/=x,c3++;
while(!(b1%x)) b1/=x,c4++;
if(c1==c2&&c3==c4&&c1<=c3) ans*=c3-c1+1;
else if((c1==c2&&((c3==c4&&c1>c3)||(c3!=c4&&c4<c1)))||
(c1!=c2&&((c3==c4&&c2>c3)||(c3!=c4&&c2!=c4)))) ans=0;
return ans;
}
int main()
{
freopen("son.in","r",stdin);
freopen("son.out","w",stdout);
scanf("%d",&T);pre();
while(T--)
{
ans=1;
scanf("%d%d%d%d",&a0,&a1,&b0,&b1);
for(int i=1,x=prime[i];x*x<=b1;i++,x=prime[i])
if(!(b1%x))
if(!calc(x))
break;
if(b1!=1) calc(b1);
printf("%d\n",ans);
}
}

3 【COGS 1265】[NOIP2012] 同余方程

求关于 x 的同余方程 ax ≡ 1 (mod b)的最小正整数解。

ax+by=1;
解exgcd,求最小的x

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
#include<cstdio>
#include<cmath>
using namespace std;
int a,b,x,y;
void exgcd(int a,int b,int &x,int &y)
{
if(!b) { x=1; y=0; }
else { exgcd(b,a%b,y,x); y-=x*(a/b); }
}
int main()
{
freopen("mod.in","r",stdin);
freopen("mod.out","w",stdout);
scanf("%d%d",&a,&b);
exgcd(a,b,x,y);
x%=b; while(x<=0) x+=b;
printf("%d",x);
}

4 【COGS 630】【bzoj 2749】[NOIP2011] 计算系数

给定一个多项式 \((ax+by)^k\) ,请求出多项式展开后 \(x^ny^m\) 项的系数mod 10007的结果。

根据多项式定理,\((ax+by)^n=\sum_{r=0}^n C_n^r (ax)^{n-r}(by)^r\)
\(\therefore ans=C_k^n\times a^n\times b^m \mod 10007\)
\(C_k^n=\frac{k!}{n!(k-n)!}=\frac{\prod_{i=1}^{i<=n} k-n+i}{n!}\)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
#include<cstdio>
using namespace std;
const int mod=10007;
int a,b,k,n,m;
int power(int x,int y)
{
int re=1;
for(;y;y>>=1,x=x*x%mod)
if(y&1) re=re*x%mod;
return re;
}
int C(int x,int y)
{
int re1=1,re2=1;
if(y-x<x) x=y-x;
for(int i=1;i<=x;i++)
{
re1=re1*(y-x+i)%mod;
re2=re2*i%mod;
}
return re1*power(re2,mod-2)%mod;
}
int main()
{
freopen("factor.in","r",stdin);
freopen("factor.out","w",stdout);
scanf("%d%d%d%d%d",&a,&b,&k,&n,&m);
printf("%d\n",C(n,k)*power(a%mod,n)%mod*power(b%mod,m)%mod);
}

5 【COGS 792】[HAOI2012] 外星人

求最小的x,使\(\varphi^x(N)=1\),即对N做\(\varphi(N)\)操作多少次后为1

由题目给的公式,对于素数x>2,phi(x)会变出好多2,而phi(2) = 1;所以求出每一个质数经过分解会变出多少个2来,2的个数就是ans。
令f[i]表示i分解出了几个2,用筛法:
i为质数,f[i]=f[i-1];
否则,f[i×prime[j]=f[i]+f[prime[j]];
tip:如果一开始n为奇数,则ans需+1,因为变出2需要第一步。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
#include<cstdio>
using namespace std;
const int N=100001;
int prime[N],tot,f[N],T,m,x,y;
long long ans;
bool notp[N],yes;
void pre()
{
notp[1]=1;f[1]=1;
for(int i=2;i<N;i++)
{
if(!notp[i]) prime[++tot]=i,f[i]=f[i-1];
for(int j=1;j<=tot&&prime[j]*i<N;j++)
{
notp[prime[j]*i]=1;
f[prime[j]*i]=f[prime[j]]+f[i];
if(!(i%prime[j])) break;
}
}
}
int main()
{
freopen("alien.in","r",stdin);
freopen("alien.out","w",stdout);
pre();
scanf("%d",&T);
while(T--)
{
scanf("%d",&m);yes=1;ans=0;
while(m--)
{
scanf("%d%d",&x,&y);ans+=(long long)f[x]*y;
if(x==2) yes=0;
}
printf("%lld\n",ans+yes);
}
}

6 【COGS 662】【bzoj 1041】[HAOI2008] 圆上的整点

平面上有一个圆,圆心坐标为(0,0),半径为n.问圆周上有多少个x,y坐标均为整数的点。

法一:见hzwer

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
#include<cstdio>
#include<cmath>
using namespace std;
long long r,ans;
long long gcd(long long x,long long y)
{
return y?gcd(y,x%y):x;
}
int main()
{
freopen("cir.in","r",stdin);
freopen("cir.out","w",stdout);
scanf("%lld",&r);
for(long long d=1;d*d<=2*r;d++)
if(!((2*r)%d))
{
for(long long a=1;a*a<=r/d;a++)
{
long long b=sqrt(2*r/d-a*a);
if(a!=b&&a*a+b*b==2*r/d&&gcd(a*a,b*b)==1)
ans++;
}
if(d*d!=2*r)
for(long long a=1;a*a<=d/2;a++)
{
long long b=sqrt(d-a*a);
if(a!=b&&a*a+b*b==d&&gcd(a*a,b*b)==1)
ans++;
}
}
printf("%lld\n",ans*4+4);
}

法二: 勾股数的性质(wiki
勾股数:\(a^2+b^2=c^2\)
素勾股数:gcd(a,b,c)=1
设m>n,gcd(m,n)=1,m、n中有一个偶数
\(a=m^2-n^2\)
\(b=2mn\)
\(c=m^2+n^2\)
所以枚举半径n的因数i,求c=i的素勾股数,即x<y且\(x^2+y^2=i^2\)且gcd(x,y)=1且x+y为奇数的个数
ans×8+4,因为枚举的x<y所以要×8,还有坐标轴上的点所以+4

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
#include<cstdio>
#include<cmath>
using namespace std;
int n;long long ans;
int gcd(int x,int y)
{
return y?gcd(y,x%y):x;
}
void calc(int z)
{
for(int x=1,y;2*x*x<z;x++)
{
y=sqrt(z-x*x);
if((long long)x*x+y*y==z&&gcd(x,y)==1&&(x+y)&1) ans++;
}
}
int main()
{
freopen("cir.in","r",stdin);
freopen("cir.out","w",stdout);
scanf("%d",&n);
for(int i=1;i*i<=n;i++)
if(!(n%i))
{
calc(i);
if(i*i!=n) calc(n/i);
}
printf("%lld\n",ans*8+4);
}

7 [【BZOJ 2705】【COGS 931】 [SDOI2012]Longge的问题]

\(\sum_{1<=i<=n} gcd(i,n)\)

枚举n的约数k,令s(k)为满足gcd(m,n)=k,(1<=m<=n)的m的个数,则\(ans=\sum_{k\mid n}k\times s(k)\)
因为gcd(m,n)=k,所以\(gcd(\frac mk,\frac nk)=1\),于是\(s(k)=\varphi(,\frac nk)\)
phi可以在根号的时间内求出

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
#include<cstdio>
using namespace std;
long long n,ans;
long long phi(long long x)
{
long long y=x;
for(long long i=2;i*i<=x;i++)
if(!(x%i))
{
y=y/i*(i-1);
while(!(x%i)) x/=i;
}
if(x>1) y=y/x*(x-1);
return y;
}
int main()
{
freopen("gcdsum.in","r",stdin);
freopen("gcdsum.out","w",stdout);
scanf("%lld",&n);
for(long long i=1;i*i<=n;i++)
if(!(n%i))
{
ans+=i*phi(n/i);
if(i*i!=n) ans+=n/i*phi(i);
}
printf("%lld\n",ans);
}

8 【COGS 2057】[ZLXOI2015]殉国

给定a,b,c,求ax+by=c(x,y>=0)的最小最大的x+y,以及解的组数

求exgcd,由\(xx=x\times \frac cd+\frac bd \times t \ge 0,yy=y\times \frac cd-\frac ad \times t\ge 0\)
得到t的区间[l,r]
方法数=r-l+1;
因为方程线性,所以最值在l、r取到
要写long double!!

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
#include<cstdio>
#include<cmath>
#define min(x,y) x<y?x:y
#define max(x,y) x>y?x:y
using namespace std;
long long a,b,c,d,x,y,xx,yy,ans1,ans2,ans;
void exgcd(long long a,long long b,long long &d,long long &x,long long &y)
{
if(!b) d=a,x=1,y=0;
else exgcd(b,a%b,d,y,x),y-=x*(a/b);
}
int main()
{
freopen("BlackHawk.in","r",stdin);
freopen("BlackHawk.out","w",stdout);
scanf("%lld%lld%lld",&a,&b,&c);
exgcd(a,b,d,x,y);
if(!(c%d))
{
xx=ceil((long double)-x/b*c);
yy=floor((long double)y/a*c);
ans=yy-xx+1;
ans1=x*c/d+y*c/d+(b-a)/d*yy;
ans2=x*c/d+y*c/d+(b-a)/d*xx;
}
if(ans<=0) printf("-1 -1\n0");
else printf("%lld %lld\n%lld",min(ans1,ans2),max(ans1,ans2),ans);
}