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_ibi??is relatively-prime with the integer?mm. And then calculate:

\displaystyle \sum_{i=1}^{p}a_{b_i}i=1∑p?abi??

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?nnand?mm. 1≤n,m≤10^8.

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.

樣例輸入復(fù)制

4 4

樣例輸出復(fù)制

14

題目來(lái)源

ACM-ICPC 2018 沈陽(yáng)賽區(qū)網(wǎng)絡(luò)預(yù)賽

題意：

給了a數(shù)組的遞推式，問(wèn)所有滿足1<=i<=n且i與m互素的ai之和

思路：

由a的遞推式可以很容易地得到an=n*(n+1)=n^2+n=?

求所有滿足條件的數(shù)不好求，我們可以用所有的減去不滿足條件的，即與m不互素的數(shù)貢獻(xiàn)的a值

根據(jù)算數(shù)基本定理將m分解，與m不互素的就是至少有其中一個(gè)因子，算所有的所以要容斥

對(duì)于每個(gè)因子積sum，會(huì)形成sum,2*sum,3*sum...[n/sum]*sum這些不互素的數(shù)，

設(shè)k=[n/sum]

則對(duì)答案的貢獻(xiàn)值是

sum*(sum+1)+2sum*(2sum+1)+3sum*(3sum+1)+...ksum*(ksum+1)

=(sum+2sum+3sum+...ksum)+(sum*sum+2sum*2sum+3sum*3sum+...ksum*ksum)

=sum ?* ?k*(1+k)/2 + sum*sum ?* ?k*(k+1)*(2k+1)/6

然后容斥一下就可以了,容斥不知道的可以看我的博客 容斥原理（二進(jìn)制枚舉）

#include<bits/stdc++.h> using namespace std;#define e exp(1) #define pi acos(-1) #define mod 1000000007 #define inf 0x3f3f3f3f #define ll long long #define ull unsigned long long #define mem(a,b) memset(a,b,sizeof(a)) int gcd(int a,int b){return b?gcd(b,a%b):a;}ll n,m,inv2,inv6; int cnt,prime[35]; ll qpow(ll a,ll b) {ll ans=1;while(b){if(b&1)ans=(ans*a)%mod;a=a*a%mod;b>>=1;}return ans%mod; } ll sum(ll x,ll k) {ll ans1=k*(1+k)%mod*inv2%mod*x%mod;ll ans2=k*(1+k)%mod*(2*k+1)%mod*inv6%mod*x%mod*x%mod;return (ans1+ans2)%mod; } ll cal(ll n) {ll ans=0;for(int i=1; i<(1<<cnt); i++){ll lcm=1;int k=0;for(int j=0; j<cnt; j++){if((1<<j)&i){k++;lcm*=prime[j];}}ll res=sum(lcm,n/lcm);if(k&1)ans=(ans+res)%mod;else ans=(ans-res+mod)%mod;}return ans; } int main() {inv2=qpow(2,mod-2);inv6=qpow(6,mod-2);while(~scanf("%lld%lld",&n,&m)){cnt=0;for(int i=2; i*i<=m; i++){if(m%i==0){prime[cnt++]=i;while(m%i==0)m/=i;}}if(m>1)prime[cnt++]=m;ll ans=sum(1,n);ans=(ans-cal(n)+mod)%mod;printf("%lld\n",ans);}return 0; }

?

總結(jié)

以上是生活随笔為你收集整理的ACM-ICPC 2018 沈阳赛区网络预赛 Spare Tire（容斥+公式推）的全部?jī)?nèi)容，希望文章能夠幫你解決所遇到的問(wèn)題。

如果覺(jué)得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。