ACM-ICPC 2018瀋陽網路賽G題(唯一分解定理、容斥原理)
A sequence of integer \lbrace a_n \rbrace{an} can be expressed as:
\displaystyle a_n = \left\{ \begin{array}{lr} 0, & n=0\\ 2, & n=1\\ \frac{3a_{n-1}-a_{n-2}}{2}+n+1, & n>1 \end{array} \right.an=⎩⎨⎧0,2,23an−1−an−2+n+1,n=0n=1n>1
Now there are two integers nn and mm. I'm a pretty girl. I want to find all b_1,b_2,b_3\cdots b_pb1,b2,b3⋯bp that 1\leq b_i \leq n1≤bi≤n and b_ibiis relatively-prime with the integer mm. And then calculate:
\displaystyle \sum_{i=1}^{p}a_{b_i}i=1∑pabi
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 1500015000 ). Each case contains two integers nn and mm. 1\leq n,m \leq 10^81≤n,m≤108.
Output
For each test case, print the answer of my question(after mod 1,000,000,0071,000,000,007).
Hint
In the all integers from 11 to 44, 11 and 33 is relatively-prime with the integer 44. So the answer is a_1+a_3=14a1+a3=14.
樣例輸入複製
4 4
樣例輸出複製
14
如果要看詳細的解釋可以看這一篇,我覺得這個公式還是挺好推的。
之所以寫這篇部落格是因為要記住一個教訓,在質數打表的時候需要開大一點的空間,比如說1e8的資料,就不能只開1e4的質數陣列,不然一直TLE 或者WA。
下面是程式碼:
#include<iostream> #include<cstdio> #include<cstring> #include<algorithm> #include<vector> #include<stack> #include<queue> #include<cmath> #include<map> #define fori(l,r) for( int i = l ; i <= r ; i++ ) #define forj(l,r) for( int j = l ; j <= r ; j++ ) #define mem(a,val) memset(a,val,sizeof a) #define inf 0x3f3f3f3f #define longinf 0x3f3f3f3f3f3f3f3f using namespace std; typedef long long ll; ll n,m; const int maxn = 1e5+5; bool isprime[maxn]; ll prime[maxn]; ll a[maxn]; int cnt,cc; const ll mod = 1e9+7; ll six; ll ans; ll down; const ll Six = 166666668; const ll Two = 500000004; //乘法逆元用的 ll fastpow( ll base,ll x ) { ll ans = 1; while( x ) { if( x&1 ) ans = (ans*base)%mod; x = x>>1; base = (base*base)%mod; } return ans; } void init() { cnt = 1; mem(isprime,true); for( int i = 4 ; i < maxn ; i+=2 ) isprime[i] = false; for( int i = 3 ; i < maxn ; i+=2 ) if( isprime[i] ) for( int j = i<<1 ; j < maxn ; j+=i ) isprime[j] = false; fori(2,maxn-1) if( isprime[i] ) prime[cnt++] = i; down = fastpow(2,mod-2); six = fastpow(6,mod-2); } void decompose( ll x ) { for( int i = 1 ; prime[i]*prime[i] <= x ; i++ ) { if( x%prime[i] == 0 ) { while( x%prime[i] == 0 ) x = x/prime[i]; a[cc++] = prime[i]; } } if( x > 1 ) a[cc++] = x; } inline ll Mod(ll a,ll b) { return (a%mod)*(b%mod)%mod; } inline ll ff(ll x) { ll k = n/x; return (Mod(x,x)*Mod(Mod(k,k+1),Mod(k+k+1,Six))%mod + Mod(Mod(1+k,k),Mod(x,Two)) )%mod; } void dfs( int start,int times,int goal,ll val ) { if( times == goal ) { if( times&1 ) ans = ( ans+ff(val) )%mod; else ans = ( ans-ff(val)+mod )%mod; return; } for( int i = start ; i < cc ; i++ ) { dfs(i+1,times+1,goal,val*a[i]); } } int main() { init(); while( scanf("%lld %lld",&n,&m) == 2 ) { ans = 0; cc = 1; decompose(m); for( int i = 1 ; i < cc ; i++ ) dfs(1,0,i,1); ans = ( ff(1)-ans+mod )%mod; printf("%lld\n",ans); } return 0; } /* 99999281 99999827 100000000 100000000 9999999 9999999 9999998 9999998 9999997 9999997 9999998 9999997 100000000 9999997 6666666 6666666 6666666 6666667 */