zoukankan      html  css  js  c++  java
  • Solution -「51nod 1584」加权约数和

    (mathcal{Description})

      Link.

      令 (sigma(n))(n) 的约数之和。求:

    [sum_{i=1}^nsum_{j=1}^nmax{i,j}sigma(ij)mod(10^9+7) ]

      多测,(nle10^6),数据组数 (le5 imes10^4)

    (mathcal{Solution})

      直 接 来 owo!

    [sum_{i=1}^nsum_{j=1}^nmax{i,j}sigma(ij)=2sum_{i=1}^nisum_{j=1}^isigma(ij)-sum_{i=1}^nisigma(i^2) ]

      先研究一下 (sigma(ij))。考虑分别枚举 (i) 的约数 (x)(j) 的约数 (y),若 (xperpfrac{j}y),则对 (sigma(ij)) 贡献一个 (xy),显然贡献不重不漏。即:

    [sigma(ij)=sum_{x|i}sum_{y|j}xy[xperpfrac{j}y] ]

      考虑原式前一项,记 (f(n)=nsum_{i=1}^nsigma(ni)),有:

    [egin{aligned} f(n)&=nsum_{i=1}^nsum_{x|n}sum_{y|i}xy[xperpfrac{i}y]\ &=nsum_{i=1}^nsum_{x|n}sum_{y|i}xysum_{d|xland d|frac{i}y}mu(d)\ &=nsum_{i=1}^nsum_{d|nland d|i}mu(d)sum_{x|nland d|x}sum_{y|iland d|y}frac{ix}y\ &=nsum_{i=1}^nsum_{d|nland d|i}mu(d)sum_{x|frac{n}d}sum_{y|frac{i}d}frac{ix}y~~~~~~~~(x,ymbox{ 同时约掉 } d)\ &=nsum_{i=1}^nsum_{d|nland d|i}mu(d)sigma(frac{n}d)sum_{y|frac{i}d}frac{i}y\ &=nsum_{i=1}^nsum_{d|nland d|i}mu(d)sigma(frac{n}d)dsum_{y|frac{i}d}frac{frac{i}y}d\ &=nsum_{i=1}^nsum_{d|nland d|i}dmu(d)sigma(frac{n}d)sigma(frac{i}d)\ &=nsum_{d|n}dmu(d)sigma(frac{n}d)sum_{i=1}^frac{n}{d}sigma(i) end{aligned} ]

      筛出 (mu,sigma),枚举 (d)(frac{n}d),可以 (mathcal O(nln n)) 算出所有 (f)

      后一项呢,就是要筛 (sigma(n^2))。和筛 (sigma(n)) 类似,记录一下当前最小素因子的等比数列求和,就可以 (mathcal O(n)) 算出来。

      综上,复杂度 (mathcal O(nln n+T))

    (mathcal{Code})

    /* Clearink */
    
    #include <cstdio>
    
    typedef long long LL;
    
    inline int rint () {
    	int x = 0, f = 1; char s = getchar ();
    	for ( ; s < '0' || '9' < s; s = getchar () ) f = s == '-' ? -f : f;
    	for ( ; '0' <= s && s <= '9'; s = getchar () ) x = x * 10 + ( s ^ '0' );
    	return x * f;
    }
    
    const int MAXN = 1e6, MOD = 1e9 + 7;
    int pn, pr[MAXN + 5], mpwr[MAXN + 5], mu[MAXN + 5], sig[MAXN + 5], sigs[MAXN + 5];
    int f[MAXN + 5], ans[MAXN + 5];
    LL dpwr[MAXN + 5], g[MAXN + 5];
    bool vis[MAXN + 5];
    
    inline void init ( const int n ) {
    	mu[1] = sig[1] = sigs[1] = g[1] = 1;
    	for ( int i = 2; i <= n; ++ i ) {
    		if ( !vis[i] ) {
    			mu[pr[++ pn] = i] = -1;
    			sig[i] = mpwr[i] = i + 1;
    			dpwr[i] = g[i] = 1ll * i * i + i + 1;
    		}
    		for ( int j = 1, t; j <= pn && ( t = i * pr[j] ) <= n; ++ j ) {
    			vis[t] = true;
    			if ( !( i % pr[j] ) ) {
    				mpwr[t] = mpwr[i] * pr[j] + 1;
    				dpwr[t] = dpwr[i] * pr[j] * pr[j] + pr[j] + 1;
    				sig[t] = sig[i] / mpwr[i] * mpwr[t];
    				g[t] = g[i] / dpwr[i] * dpwr[t];
    				break;
    			}
    			mu[t] = -mu[i];
    			mpwr[t] = mpwr[pr[j]];
    			dpwr[t] = dpwr[pr[j]];
    			sig[t] = sig[i] * sig[pr[j]];
    			g[t] = g[i] * g[pr[j]];
    		}
    		sigs[i] = ( sigs[i - 1] + sig[i] ) % MOD;
    	}
    	for ( int i = 1; i <= n; ++ i ) g[i] = i * g[i] % MOD;
    	for ( int i = 1; i <= n; ++ i ) {
    		for ( int j = 1, t = n / i; j <= t; ++ j ) {
    			f[i * j] = ( f[i * j] + 1ll * i * mu[i] * sig[j] % MOD * sigs[j] ) % MOD;
    		}
    	}
    	for ( int i = 1; i <= n; ++ i ) {
    		f[i] = 1ll * i * ( f[i] + MOD ) % MOD;
    		ans[i] = ( ( ans[i - 1] + 2ll * f[i] - g[i] ) % MOD + MOD ) % MOD;
    	}
    }
    
    int main () {
    	init ( MAXN );
    	for ( int T = rint (), i = 1; i <= T; ++ i ) {
    		printf ( "Case #%d: %d
    ", i, ans[rint ()] );
    	}
    	return 0;
    }
    

    (mathcal{Details})

      突然觉得推式子好养生啊。(

  • 相关阅读:
    SkylineGlobe 如何实现绘制圆形Polygon和对图层的圆形范围选择查询
    如何修改Oracle服务IP地址
    TerraGate软件安装后,不能启动的解决办法
    SkylineGlobe SFS发布的WFS和WMS服务测试
    SkylineGlobe 如何实现二次开发加载KML文件
    如何屏蔽SkylineGlobe提供的三维地图控件上的快捷键
    SkylineGlobe 如何实现FlyTo定位到目标点之后触发的事件函数
    转子百度知道:LTE上行采用DC子载波,而下行置零,为什么?
    转载:LTE小区搜索过程
    dB dBc dBi dBd dBm dBW 定义
  • 原文地址:https://www.cnblogs.com/rainybunny/p/13685651.html
Copyright © 2011-2022 走看看