ACM-ICPC 2018 沈阳赛区网络预赛-G Spare Tire 素因子分解+容斥

A sequence of integer {an} can be expressed as:

Now there are two integers n and m. I'm a pretty girl. I want to find all b1,b2,b3⋯bpthat 1≤bi​≤n and bi is relatively-prime with the integer m. And then calculate:

       

But I have no time to solve this problem because I am going to date my boyfriend soon. So can you help me?

 

Input

Input contains multiple test cases ( about 15000 ). Each case contains two integers nnn and mmm. 1≤n,m≤108.

Output

For each test case, print the answer of my question(after mod 1,000,000,007).

Hint

In the all integers from 1 to 4, 1 and 3 is relatively-prime with the integer 4. So the answer is a1+a3=14.

样例输入

4 4

样例输出

14

---

题意：

       给你一个可以求a数组的公式（如上），现在给你两个数字n和m，先找出1-n之中与m互质的数b1,b2,b3..bp,然后计算。

 

做法：

       其实公式是好求的，ai=i*(i+1)。就是在互质的处理上容易出问题。这就用到了容斥。

       我们要求的是1-n中与m不互质的数xi，设

                                  

       如果xi与p不互质，那么必有公共的素因子，

       如果枚举到一个公共素因子p1的情况（也就是奇数个），那么[1,n]中有n/p1个以p1为因数的数，（如枚举到素因子2时，假设n为8，那么就会有8/2=4个有公共素因子的数）。所以要将其删除。而                                                                                                                             

      所以可以将这些数直接筛掉（具体见代码，因为存在倍数关系所以要乘上除数）。

      如果枚举到两个公共素因子（也就是偶数个），因为在1个的时候两个都各自被处理过，所以会有重复的部分，所以要加上来。

       为什么这样的暴力枚举可行呢，因为2*3*5*7*11*13*17*19*23=223092870所以素因子的个数最多不会超过9个，所以可以用的复杂度直接计算。

 

---

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod=1e9+7;
ll a[20005],inv6=166666668,inv2=500000004,n,m;
ll Calsum(ll n,ll now){
    n/=now;
    return (n*(n+1)%mod*(2*n+1)%mod*inv6%mod*now%mod*now%mod+n*(n+1)%mod*inv2%mod*now%mod)%mod;
}
int init(ll m){
    int num=0;
    for(int i=2;i*i<=m;i++){
        if(m%i==0){
            a[++num]=i;
            while(m%i==0) m/=i;
        }
    }
    if(m!=1) a[++num]=m;
    return num;
}
int main(){
    while(scanf("%lld%lld",&n,&m)!=EOF){
        int num=init(m);
        ll ans1=Calsum(n,1),ans2=0;
        for(int i=1;i<(1<<num);i++){
            ll tmp=1,mark=0;
            for(int j=1;j<=num;j++){
                if(i&(1<<(j-1))){
                    mark++;
                    tmp*=a[j];
                }
            }
            tmp=Calsum(n,tmp);
            if(mark%2) ans2=(ans2+tmp)%mod;
            else ans2=(ans2-tmp)%mod;
        }
        printf("%lld\n",(ans1-ans2+mod)%mod);
    }
    return 0;
}