zoukankan      html  css  js  c++  java
  • BZOJ4176 Lucas的数论 【莫比乌斯反演 + 杜教筛】

    题目

    去年的Lucas非常喜欢数论题,但是一年以后的Lucas却不那么喜欢了。

    在整理以前的试题时,发现了这样一道题目“求Sigma(f(i)),其中1<=i<=N”,其中 表示i的约数个数。他现在长大了,题目也变难了。
    求如下表达式的值:

    其中 表示ij的约数个数。
    他发现答案有点大,只需要输出模1000000007的值。

    输入格式

    第一行一个整数n。

    输出格式

    一行一个整数ans,表示答案模1000000007的值。

    输入样例

    2

    输出样例

    8

    提示

    对于100%的数据n <= 10^9。

    题解

    这题推导和SDOI2015约数个数和那道题是一样的
    只不过计算的方式有差别

    这道题没有多组询问,而且n特别大【不能O(n)实现】,要用杜教筛

    最后推出式子:

    [ans=sumlimits_{d=1}^{n} mu(d) (sumlimits_{i=1}^{lfloor frac{n}{d} floor} lfloor frac{lfloor frac{n}{d} floor}{i} floor)^2 ]

    如果我们记

    [sum(n) = sumlimits_{i=1}^{n}lfloor frac{n}{i} floor ]

    那么式子可以写成:

    [ans=sumlimits_{d=1}^{n} mu(d) sum(lfloor frac{n}{d} floor)^2 ]

    显然可以分块计算
    因为(n<=10^9),所以对于(mu)的前缀和我们采用杜教筛,时间复杂度(O(n^{frac{2}{3}}logn))
    对于(sum(n)),我们内部也分块计算,时间复杂度(O(int_{0}^{sqrt{n}} x^{frac{1}{2}} dx) = O(frac{2}{3} n^{frac{3}{4}}) = O(n^{frac{3}{4}}))

    所以总的复杂度(O(n^{frac{2}{3}}logn + n^{frac{3}{4}}))

    #include<iostream>
    #include<cstdio>
    #include<cmath>
    #include<map>
    #include<cstring>
    #include<algorithm>
    #define LL long long int
    #define Redge(u) for (int k = h[u],to; k; k = ed[k].nxt)
    #define REP(i,n) for (int i = 1; i <= (n); i++)
    #define BUG(s,n) for (int i = 1; i <= (n); i++) cout<<s[i]<<' '; puts("");
    using namespace std;
    const int maxn = 1000005,maxm = 100005,INF = 1000000000,P = 1000000007;
    typedef map<LL,LL> Map;
    Map _mu;
    Map::iterator it;
    LL p[maxn],pi,mu[maxn],N,n;
    int isn[maxn];
    void init(LL n){
    	N = (LL)pow(n,2.0 / 3.0);
    	mu[1] = 1;
    	for (int i = 2; i < N; i++){
    		if (!isn[i]) p[++pi] = i,mu[i] = -1;
    		for (int j = 1; j <= pi && i * p[j] < N; j++){
    			isn[i * p[j]] = true;
    			if (i % p[j] == 0){
    				mu[i * p[j]] = 0;
    				break;
    			}
    			mu[i * p[j]] = -mu[i];
    		}
    	}
    	for (int i = 1; i < N; i++) mu[i] = (mu[i - 1] + mu[i]) % P;
    }
    LL sum(LL x){
    	LL ans = 0;
    	for (int i = 1,nxt; i <= x; i = nxt + 1){
    		nxt = x / (x / i);
    		ans = (ans + (nxt - i + 1) * (x / i) % P) % P;
    	}
    	return ans;
    }
    LL S(LL n){
    	if (n < N) return mu[n];
    	if ((it = _mu.find(n)) != _mu.end())
    		return it->second;
    	LL ans = 1;
    	for (int i = 2,nxt; i <= n; i = nxt + 1){
    		nxt = n / (n / i);
    		ans = (ans - (nxt - i + 1) * S(n / i) % P) % P;
    	}
    	return _mu[n] = ans;
    }
    int main(){
    	cin >> n;
    	init(n);
    	LL ans = 0;
    	for (int i = 1,nxt; i <= n; i = nxt + 1){
    		nxt = n / (n / i);
    		LL tmp = sum(n / i);
    		tmp = tmp * tmp % P;
    		ans = ans + (S(nxt) - S(i - 1)) % P * tmp % P;
    	}
    	ans = (ans % P + P) % P;
    	cout << ans << endl;
    	return 0;
    }
    
    
  • 相关阅读:
    个人网站上线
    从零开始学C++之虚函数与多态(二):纯虚函数、抽象类、虚析构函数
    Triangle
    Hard problem
    [置顶] Android AlarmManager实现不间断轮询服务
    Python学习入门基础教程(learning Python)--5.1 Python下文件处理基本过程
    编程之美 裴波那楔数列
    130712周赛(CF)
    Python学习入门基础教程(learning Python)--5.2 Python读文件基础
    phantomjs
  • 原文地址:https://www.cnblogs.com/Mychael/p/8746620.html
Copyright © 2011-2022 走看看