Description
给出区间[l,r],问[l,r]中能整除7且模pi不等于ai的数的个数
Input
第一行一整数T表示用例组数,每组用例首先输入三个整数n,l,r,之后n行每行两个整数pi,ai
(T<=20,0<=n<=15,0< l< r< 10^18,pi互不相同且pi!=7,p1*…*pn<=10^18,0<=ai< pi<=10^5)
Output
对于每组用例,输出[l,r]中能够整除7且模pi不等于ai的数的个数
Sample Input
2
2 1 100
3 2
5 3
0 1 100
Sample Output
Case #1: 7
Case #2: 14
Solution
令S为[1,N]中能够整除7的数的集合,A为[1,N]中能够整除7且模pi不等于ai的数的集合,Ai为[1,N]中能够整数7且模pi等于ai的数的集合,由容斥定理有
后面若干项Ai的并每一项都可以用一次中国剩余定理解决,注意中间有一处乘法可能会爆longlong所以要用快速乘
Code
#include<cstdio>
#include<iostream>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
#define maxn 25
ll a[maxn],p[maxn],c[maxn],m[maxn];
ll extend_gcd(ll a,ll b,ll &x,ll &y)//扩展欧几里得,求ax+by=gcd(a,b)
{
ll d=a;
if(b!=0)
{
d=extend_gcd(b,a%b,y,x);
y-=(a/b)*x;
}
else
{
x=1;
y=0;
}
return d;
}
ll mul(ll a,ll b,ll p)//快速乘,求a*b(mod p)
{
a=(a%p+p)%p,b=(b%p+p)%p;
ll ans=0;
if(a<b)swap(a,b);
while(b)
{
if(b&1)ans=(ans+a)%p;
a=(a+a)%p;
b>>=1;
}
return ans;
}
ll CRT(int n,ll N)//a数组表示余数,m数组表示模数,n表示方程个数,N表示解的上限
{
ll M=1,Mi,ans=0;
ll x,y,d;
for(int i=0;i<n;i++) M*=m[i];
for(int i=0;i<n;i++)
{
Mi=M/m[i];
d=extend_gcd(Mi,m[i],x,y);
ans=(ans+mul(mul(Mi,x,M),a[i],M))%M;
}
if(ans<0) ans=(ans%M+M)%M;
if(ans>N)return 0;
return (N-ans)/M+1;
}
ll Solve(int n,ll N)
{
if(N==0)return 0;
ll ans=N/7;
for(int i=1;i<(1<<n);i++)
{
int res=0;
for(int j=0;j<n;j++)
if(i&(1<<j))
a[res]=c[j],m[res++]=p[j];
a[res]=0,m[res++]=7;
ll cnt=CRT(res,N);
if(res%2)ans+=cnt;
else ans-=cnt;
}
return ans;
}
int main()
{
int T,n,Case=1;
ll x,y;
scanf("%d",&T);
while(T--)
{
scanf("%d%I64d%I64d",&n,&x,&y);
for(int i=0;i<n;i++)
scanf("%I64d%I64d",&p[i],&c[i]);
ll ans=Solve(n,y)-Solve(n,x-1);
printf("Case #%d: %I64d\n",Case++,ans);
}
return 0;
}