CodeChef Sereja and LCM(矩阵快速幂)

Sereja and LCM

Problem code: SEALCM

All submissions for this problem are available.

Read problems statements in Mandarin Chinese and Russian.

In this problem Sereja is interested in number of arrays A[1], A[2], ..., A[N] (1 ≤ A[i] ≤ M, A[i] - integer)
such that least common multiple (LCM) of all its elements is divisible by number D.

Please, find sum of answers to the problem with D = L, D = L+1, ..., D = R. As answer could be large, please output it modulo
(10⁹ + 7).

Input

First line of input contain integer T - number of test cases.
For each test case, only line of input contain four integers N, M, L, R.

Output

For each test case, output required sum modulo (10⁹ + 7).

Constraints

1 ≤ T ≤ 10
1 ≤ N ≤ 5*10⁶
1 ≤ L ≤ R ≤ M ≤ 1000

Subtasks

Subtask #1: 1 ≤ N, M ≤ 10 (25 points)
Subtask #2: 1 ≤ N, M ≤ 1000 (25 points)
Subtask #3: original constraints (50 points)

Example

Input:
2
5 5 1 5
5 5 4 5

Output:
12310
4202

给出 N , M , L , R 。

求用 1~M ，组成N个数，它们的LCM能够整除d的序列个数 , d = L ~ R 。

枚举d , 将d分解，d = p1^k1*p2^k2*....*pn^kn .

当一个数与d有相同质因子 pi 且它的ki 大于等于 d的ki时 ,就能够使得LCM起到作用。

M只有1000 ，直接用二进制表示个个质因子的指数项是否大于等于d 的 ki 。

然后弄好转移矩阵，直接快速幂就OK了。

#include <bits/stdc++.h>
using namespace std;
typedef long long LL;
const int mod = 1e9+7;
const int N = 1010;
int fac[20] , tot , all ;

struct Matrix {
    LL v[50][50] ;
    Matrix(){ memset( v,0,sizeof v); }
    Matrix operator * ( const Matrix &a ) const {
        Matrix res ;
        for( int i = 0 ; i < all ; ++i )
            for( int k = 0 ; k < all ; ++k ) {
                if( v[i][k] == 0 ) continue ;
                for( int j = 0 ; j < all ; ++j ) {
                    res.v[i][j] = ( res.v[i][j] + v[i][k] * a.v[k][j] % mod )% mod;
                }
            }
        return res;
    }
};

Matrix q_pow( Matrix a , int b ) {
    Matrix res = a ;
    b--;
    while( b ) {
        if( b&1 ) res = res * a ;
        a = a * a ;
        b >>= 1 ;
    }
    return res ;
}

void Work( int d ) {
    tot = 0 ;
    for( int i = 2 ; i * i <= d ; ++i ) {
        if( d % i == 0 ) {
            fac[tot++] = 1 ;
            while( d % i == 0 ) fac[tot-1] *= i , d /= i ;
        }
    }
    if( d > 1 ) fac[tot++] = d ;
}

int main() {
//    freopen("in.txt","r",stdin);
    int _ ; cin >> _ ;
    while( _-- ) {
        int n , m , l , r ;
        LL ans = 0 ;
        cin >> n >> m >> l >> r ;
        for( int d = l ; d <= r ; ++d ) {
            Work(d);  Matrix a ;
            all = 1<<tot ;
            for( int i = 0 ; i < all ; ++i ) {
                for( int j = 1 ; j <= m ; ++j ) {
                    int st = i ;
                    for( int k = 0 ; k < tot ; ++k ) {
                        if( j % fac[k] == 0 ) st |= (1<<k) ;
                    }
                    a.v[i][st]++;
                }
            }
            a = q_pow(a,n);
            ans = ( ans + a.v[0][all-1] ) % mod ;
        }
        cout << ans << endl ;
    }
}

View Code