题解 P3768 【简单的数学题】

· · 题解

\sum_{i=1}^n\sum_{j=1}^nijgcd(i,j) \sum_{k=1}^n\sum_{i=1}^n\sum_{j=1}^nij[gcd(i,j)=k] \sum_{k=1}^nk^3\sum_{i=1}^{\lfloor \frac{n}{k} \rfloor}\sum_{j=1}^{}ij[gcd(i,j)=1] \sum_{k=1}^nk^3\sum_{d=1}^{\lfloor \frac{n}{k} \rfloor}\mu(d)~d^2\sum_{i=1}^{\lfloor \frac{n}{kd} \rfloor}\sum_{j=1}^{\lfloor \frac{n}{kd} \rfloor}ij

sum^2(n)=\sum_{i=1}^ni^2 , sum^3(n)=\sum_{i=1}^ni^3

利用小学奥数知识有:

sum^2(n)=\frac{n(n+1)(2n+1)}{6} ~~~ , ~~~ sum^3(n)=\frac{n^2(n+1)^2}{4}

上式为:

\sum_{k=1}^nk^3\sum_{d=1}^{\lfloor \frac{n}{k} \rfloor}\mu(d)~d^2sum^3(\frac{n}{kd})

T=kd , 则原式可继续化简为

\sum_{T=1}^n \sum_{k|T}k^3 \mu(\frac{T}{k})(\frac{T}{k})^2sum^3(\frac{n}{T}) \sum_{T=1}^n sum^3(\frac{n}{T})~T^2 \sum_{k|T}k \mu(\frac{T}{k})

后面是一个卷积形式: id*\mu=\varphi

\sum_{T=1}^n sum^3(\frac{n}{T})~T^2 \varphi(T)

现在只需求得 \varphi \times id \times id 的前缀和就可以整除分块了。

f(n)=\varphi(n) \times n \times n g(n)=n \times n h(n)=\sum_{d|n}f(d)* g(\frac{n}{d}) ~~ ~~~~~~~~~~~~~~~~~~~~~~~ =\sum_{d|n}\varphi(d) \times d \times d * \frac{n}{d} * \frac{n}{d} ~~ ~~~~ =n^2\sum_{d|n}\varphi(d)=n^3 g(1)S(n)=\sum_{i=1}^{n}h(i)-\sum_{d=2}^ng(d)S(\lfloor \frac{n}{d} \rfloor) S(n)=sum^3(n)-\sum_{d=2}^ng(d)S(\lfloor \frac{n}{d} \rfloor)

这道题就完了。

#include <map>
#include <cstdio>
using namespace std;

const int MAXN = 8000000;
long long n;
int p , inv4 , inv6 , k , prime[ MAXN + 5 ] , phi[ MAXN + 5 ] , f[ MAXN + 5 ];
bool vis[ MAXN + 5 ];

void sieve( ) {
    f[ 1 ] = phi[ 1 ] = 1;
    for( int i = 2 ; i <= MAXN ; i ++ ) {
        if( !vis[ i ] ) {
            prime[ ++ k ] = i;
            phi[ i ] = i - 1;
        }
        for( int j = 1 ; j <= k && 1ll * i * prime[ j ] <= MAXN ; j ++ ) {
            vis[ i * prime[ j ] ] = 1;
            if( i % prime[ j ] == 0 ) {
                phi[ i * prime[ j ] ] = phi[ i ] * prime[ j ];
                break;
            }
            else
                phi[ i * prime[ j ] ] = phi[ i ] * ( prime[ j ] - 1 );  
        }
        f[ i ] = ( f[ i - 1 ] + 1ll * i * i % p * phi[ i ] % p ) % p;
    }
}

map< long long , int > Map;
int sum2( long long x ) {
    x %= p;
    return 1ll * x * ( x + 1 ) % p * ( 2 * x + 1 ) % p * inv6 % p;
}
int sum3( long long x ) {
    x %= p;
    return 1ll * x * x % p * ( x + 1 ) % p * ( x + 1 ) % p * inv4 % p; 
} 
int Sumphii( long long n ) {
    if( n <= MAXN ) return f[ n ];
    if( Map[ n ] ) return Map[ n ];

    int Ans = sum3( n );
    for( long long l = 2 , r ; l <= n ; l = r + 1 ) {
        r = n / ( n / l );
        Ans = ( Ans - 1ll * ( sum2( r ) - sum2( l - 1 ) ) % p * Sumphii( n / l ) % p ) % p; 
    }
    return Map[ n ] = ( Ans + p ) % p;
}
int solve( long long n ) {
    int Ans = 0;
    for( long long l = 1 , r ; l <= n ; l = r + 1 ) {
        r = n / ( n / l );
        Ans = ( Ans + 1ll * ( Sumphii( r ) - Sumphii( l - 1 ) ) * sum3( n / l ) % p ) % p;
    }
    return ( Ans + p ) % p;
}

int Quick_pow( int x , int po ) {
    int Ans = 1;
    while( po ) {
        if( po & 1 ) Ans = 1ll * Ans * x % p;
        x = 1ll * x * x % p;
        po /= 2;
    }
    return Ans;
}
int Inv( int x ) {
    return Quick_pow( x , p - 2 );
}

signed main( ) {
    scanf("%d %lld",&p,&n);
    sieve( ); 
    inv4 = Inv( 4 ) , inv6 = Inv( 6 );
    printf("%d", solve( n ) );
    return 0;
}