ICode9

精准搜索请尝试: 精确搜索
首页 > 其他分享> 文章详细

[SDOI 2015] 约数个数和

2021-01-13 16:01:14  阅读:203  来源: 互联网

标签:约数 lfloor right frac SDOI sum rfloor times 2015


\(\text{Description}\)

传送门

\(\text{Solution}\)

首先有:

\[\sigma(x\times y)=\sum_{i|x}\sum_{j|y} [\gcd(i,j)=1] \]

具体证明戳这

柿子变成了:

\[\sum_{i=1}^n\sum_{j=1}^m\sum_{x|i}\sum_{y|j} [\gcd(x,y)=1] \]

而常见形式是:

\[\sum_{i=1}^n\sum_{j=1}^m [\gcd(i,j)=1] \]

尝试将 \(x,y\) 提出去,向常见形式靠拢。

其实可以发现 \(\sum_{i=1}^n\sum_{j=1}^m\) 就是个限制条件,它决定了次数。

\[\sum_{x=1}^n\sum_{y=1}^m[\gcd(x,y)=1]\times \left \lfloor\frac{n}{x}\right \rfloor\times \left \lfloor\frac{m}{y}\right \rfloor \]

设,

\[g(k)=\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=k]\times \left \lfloor\frac{n}{i}\right \rfloor\times \left \lfloor\frac{m}{j}\right \rfloor \]

\[f(k)=\sum_{k|d}g(d) \]

所以,

\[\text{Ans}=g(1)=\sum_{i=1}^{\min\{n,m\}}\mu(i)\times f(i) \]

我们需要快速求出 \(f(i)\)(根据经验,应该是可以整除分块的)。

不过这次似乎无法像 \(\sum_{i=1}^n\sum_{j=1}^m [\gcd(i,j)=1],\sum_{i=1}^n\sum_{j=1}^m [\gcd(i,j)=1]\times i\times j\) 这样的柿子猛拆二重循环,因为有 \(\left \lfloor\frac{n}{i}\right \rfloor\times \left \lfloor\frac{m}{j}\right \rfloor\) 碍事的。

尝试直接化柿子:

\[f(k)=\sum_{i=1}^n\sum_{j=1}^m[k|\gcd(i,j)]\times \left \lfloor\frac{n}{i}\right \rfloor\times \left \lfloor\frac{m}{j}\right \rfloor \]

显然我们无法简化 \(i,j\),显然不能枚举 \(\gcd\)。

其实可以将 \(k\) 提出来,这样 \(\gcd\) 肯定是 \(k\) 的倍数:

\[=\sum_{i=1}^{\left \lfloor\frac{n}{k}\right \rfloor}\sum_{j=1}^{\left \lfloor\frac{m}{k}\right \rfloor}\left \lfloor\frac{n}{i\times k}\right \rfloor\times \left \lfloor\frac{m}{j\times k}\right \rfloor \]

预处理 \(\sum_{i=1}^n \left \lfloor\frac{n}{i}\right \rfloor\) 就行了。

时间复杂度 \(\mathcal O((T+n)\times \sqrt n)\)。

\(\text{Code}\)

#include <cstdio>

#define rep(i,_l,_r) for(register signed i=(_l),_end=(_r);i<=_end;++i)
#define fep(i,_l,_r) for(register signed i=(_l),_end=(_r);i>=_end;--i)
#define erep(i,u) for(signed i=head[u],v=to[i];i;i=nxt[i],v=to[i])
#define efep(i,u) for(signed i=Head[u],v=to[i];i;i=nxt[i],v=to[i])
#define print(x,y) write(x),putchar(y)

template <class T> inline T read(const T sample) {
    T x=0; int f=1; char s;
    while((s=getchar())>'9'||s<'0') if(s=='-') f=-1;
    while(s>='0'&&s<='9') x=(x<<1)+(x<<3)+(s^48),s=getchar();
    return x*f;
}
template <class T> inline void write(const T x) {
    if(x<0) return (void) (putchar('-'),write(-x));
    if(x>9) write(x/10);
    putchar(x%10^48);
}
template <class T> inline T Max(const T x,const T y) {if(x>y) return x; return y;}
template <class T> inline T Min(const T x,const T y) {if(x<y) return x; return y;}
template <class T> inline T fab(const T x) {return x>0?x:-x;}
template <class T> inline T gcd(const T x,const T y) {return y?gcd(y,x%y):x;}
template <class T> inline T lcm(const T x,const T y) {return x/gcd(x,y)*y;}

typedef long long ll;

const int maxn=5e4+5;

int n,m,u[maxn],p[maxn],pc;
bool is[maxn];
ll f[maxn];

void init() {
	u[1]=1;
	rep(i,2,maxn-5) {
		if(!is[i]) p[++pc]=i,u[i]=-1;
		rep(j,1,pc) {
			if(p[j]*i>maxn-5) break;
			is[p[j]*i]=1;
			if(i%p[j]==0) break;
			u[i*p[j]]=-u[i];
		}
	} 
	rep(i,2,maxn-5) u[i]+=u[i-1];
	int r; ll ret;
	rep(i,1,maxn-5) {
		ret=0;
		for(int l=1;l<=i;l=r+1) r=i/(i/l),ret+=1ll*(r-l+1)*(i/l);
		f[i]=ret;
	}
}

ll Query() {
	int lim=Min(n,m),r; ll ret=0;
	for(int l=1;l<=lim;l=r+1) {
		r=Min(n/(n/l),m/(m/l));
		ret+=1ll*(u[r]-u[l-1])*f[n/l]*f[m/l];
	}
	return ret;
}

int main() {
	init();
	for(int T=read(9);T;--T) {
		n=read(9),m=read(9);
		print(Query(),'\n');
	}
	return 0;
} 

标签:约数,lfloor,right,frac,SDOI,sum,rfloor,times,2015
来源: https://www.cnblogs.com/AWhiteWall/p/14269865.html

本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享;
2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关;
3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关;
4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除;
5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

专注分享技术,共同学习,共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有