【BZOJ】2956:模积和

Time Limit: 10 Sec  Memory Limit: 128 MB

Description

求∑∑((n mod i)*(m mod j))其中1<=i<=n,1<=j<=m,i≠j。

Input

第一行两个数n,m。

Output

一个整数表示答案mod 19940417的值

Sample Input

3 4

Sample Output

1


样例说明
  
答案为(3 mod 1)*(4 mod 2)+(3 mod 1) * (4 mod 3)+(3 mod 1) * (4 mod 4) + (3 mod 2) * (4 mod 1) +

(3 mod 2) * (4 mod 3) + (3 mod 2) * (4 mod 4) + (3 mod 3) * (4 mod 1) + (3 mod 3) * (4 mod 2) +

(3 mod 3) * (4 mod 4) = 1

数据规模和约定
  对于100%的数据n,m<=10^9

 

这里是链接:【BZOJ】2956:模积和

这里是题解:

首先,暴力枚举将会很凄惨:O(nm)。早就 GG ( Time Limit Exceeded )了。

所以从公式入手:原公式是:TIM图片20171020172807_thumb[11]

展开为:TIM图片20171020173411_thumb[2]

观察式子:[n mod i](同理 m mod j)

根据mod的定义可以将上式写成这个样子:TIM图片20171020173710_thumb[1]

[n/i]向下取整,就是C++中的整型(n/i),然后再乘以 i 就相当于下图灰色区域

再用n减掉就能得到mod后的值下图模拟mod的转化

TIM图片20171020174120_thumb[3]

所以式子就可以简化为:TIM图片20171020174328_thumb[3]

注意:因为题目中i!=j,所以当i、j 相同就直接减掉】

TIM图片20171020174618_thumb[2]

然而这样还是O(n^2)的复杂度。

所以继续化简:将第一个式子的Σ移动,使时间复杂度变为O(n)

这里是最终式子:

TIM图片20171020175357_thumb[2]

能够移动Σ的证明设n-[n/i]*i为Xi,m-[m/j]*j为Yj。

(下图有点错误:n、m的值应该是不同的,但是n、m不同的证明也是这样的。)

   这个是原式子展开的ans:

               TIM图片20171020175536_thumb[2]

  这个是化简后展开的ans:

TIM图片20171020175656_thumb[2]

显然它们的ans值是相等的。那么,第一步化简式子已经完成了。

虽然移动Σ已经将复杂度降低到O(n),但很不幸的是依然过不了

考虑如何优化

低于O(n)的复杂度一般就三种:O(1)、O(logn)、O(√n)

注意最终式子,都有一个式子like this:[n/i](其中n为一个定值,i是从1到n的一个变量。)

但是这里有个很美妙的事情就是:

假如n==1000时,i在91~100之间,n/i的值都是为10的;

假如n==100000000时,i在9090910~10000000之间,n/i的值都是为10的。

这里因为在某一个区间中的值都是相等的。所以我们可以很愉快地利用分块的思想

那么怎么算这个块的大小呢?

假设有一个块里面的[n/i]的值都是为k,那么其区间就是:[(n/(k+1))+1,n/k].

推导:因为[n/i]是向下取整的,所以k*i<=n,n/k为定值,所以i<=n/k,但i一定也有下界,所以i>n/(k+1),

即i>=n/(k+1)+1.

注意:分块求值套用公式的时候需要除法,并不能先取模,然而不先取模会爆long long.

但好在除的数是固定的6,所以就直接在求平方和的时候,MOD开大6倍,最后再模回去就行了。

(其实反过来也就是网上普遍流传的3323403,是[mod/6]的值)

这里是代码:

 

#include<iostream>
#include<cstring>
#include<cstdio>
#include<algorithm>
#define LL long long 
#define mod 19940417
#define MOD 119642502
using namespace std;
LL n,m,tmp1,tmp2,ans,ine;

LL sum(LL x){
    return x*(x+1)/2%mod;
}

LL SUM(LL x){
    return x*(x+1)%MOD*(2*x+1)%MOD/6;
}

LL calc(LL x){
    LL ans=x*x%mod;
    for(LL i=1;i<=x;i=tmp1+1){
        tmp1=x/(x/i);//(x/i)求k的值,n/k为块的上界
         ans=(ans-(x/i)*(sum(tmp1)-sum(i-1)+mod)%mod)%mod; 
    }
    return ans%mod;
}

int main(){
    scanf("%lld %lld",&n,&m);
    if(n>m) swap(n,m);
    ans=calc(n)*calc(m)%mod;
    //处理i,j相同情况
    for(LL i=1;i<=n;i=tmp2+1){
        tmp2=min(n/(n/i),m/(m/i));
        ans=(ans-n*m%mod*(tmp2-i+1)%mod
                +m*(n/i)%mod*(sum(tmp2)-sum(i-1)+mod)%mod
                +n*(m/i)%mod*(sum(tmp2)-sum(i-1)+mod)%mod
                -(n/i)*(m/i)%mod*(SUM(tmp2)-SUM(i-1)+mod)%mod+2*mod)%mod;
    }
    printf("%lld",ans%mod);
    return 0;
}
【BZOJ】2956

 

 

梦想总是要有的,万一实现了呢?

 

posted @ 2017-10-20 20:54  Σ1102  阅读(398)  评论(0编辑  收藏  举报