ACM-ICPC 2018 沈阳赛区网络预赛 G Spare Tire (素因子分解+容斥)
程序员文章站
2022-03-04 22:37:10
...
.
样例输入复制
4 4
样例输出复制
14
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll MOD=1e9+7,inv2=500000004,inv6=166666668;
ll n,m;
const int maxn=100000;
ll A(ll x) {
x %= MOD;
return ((x * x) + x) % MOD;
}
ll mul(ll a, ll b) {
a %= MOD, b %= MOD;
return (a *b) % MOD;
}
ll kA(ll t,ll k) {
ll res= (t+k*t)%MOD;
res = (res*k)%MOD;
res = (res*inv2)%MOD;
ll h = (mul(mul(k,k+1),2*k+1)*inv6)%MOD;
h = (h*t)%MOD;
h = (h*t)%MOD;
res = (res+h)%MOD;
return res;
}
ll a[maxn];
ll b[maxn];
void factor(ll n,ll &tot) {
ll temp,i,now;
temp=(ll)((double)sqrt(n)+1);
tot=0;
now=n;
for(i=2; i<=temp; ++i)if(now%i==0) {
a[++tot]=i;
b[tot]=0;
while(now%i==0) {
++b[tot];
now/=i;
}
}
if(now!=1) {
a[++tot]=now;
b[tot]=1;
}
}
int main() {
while(scanf("%d%d",&n,&m)!=EOF) {
ll all;
factor(m,all);
ll g=0;
for(int i=1; i<(1<<all); i++) {
int cnt=0,t=1;
for(int j=0; j<all; j++) {
if((i>>j)&1) {
t*=a[j+1];
cnt++;
}
}
if(cnt&1) {
if(n/t==0||t>n) continue;
g = (g+kA(t,n/t))%MOD;
} else {
if(n/t==0||t>n) continue;
g = (g-kA(t,n/t)+MOD)%MOD;
}
}
ll ans=g;
ll other = (mul(n, n + 1)*inv2 + mul(mul(n,n+1),2*n+1)*inv6)%MOD;
ans = (other - ans + MOD)% MOD;
cout<<ans<<endl;
}
return 0;
}