[洛谷 P6156] 简单题 及其加强

题目

点这里看原版题目。
点这里看加强题目。

分析

简单版

不难发现 f(n)=μ2(n)

(这里定义 fk(n)=(f(n))k

然后开始娴熟地推式子:

i=1nj=1nμ2(gcd(i,j))gcd(i,j)(i+j)k=d=1nμ2(d)di=1ndj=1nd(id+jd)k[gcd(i,j)=1]=d=1nμ2(d)dk+1t=1ndμ(t)tki=1ndtj=1ndt(i+j)k

为了方便简洁,我们定义 S(n)=i=1nj=1n(i+j)k ,然后就有:

d=1nμ2(d)dk+1t=1ndμ(t)tki=1ndtj=1ndt(i+j)k=d=1nμ2(d)dk+1t=1ndμ(t)tkS(ndt)

现在,只要我们可以快速地求出 S(n) ,我们就可以两次整除分块来解决这个问题。

首先不难发现 S(n) 有一个相似的定义:

S(n)=i=12nik×min{i1,2ni+1}

注意到 S(n) 取值只有 O(n) 种,因此我们可以对于它的每一种取值 m ,用 O(m) 的时间算出 S(n) 。于是我们就可以用 o(nlnn) 的时间求出所有需要的 S

然后我们就可以愉快地计算了。

加强版

为了更方便地计算,我们需要继续推式子。下面令 T=dt

d=1nμ2(d)dk+1t=1ndμ(t)S(ndt)=T=1nS(nT)d|Tμ2(d)dk+1μ(Td)(Td)k=T=1nS(nT)Tkd|Tμ2(d)μ(Td)d

式子已经化到底了。现在,如果我们可以快速地求出来 f(n)=d|nμ2(d)μ(nd)dS(n) ,我们就可以运用整除分块,用 O(n) 的时间回答一次查询。

首先考虑 S 。不难发现 S 实际上是一个三角形系数乘上自然数幂的和(类似于前三角加后三角)。
简单推一推就可以发现,如果令 g(n)=i=1n(ni+1)×ik+1 ,那么 S(n)=g(2n)2g(n)

现在我们只需要考虑求 f(n) 了。注意到 f(n) 实际上是 μ2(n)nμ(n) 的狄利克雷卷积,所以它是积性函数

积性函数,我们就可以使用欧拉筛来求值。考虑新加入一个质因子 p ,我们只需要考虑 f(pk)

  1. f(p0)=1
  2. f(p1)=μ2(1)μ(p)+pμ2(p)μ(1)=p1
  3. f(p2)=μ2(1)μ(p2)+pμ2(p)μ(p)+p2μ2(p2)μ(1)=p
  4. f(pk)(k>2)=0 ,因为无论怎么分, μ2(d)μ(nd) 中都会有一个的幂次大于 2 。

于是我们就可以在欧拉筛中,加上判断幂次的过程,并求得 f 的值。

注意,本题有那么一点点地卡空间

一些有价值的点:

  1. 求出 f(n) 的方法比较通用。 n 比较小的时候,我们通常可以使用欧拉筛来计算各处的值。(欧拉筛可以被分为增加幂次增加新质因子两个部分,并且进行分类讨论)

代码

简单版:

#include <cstdio>

typedef long long LL;

const int mod = 998244353;
const int MAXN = 1e7 + 5;

template<typename _T>
void read( _T &x )
{
	x = 0;char s = getchar();int f = 1;
	while( s > '9' || s < '0' ){if( s == '-' ) f = -1; s = getchar();}
	while( s >= '0' && s <= '9' ){x = ( x << 3 ) + ( x << 1 ) + ( s - '0' ), s = getchar();}
	x *= f;
}

template<typename _T>
void write( _T x )
{
	if( x < 0 ){ putchar( '-' ); x = ( ~ x ) + 1; }
	if( 9 < x ){ write( x / 10 ); }
	putchar( x % 10 + '0' );
}

template<typename _T>
_T MIN( const _T a, const _T b )
{
	return a < b ? a : b;
}

int su1[MAXN], su2[MAXN];
int pw[MAXN], mu[MAXN];
int prime[MAXN], pn;
bool isPrime[MAXN];

int S[MAXN];
int N; LL K;

int Sub( int x, int v ) { return x < v ? x + mod - v : x - v; }
int Mul( LL x, int v ) { x *= v; if( x >= mod ) x %= mod; return x; }
int Add( int x, int v ) { return x + v >= mod ? x + v - mod : x + v; }

int Qkpow( int base, LL indx )
{
	int ret = 1;
	while( indx )
	{
		if( indx & 1 ) ret = Mul( ret, base );
		base = Mul( base, base ), indx >>= 1;
	}
	return ret;
}

void EulerSieve( const int siz )
{
	su1[1] = su2[1] = pw[1] = 1;
	for( int i = 2 ; i <= siz ; i ++ )
	{
		if( ! isPrime[i] ) mu[prime[++ pn] = i] = mod - 1, pw[i] = Qkpow( i, K );
		for( int j = 1 ; j <= pn && 1ll * i * prime[j] <= siz ; j ++ )
		{
			isPrime[i * prime[j]] = true, pw[i * prime[j]] = Mul( pw[i], pw[prime[j]] );
			if( ! ( i % prime[j] ) ) break; mu[i * prime[j]] = mod - mu[i];
		}
		su1[i] = Add( su1[i - 1], Mul( Mul( mu[i], mu[i] ), Mul( pw[i], i ) ) );
		su2[i] = Add( su2[i - 1], Mul( mu[i], pw[i] ) );
	}
}

int GetSum( const int lim )
{
	int ret = 0;
	for( int k = 1 ; k <= 2 * lim ; k ++ )
		ret = Add( ret, Mul( pw[k], MIN( k - 1, 2 * lim - k + 1 ) ) );
	return ret;
}

void Init()
{
	EulerSieve( N << 1 );
	for( int l = 1, r ; l <= N ; l = r + 1 )
	{
		r = N / ( N / l );
		S[N / l] = GetSum( N / l );
	}
}

int f( const int n )
{
	int ret = 0;
	for( int l = 1, r ; l <= n ; l = r + 1 )
	{
		r = n / ( n / l );
		ret = Add( ret, Mul( Sub( su2[r], su2[l - 1] ), S[n / l] ) );
	}
	return ret;
}

int main()
{
	read( N ), read( K ); 
	Init(); int ans = 0;
 	for( int l = 1, r ; l <= N ; l = r + 1 )
	{
		r = N / ( N / l );
		ans = Add( ans, Mul( Sub( su1[r], su1[l - 1] ), f( N / l ) ) );
	}
	write( ans ), putchar( '\n' );
	return 0;
}

加强版:

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

typedef long long LL;
typedef unsigned int ui;

const int MAXN = 2e7 + 5;

template<typename _T>
void read( _T &x )
{
	x = 0;char s = getchar();int f = 1;
	while( s > '9' || s < '0' ){if( s == '-' ) f = -1; s = getchar();}
	while( s >= '0' && s <= '9' ){x = ( x << 3 ) + ( x << 1 ) + ( s - '0' ), s = getchar();}
	x *= f;
}

template<typename _T>
void write( _T x )
{
	if( x < 0 ){ putchar( '-' ); x = ( ~ x ) + 1; }
	if( 9 < x ){ write( x / 10 ); }
	putchar( x % 10 + '0' );
}

template<typename _T>
_T MIN( const _T a, const _T b )
{
	return a < b ? a : b;
}

ui pw[MAXN], f[MAXN];
int prime[MAXN / 10], pn;
bitset<MAXN> isPrime;

int N; LL K;

ui Qkpow( ui base, int indx )
{
	ui ret = 1;
	while( indx )
	{
		if( indx & 1 ) ret *= base;
		base *= base, indx >>= 1;
	}
	return ret;
}

void EulerSieve( const int siz )
{
	f[1] = pw[1] = 1;
	for( int i = 2 ; i <= siz ; i ++ )
	{
		if( ! isPrime[i] ) pw[prime[++ pn] = i] = Qkpow( i, K ), f[i] = i - 1;
		for( int j = 1 ; j <= pn && 1ll * i * prime[j] <= siz ; j ++ )
		{
			isPrime[i * prime[j]] = true, pw[i * prime[j]] = pw[i] * pw[prime[j]];
			if( ! ( i % prime[j] ) ) 
			{
				int q = i / prime[j];
				if( q % prime[j] ) f[i * prime[j]] = - prime[j] * f[q];
				break; 
			}
			f[i * prime[j]] = f[i] * f[prime[j]];
		}
	}
}

void Init()
{
	EulerSieve( N << 1 );
	for( int i = 1 ; i <= N << 1 ; i ++ ) 
		f[i] = f[i - 1] + pw[i] * f[i], pw[i] = pw[i - 1] + pw[i];
	for( int i = 1 ; i <= N << 1 ; i ++ ) pw[i] = pw[i - 1] + pw[i];
}

int S( const int n ) { return pw[n << 1] - ( pw[n] << 1 ); }

int main()
{
	int T, n;
	read( T ), read( N ), read( K );
	Init();
	while( T -- )
	{
		read( n ); ui ans = 0;
		for( int l = 1, r ; l <= n ; l = r + 1 )
		{
			r = n / ( n / l );
			ans += ( f[r] - f[l - 1] ) * S( n / l );
		}
		write( ans ), putchar( '\n' );
	}
	return 0;
}
posted @   crashed  阅读(156)  评论(0编辑  收藏  举报
编辑推荐:
· 为什么说在企业级应用开发中,后端往往是效率杀手?
· 用 C# 插值字符串处理器写一个 sscanf
· Java 中堆内存和栈内存上的数据分布和特点
· 开发中对象命名的一点思考
· .NET Core内存结构体系(Windows环境)底层原理浅谈
阅读排行:
· 为什么说在企业级应用开发中,后端往往是效率杀手?
· DeepSeek 解答了困扰我五年的技术问题。时代确实变了!
· 本地部署DeepSeek后,没有好看的交互界面怎么行!
· 趁着过年的时候手搓了一个低代码框架
· 推荐一个DeepSeek 大模型的免费 API 项目!兼容OpenAI接口!
点击右上角即可分享
微信分享提示