[2016ACM多校] HDU5768 容斥原理 中国剩余定理

题意

给出数组p和a，是7的倍数且除以任何一个pi余数不为ai的数是幸运数，问[x,y]中有多少幸运数。

思路

因为幸运数n要满足使得任意一个n%pi=ai不成立，鉴于中国剩余定理可以方便的解出{x%pi=ai}这样的模方程组，所以用容斥原理，先求出7的倍数在减去满足一个方程的，+满足两个方程的-。。。就可以得到结果了。当然，在每次解方程组时都要加上x%7=0这个方程。
注意！pi累乘<10……18且ai<10^5，会爆long long，有取模所以要用快速乘法。赛场上因为计数有点整不清楚，在这wa了好几炮。

AC代码 C++

#include <stdio.h>
#include <string.h>
#include <math.h>
#include <algorithm>

using namespace std;

#define MAXN 20
long long p[MAXN];
long long a[MAXN];
long long tp[MAXN] = {7};
long long ta[MAXN];

//类似快速幂方法解决超过long long范围的long long取模乘法
long long mul(long long a, long long b, long long c)
{
    if(a < 0) a+= c;
    if(b < 0) b+= c;
    long long ans = 0;
    while(b)
    {
        if(b&1) ans = (ans+a)%c;
        a = a*2%c;
        b >>= 1;
    }
    return ans;
}

long long extend_gcd(long long a, long long b, long long& x, long long& y)
{
    if(!b)
    {
        x = 1;
        y = 0;
        return a;
    }
    long long r = extend_gcd(b, a%b, y, x);
    y -= x * (a / b);
    return r;
}

long long M;

long long crt(long long n, long long* a, long long* m)
{
    long long ret=0, x, y, tm, i;
    for(i=0, M=1; i<n; i++)
        M *= m[i];
    for(i=0; i<n; i++)
    {
        tm = M / m[i];
        extend_gcd(tm, m[i], x, y);
        ret = (ret + mul(a[i], mul(tm, x%M, M), M))%M;
    }
    return (ret + M) % M;
}

int main()
{
    long long T, t, n, i, j, cnt, ans;
    long long x, y, r;
    scanf("%I64d", &T);
    for(t=1; t<=T; t++)
    {
        scanf("%I64d%I64d%I64d", &n, &x, &y);
        for(i=0; i<n; i++)
            scanf("%I64d%I64d", p+i, a+i);
        for(i=ans=0; i<(1<<n); i++)
        {
            tp[0] = 7, ta[0] = 0;
            for(j=cnt=0; j<n; j++)
            {
                if(i & (1 << j))
                {
                    tp[++cnt] = p[j];
                    ta[cnt] = a[j];
                }
            }
            r = crt(cnt + 1, ta, tp);
            if(x < r && y < r) continue;
            long long st = (x-r)%M == 0 ? (x-r)/M :  ceil((long double)(1.0*(x-r))/M);
            long long ed = (y-r)%M == 0 ? (y-r)/M : floor((long double)(1.0*(y-r))/M);
            long long sum = max(0LL, ed - st + 1);
            ans += (cnt & 1) ? -sum : sum;
        }
        printf("Case #%I64d: %I64d\n", t, ans);
    }
    return 0;
}