SDOI2014 数表

wjh15101051's space

卧薪尝胆，厚积薄发。

SDOI2014 数表

Date: Wed Jul 04 11:17:15 CST 2018

In Category: NoCategory

Description：

求 $\begin{align}\sum_{i=1}^n\sum_{j=1}^m\sum_{d|i,d|j}d[\sum_{d|i,d|j}d\le a]\end{align}$

多组询问，每次给出 $n,m,a$ 。

Solution：

限制不好做，还是不考虑 $a$ 。

要求的就是： $$ \begin{align}\sum_{i=1}^n\sum_{j=1}^m\sigma(gcd(i,j))\end{align} $$ 下一步不好化，设 $gcd(i,j)=d$

则原式： $$ \begin{align}=\sum_d\sum_{i=1}^n\sum_{j=1}^m\sigma(d)[gcd(i,j)==d]\end{align} $$ 后面的式子可以反演：

$$ \begin{align}=\sum_d\sigma(d)\sum_{i=1}^{\lfloor \frac n d\rfloor}\sum_{j=1}^{\lfloor \frac m d\rfloor}\sum_{k|i,k|j}\mu(k)\end{align} $$ 交换求和号，得：

$$ \begin{align}=\sum_d\sigma(d)\sum_k\mu(k)\sum_{i=1,k|i}^{\lfloor \frac n d\rfloor}\sum_{j=1,k|j}^{\lfloor \frac m d\rfloor}1=\sum_d\sigma(d)\sum_k\mu(k)\lfloor \frac n {dk}\rfloor\lfloor \frac m {dk}\rfloor\end{align} $$ 按惯例设一个 $T=d\times k$ 并在外层枚举。

原式： $$ \begin{align}=\sum_{T=1}^n\lfloor\frac n T\rfloor\lfloor\frac m T\rfloor\sum_{d|T}\sigma(d)\mu(\frac T d)\end{align} $$ 前面是一个除法分块，后面是一个狄利克雷卷积，设其为 $g$ ，则： $$ \begin{align}=\sum_{T=1}^n\lfloor\frac n T\rfloor\lfloor\frac m T\rfloor g(T)\end{align} $$ 由于 $\sigma$ 和 $\mu$ 都积性，所以 $g$ 也积性。 $\mu$ 线性筛， $\sigma$ 枚举每个数把他的倍数加上这个数，由调和级数得时间复杂度 $O(log n)$

考虑 $a$ 的影响，由于除法分块要求 $g$ 得前缀和，即使用数据结构维护在线也非常不好做。

把所有询问离线，按 $a$ 排序，将所有 $\sigma$ 值求出来排序。

用树状数组维护 $g$ ，对于每个询问，把 $\sigma(x)$ 比 $a$ 小的数乘上 $\mu$ 插到 $x$ 的倍数的位置，这样树状数组中保存的就是合法的 $g$ 值可以直接求前缀和了。

Code：


#include<algorithm>

#include<iostream>

#include<cstdlib>

#include<cstdio>

#include<cmath>

#include<cstring>

using namespace std;

int T;

#define MAXT 20010

#define MOD (1ll << 31)

typedef long long ll;

struct query

{

	int n,m,a,id;

	ll ans;

}q[MAXT];

bool cmpq(query x,query y){return x.a < y.a;}

bool cmpid(query x,query y){return x.id < y.id;}

#define MAXN 100010

pair<ll,int> s[MAXN];

int mu[MAXN];

bool isprime[MAXN];

int prime[MAXN],tot = 0;

void init()

{

	for(int i = 2;i < MAXN;++i)isprime[i] = true;

	mu[1] = 1;

	for(int i = 2;i < MAXN;++i)

	{

		if(isprime[i])

		{

			prime[++tot] = i;

			mu[i] = -1;

		}

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

		{

			isprime[i * prime[j]] = false;

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

			{

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

				break;

			}

			else

			{

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

			}

		}

	}

	for(int i = 1;i < MAXN;++i)

	{

		s[i].second = i;

		for(int j = 1;j * i < MAXN;++j)

		{

			s[i * j].first += i;

		}

	}

	sort(s + 1,s + MAXN);

	return;

}

ll c[MAXN];

int lowbit(int x){return x & (-x);}

void add(int x,ll k){for(int i = x;i < MAXN;i += lowbit(i))c[i] = (c[i] + k) % MOD;return;}

ll sum(int x){ll res = 0;for(int i = x;i >= 1;i -= lowbit(i))res = (res + c[i]) % MOD;return res;}

int main()

{

	scanf("%d",&T);

	for(int i = 1;i <= T;++i)

	{

		scanf("%d%d%d",&q[i].n,&q[i].m,&q[i].a);

		q[i].id = i;

	}

	sort(q + 1,q + 1 + T,cmpq);

	init();

	for(int t = 1,j = 0;t <= T;++t)

	{

		while(j + 1 < MAXN && s[j + 1].first <= q[t].a)

		{

			for(int k = 1;s[j + 1].second * k < MAXN;++k)

			{

				add(s[j + 1].second * k,s[j + 1].first * mu[k]);

			}

			++j;

		}

		int n = q[t].n,m = q[t].m;

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

		ll ans = 0;

		for(int l = 1,r = 0;l <= n;l = r + 1)

		{

			r = min(n / (n / l),m / (m / l));

			ans += (sum(r) - sum(l - 1) + MOD) % MOD * (ll)(n / l) % MOD * (ll)(m / l) % MOD;

			ans %= MOD;

		}

		q[t].ans = ans;

	}

	sort(q + 1,q + 1 + T,cmpid);

	for(int i = 1;i <= T;++i)

	{

		printf("%lld\n",q[i].ans);

	}

	return 0;

}

In tag: 数学-莫比乌斯反演

wjh15101051