原題鏈接

?

最近被51NOD的數論題各種刷……（NOI快到了我在干什么啊！

然后發現這題在網上找不到題解……那么既然A了就來騙一波訪問量吧……

（然而并不怎么會用什么公式編輯器，寫得丑也湊合著看吧……

?

?

$$
ANS=\sum_{i=1}^{n}\sum_{j=1}^{n} \frac{i*j}{gcd(i,j)}
$$
$$
=\sum_{d=1}^{n} d*\sum_{i=1}^{\left\lfloor\frac{n}ze8trgl8bvbq\right\rfloor}?\sum_{j=1}^{\left\lfloor\frac{n}ze8trgl8bvbq\right\rfloor} i*j*[gcd(i,j)==1]
$$
$$
=\sum_{d=1}^{n} d*\sum_{d'=1}^{\left\lfloor\frac{n}ze8trgl8bvbq\right\rfloor} f(\left\lfloor\frac{n}{d*d'}\right\rfloor)*d'^2*μ(d')
$$
$$
=\sum_{d'=1}^{n} d'^2*μ(d')*\sum_{d=1}^{\left\lfloor\frac{n}{d'}\right\rfloor} d*f(\left\lfloor\frac{n}{d*d'}\right\rfloor)
$$

?

這里$f(x)=(\frac{x*(x+1)}{2})^2$，即1到x中數字之間兩兩乘積之和。

?

可以看出不同的$\left\lfloor\frac{n}{d'}\right\rfloor$只有根號n個，對于某一個$\left\lfloor\frac{n}{d'}\right\rfloor$，不同的$\left\lfloor\frac{n}{d*d'}\right\rfloor$也只有根號n個，那么把它們壓在一起處理就好了。

但是要做到這一點需要快速求出$d'^2*μ(d')$的前綴和，這時就要用上杜教篩了。

記$g(x)=x^2*μ(x)$，$S(x)=\sum_{i=1}^{x} g(i)$

那么

$ ?\sum_{i=1}^{n} \sum_{x|i} g(x)*(\frac{i}{x})^2$

$=\sum_{i=1}^{n} ?i*\sum_{x|i}?μ(x)$

$=\sum_{i=1}^{n} [i==1] $

$= 1$

同時

$ ?\sum_{i=1}^{n} \sum_{x|i} g(x)*(\frac{i}{x})^2$

$=\sum_{i=1}^{n}?\sum_{x=1}^{[\frac{n}{i}]} g(x)*i^2$

$=\sum_{i=1}^{n}?i^2*S[\frac{n}{i}]$

可以得到$S(n)=1-\sum_{i=2}^{n} i^2*S[\frac{n}{i}]$，預處理+哈希維護一波就行了

但是不知道是我常數太大，還是方法沒別人優越，卡了很久常數才卡過去。

?

#include<ctime> #include<cstdio> #include<algorithm> #define ll long long #define N 1000001 #define K 500000004 using namespace std; ll read_p,read_ca,read_f; inline ll read(){read_p=0;read_ca=getchar();read_f=1;while(read_ca<'0'||read_ca>'9') {if (read_ca=='-') read_f=-1;read_ca=getchar();}while(read_ca>='0'&&read_ca<='9') read_p=read_p*10+read_ca-48,read_ca=getchar();return read_p*read_f; } const int HA=1e6+7; const int MOD=1e9+7; int MMH=0,p[N],num=0,mu[N],_ff[N],l[HA],Num,f[N],s[N],U=0; ll n,P; bool bo[N]; struct na{int z,ne;ll y;}b[HA]; inline void M(int &x){while (x>=MOD) x-=MOD;while(x<0)x+=MOD;} inline ll MM(ll x){return x<MOD&&x>-MOD?x:x%MOD;} inline void in(int x,ll y,int z){b[++Num].y=y;b[Num].z=z;b[Num].ne=l[x];l[x]=Num;} inline int _f(ll x){if (x>MOD) x%=MOD;if (x<N) return _ff[x];return MM(1LL*x*(x+1)>>1); } inline ll sqr(ll x){return x*x%MOD;} inline int F(ll x){int MMH=0;ll P;for (register ll i=1,j;i<=x;i=j+1) j=x/(P=x/i),M(MMH+=MM(sqr(_f(P))*(_f(j)-_f(i-1))));return MMH; } inline int Q(ll x){if (x>MOD) x%=MOD;return MM(MM(x*(x+1))*(x+x+1))*166666668%MOD;} inline int mmh(ll x){if (x<N) return s[x];register ll i,j;ll P;int o;for (i=l[o=x%HA];i;i=b[i].ne) if (b[i].y==x) return b[i].z;int MMH=0;for (i=2;i<=x;i=j+1) j=x/(P=x/i),M(MMH=MM(1LL*(Q(j)-Q(i-1))*mmh(P)+MMH));M(MMH=1-MMH);in(o,x,MMH);return MMH; } int main(){register int i,j;mu[1]=1;for (i=2;i<N;i++){if (!bo[i]) p[++num]=i,mu[i]=-1;for (j=1;j<=num&&i*p[j]<N;j++){bo[i*p[j]]=1;if (i%p[j]) mu[i*p[j]]=-mu[i];else {mu[i*p[j]]=0;break;}}}for (i=1;i<N;i++) M(f[i]=1LL*i*i*mu[i]%MOD+MOD),M(s[i]=s[i-1]+f[i]),_ff[i]=MM(1LL*i*(i+1)>>1);n=read();for (ll i=1,j;i<=n;i=j+1) j=n/(P=n/i),M(MMH+=MM(1LL*(mmh(j)-mmh(i-1))*F(P)));printf("%d\n",MMH); } View Code

?

?

?

?

轉載于:https://www.cnblogs.com/Enceladus/p/5655267.html

總結

以上是生活随笔為你收集整理的51 NOD 1238 最小公倍数之和 V3的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。