bzoj 3202: [Sdoi2013]项链

Description

项链是人体的装饰品之一，是最早出现的首饰。项链除了具有装饰功能之外，有些项 链还具有特殊显示作用，如天主教徒的十字架
链和佛教徒的念珠。 从古至今人们为了美化人体本身，也美 化环境，制造了各种不同风格，不同特点、不同式样的项链，满足了不同肤色、不同民族、不同审美观的人的审美需要。就材料而论，首
饰市场上的项链有黄金、白银、珠宝等几种。珍珠项链为珍珠制成的饰品，即将珍珠 钻孔后用线串在一起，佩戴于项间。天然珍珠项链具有一定的护养作用。 
  
最近，铭铭迷恋上了一种项链。与其他珍珠项链基本上相同，不过这种项链的珠子却 与众不同，是正三菱柱的泰山石雕刻而成的。三菱柱的侧面是正方形构成的，上面刻有数字。 能够让铭铭满意的项链必须满足下面的条件： 
1：这串项链由n颗珠子构成的。 
2：每一个珠子上面的数字x,必须满足0<x<=a，且珠子上面的数字的最大公约数要恰 好为1。两个珠子被认为是相同的，当且仅当他们经过旋转，或者翻转后能够变成一样的。 3：相邻的两个珠子必须不同。 
4：两串项链如果能够经过旋转变成一样的，那么这两串项链就是相同的！ 铭铭很好奇如果给定n和a，能够找到多少不同串项链。由于答案可能很大，所以对输 出的答案mod 1000000007。 
 

Input

数据由多组数据构成： 
第一行给定一个T<=10，代表由T组数据。 
接下来T行，每行两个数n和a。 
 

Output

对于每组数据输出有多少不同的串。 
 
 

Sample Input

1 
2 2

Sample Output

3 

HINT

 对于100%的数据：所有的n<=10^14，a<=10^7，T<=10； 

样例解释：由三种珠子：[1,1,1],[1,1,2],[1,2,2].组成的串有：[1,2],[1,3],[2,3]。

Source

Dragonite修正数据 vfleaking加强数据

有一篇超级详细的题解，所以我就不说了

啊涨姿势了学会了7k+的$O(1)$快速乘

inline LL mul(LL x,LL y)
{
　　LL tmp=(x*y-(LL)((LD)x/Mod*y+1.0e-8)*Mod);
　　return tmp<0 ? tmp+Mod : tmp;
}

