【bzoj3561】DZY Loves Math VI 莫比乌斯反演

题目描述

给定正整数n,m。求

输入

一行两个整数n,m。

输出

一个整数，为答案模1000000007后的值。

样例输入

5 4

样例输出

424

题解

莫比乌斯反演

（为了方便，以下公式默认$n\le m$）

$\ \ \ \ \sum\limits_{i=1}^n\sum\limits_{j=1}^mlcm(i,j)^{gcd(i,j)}\\=\sum\limits_{d=1}^n\sum\limits_{i=1}^n\sum\limits_{j=1}^m[gcd(i,j)=d](\frac{ij}d)^d\\=\sum\limits_{d=1}^n\sum\limits_{i=1}^{\lfloor\frac nd\rfloor}\sum\limits_{j=1}^{\lfloor\frac md\rfloor}[gcd(i,j)=1](ijd)^d\\=\sum\limits_{d=1}^nd^d\sum\limits_{i=1}^{\lfloor\frac nd\rfloor}\sum\limits_{j=1}^{\lfloor\frac md\rfloor}[gcd(i,j)=1](ij)^d\\=\sum\limits_{d=1}^nd^d\sum\limits_{i=1}^{\lfloor\frac nd\rfloor}\sum\limits_{j=1}^{\lfloor\frac md\rfloor}(ij)^d\sum\limits_{p|gcd(i,j)}\mu(p)\\=\sum\limits_{d=1}^nd^d\sum\limits_{p=1}^{\lfloor\frac nd\rfloor}\mu(p)\sum\limits_{i=1}^{\lfloor\frac n{dp}\rfloor}\sum\limits_{j=1}^{\lfloor\frac m{dp}\rfloor}(ijp^2)^d\\=\sum\limits_{d=1}^nd^d\sum\limits_{p=1}^{\lfloor\frac nd\rfloor}p^{2d}\mu(p)\sum\limits_{i=1}^{\lfloor\frac n{dp}\rfloor}\sum\limits_{j=1}^{\lfloor\frac m{dp}\rfloor}i^dj^d\\=\sum\limits_{d=1}^nd^d\sum\limits_{p=1}^{\lfloor\frac nd\rfloor}p^{2d}\mu(p)\sum\limits_{i=1}^{\lfloor\frac n{dp}\rfloor}i^d\sum\limits_{j=1}^{\lfloor\frac m{dp}\rfloor}j^d$

此时暴力即可。。。本题就做完了。。。复杂度竟然是对的。。。

快筛$\mu$，再先枚举$d$，然后处理出$1\sim\lfloor\frac md\rfloor$中每个数的$d$次方及其前缀和。由于$d$是从小到大枚举的，因此可以递推出来。

然后只需要先用快速幂求$d^d$，其余的对于每一个$p$均可$O(1)$解出（其中$p^2d=(p^d)^2$），因此总的时间复杂度是每次快速幂的时间复杂度$O(n\log n)$+调和级数时间复杂度$O(n\ln n)$=$O(n\log n)$。

#include <cstdio>

#include <algorithm>

#define N 500010

#define mod 1000000007

using namespace std;

typedef long long ll;

int mu[N] , np[N] , prime[N] , tot;

ll b[N] , sum[N];

inline ll pow(ll x , int y)

{

	ll ans = 1;

	while(y)

	{

		if(y & 1) ans = ans * x % mod;

		x = x * x % mod , y >>= 1;

	}

	return ans;

}

int main()

{

	int n , m , i , j;

	ll d , ans = 0;

	scanf("%d%d" , &n , &m);

	if(n > m) swap(n , m);

	mu[1] = 1;

	for(i = 2 ; i <= n ; i ++ )

	{

		if(!np[i]) mu[i] = -1 , prime[++tot] = i;

		for(j = 1 ; j <= tot && i * prime[j] <= n ; j ++ )

		{

			np[i * prime[j]] = 1;

			if(i % prime[j] == 0)

			{

				mu[i * prime[j]] = 0;

				break;

			}

			else mu[i * prime[j]] = -mu[i];

		}

	}

	for(i = 1 ; i <= m ; i ++ ) b[i] = 1;

	for(i = 1 ; i <= n ; i ++ )

	{

		for(j = 1 ; j * i <= m ; j ++ )

			b[j] = b[j] * j % mod , sum[j] = (sum[j - 1] + b[j]) % mod;

		d = pow(i , i);

		for(j = 1 ; j * i <= n ; j ++ )

			ans = (ans + mu[j] * b[j] * b[j] % mod * d % mod * sum[n / i / j] % mod * sum[m / i / j] % mod + mod) % mod;

	}

	printf("%lld\n" , ans);

	return 0;

}

巴特西

【bzoj3561】DZY Loves Math VI 莫比乌斯反演

最新文章

热门文章