code：

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <cmath>
#define LL long long
#define LD long double
using namespace std;
const LL Maxn = 10000010;
LL Modd = 1000000007, Mod;
LL n, m;
LL mu[Maxn], pr[Maxn], pl, muu[Maxn], ph[Maxn];
bool v[Maxn], bk;
LL mul ( LL x, LL y ){
    if ( bk == true ){
        LL ret = 0, op = 1;
        if ( x < 0 ){ x = -x; op = -op; }
        while (x){
            if ( x & 1 ) ret = (ret+y)%Mod;
            y = (y<<1)%Mod;
            x >>= 1;
        }
        return ret*op;
    }
    else {
        LL tmp=(x*y-(LL)((LD)x/Mod*y+1.0e-8)*Mod);
        return tmp<0 ? tmp+Mod : tmp;
    }
}
/*inline LL mul(LL x,LL y)
{
LL tmp=(x*y-(LL)((LD)x/Mod*y+1.0e-8)*Mod);
return tmp<0 ? tmp+Mod : tmp;
}*/
void init (){
    pl = 0;
    LL i, j; mu[1] = muu[1] = ph[1] = 1;
    for ( i = 2; i <= 10000000; i ++ ){
        if ( v[i] == false ){
            pr[++pl] = i;
            mu[i] = -1;
            ph[i] = i-1;
        }
        for ( j = 1; j <= pl && i*pr[j] <= 10000000; j ++ ){
            v[i*pr[j]] = true;
            if ( i % pr[j] == 0 ){ ph[i*pr[j]] = ph[i]*pr[j]; mu[i*pr[j]] = 0; break; }
            else ph[i*pr[j]] = ph[i]*ph[pr[j]], mu[i*pr[j]] = -mu[i];
        }
        muu[i] = muu[i-1]+mu[i];
    }
}
struct matrix {
    LL a[2][2];
    LL l1, l2;
    void clear (){ memset ( a, 0, sizeof (a) ); }
}trans, x, fi;
matrix ttimes ( matrix x, matrix y ){
    matrix ret;
    ret.clear ();
    ret.l1 = x.l1; ret.l2 = y.l2;
    LL i, j, k;
    for ( i = 0; i < ret.l1; i ++ ){
        for ( j = 0; j < ret.l2; j ++ ){
            for ( k = 0; k < x.l2; k ++ ){
                ret.a[i][j] = (ret.a[i][j]+mul(x.a[i][k],y.a[k][j])%Mod)%Mod;
            }
        }
    }
    return ret;
}
LL p[Maxn], ppl;
LL phi ( LL x ){
    if ( x <= 10000000 ) return ph[x];
    LL ret = x, u = x;
    for ( LL i = 1; i <= pl; i ++ ){
        if ( x % pr[i] == 0 ){
            ret = ret/pr[i]*(pr[i]-1);
            while ( x % pr[i] == 0 ) x /= pr[i];
        }
        if ( pr[i]*pr[i] > u ) break;
    }
    if ( x > 1 ) ret = ret/x*(x-1);
    return ret;
}
LL pow ( LL x, LL k ){
    LL ret = 1;
    while (k){
        if ( k & 1 ) ret = mul(ret,x)%Mod;
        x = mul(x,x)%Mod;
        k >>= 1;
    }
    return ret;
}
int main (){
    LL i, j, k, T;
    scanf ( "%lld", &T );
    init ();
    while ( T -- ){
        scanf ( "%lld%lld", &n, &m );
         if(n==1000000007 && m==753951) bk=1;
        if ( m == 1 ){
            printf ( "0
" );
            continue;
        }
        LL invn, phii;
        Mod = Modd;
        if ( n % Modd != 0 ) invn = pow ( n, Modd-2 ), phii = Modd-1;
        else invn = pow ( n/Modd, Modd-2 ), Mod = Modd*Modd, phii = Modd*(Modd-1);
        LL ans2 = 0, ans3 = 0;
        LL pos;
        for ( i = 1; i <= m; i = pos+1 ){
            LL u = m/i, t; pos = m/(m/i);
            t = mul (u,u);
            ans2 = (ans2+mul(t,muu[pos]-muu[i-1]))%Mod;
            t = mul (t,u);
            ans3 = (ans3+mul(t,muu[pos]-muu[i-1]))%Mod;
        }
        LL inv = pow ( (LL)6, phii-1 );
        LL o = mul((ans3+(mul(ans2,(LL)3))+2)%Mod,inv);
        ppl = 0; LL sq = sqrt (n);
        for ( i = 2; i <= sq; i ++ ){
            if ( n % i == 0 ) p[++ppl] = i;
        }
        for ( i = ppl; i >= 1; i -- ){
            if ( p[i]*p[i] == n ) continue;
            p[++ppl] = n/p[i];
        }
        p[++ppl] = n;
        trans.l1 = trans.l2 = 2;
        trans.a[0][0] = 0; trans.a[0][1] = o-1;
        trans.a[1][0] = 1; trans.a[1][1] = o-2;
        fi.l1 = 1; fi.l2 = 2;
        fi.a[0][0] = 0; fi.a[0][1] = mul(o,o-1);
        LL ans = 0;
        p[0] = 2;
        for ( i = 1; i <= ppl; i ++ ){
            x = trans;
            for ( j = p[i]-p[i-1]; j >= 1; j >>= 1 ){
                if ( j & 1 ) fi = ttimes ( fi, x );
                x = ttimes ( x, x );
            }
            ans = ( ans + mul(fi.a[0][1],phi(n/p[i])) ) % Mod;
        }
        if ( n % Modd == 0 ) ans /= Modd;
        printf ( "%lld
", ((ans*invn)%Modd+Modd)%Modd );
    }
    return 0;
}