zoukankan      html  css  js  c++  java
  • 杜教筛学习笔记

    前置知识

    积性函数

    积性函数分为积性函数和完全积性函数。

    • 积性函数:(forall{aperp{b}},f(a*b)=f(a)*f(b))。常见的有(varphi.mu,sigma,d)
    • 完全积性函数:(forall{a,b},f(a*b)=f(a)*f(b))。常见的有(epsilon,I,id)

    其中(epsilon(n)=[n=1],I(n)=1,id(n)=n)

    狄利克雷卷积

    ((f*g)(n)=sumlimits_{d|n}f(d)g(frac{n}{d})),满足交换律和结合律。

    定义(epsilon)为单位元((f*epsilon=f),根据定义易得)

    (egin{cases}mu*I=epsilon\varphi*I=id\mu*id=varphiend{cases})

    莫比乌斯反演

    (Link)

    杜教筛

    原理

    (sumlimits_{i=1}^nf(i)=S(n))

    再找一个积性函数(g),考虑((f*g))的前缀和。

    [sumlimits_{i=1}^{n}(f*g)(i) ]

    [=sumlimits_{i=1}^nsumlimits_{d|i}f(d)g(frac{i}{d}) ]

    [=sumlimits_{d=1}^ng(d)sumlimits_{i=1}^{lfloorfrac{n}{d} floor}f(i) ]

    [=sumlimits_{d=1}^ng(d)S(lfloorfrac{n}{d} floor) ]

    又因为

    [g(1)S(n) ]

    [=sumlimits_{d=1}^ng(d)S(lfloorfrac{n}{d} floor)-sumlimits_{d=2}^ng(d)S(lfloorfrac{n}{d} floor) ]

    [=sumlimits_{i=1}^{n}(f*g)(i)-sumlimits_{d=2}^ng(d)S(lfloorfrac{n}{d} floor) ]

    如果我们能找到合适的(g)来快速算出(sumlimits_{i=1}^{n}(f*g)(i)),剩下的就能够整除分块了,从而可以递归求解。

    复杂度为(O(n^{frac{3}{4}})),预处理前(m)项使复杂度达到(O(frac{n}{sqrt{m}}))(m)(n^{frac{2}{3}})时达到最优(O(n^{frac{2}{3}}))

    实践

    (f(n)=mu(n)),取(g(n)=I),则((f*g)(n)=epsilon(n))

    (f(n)=varphi(n)),取(g(n)=I),则((f*g)(n)=id(n))

    注意若(n=2^{31}-1)时,整除分块中(l=r+1)可能会爆(int),改为(unsigned int)即可。

    下面是模板代码:

    #include <bits/stdc++.h>
    
    using namespace std;
    
    #define ll long long
    
    const int lim = (1 << 22), mx = 2147483647; 
    
    int t, n, tot, mu[4200005], phi[4200005], vis[4200005], pr[1250005], prs_mu[4200005];
    
    ll prs_phi[4200005];
    
    map < int, int > fmu;
    map < int, ll > fphi;
    
    int read()
    {
    	int x = 0, fl = 1; char ch = getchar();
    	while (ch < '0' || ch > '9') { if (ch == '-') fl = -1; ch = getchar();}
    	while (ch >= '0' && ch <= '9') {x = x * 10ll + ch - '0'; ch = getchar();}
    	return x * fl;
    }
    
    void init()
    {
    	phi[1] = mu[1] = 1;
    	for (int i = 2; i <= lim; i ++ )
    	{
    		if (!vis[i])
    		{
    			vis[i] = i;
    			pr[ ++ tot] = i;
    			phi[i] = i - 1;
    			mu[i] = -1;
    		}
    		for (int j = 1; j <= tot; j ++ )
    		{
    			if (i * pr[j] > lim || pr[j] > vis[i]) break;
    			vis[i * pr[j]] = pr[j];
    			if (i % pr[j])
    			{
    				mu[i * pr[j]] = -mu[i];
    				phi[i * pr[j]] = phi[i] * phi[pr[j]];
    			}
    			else
    			{
    				mu[i * pr[j]] = 0;
    				phi[i * pr[j]] = phi[i] * pr[j];
    			}
    		}
    	}
    	for (int i = 1; i <= lim; i ++ )
    	{
    		prs_mu[i] = prs_mu[i - 1] + mu[i];
    		prs_phi[i] = prs_phi[i - 1] + phi[i];
    	}
    	return;
    }
    
    int d_mu(unsigned int x)
    {
    	if (x <= lim) return prs_mu[x];
    	if (fmu[x]) return fmu[x];
    	int cnt = 1;
    	for (unsigned int l = 2, r; l <= x; l = r + 1)
    	{
    		r = min(x, x / (x / l));
    		cnt = cnt - (r - l + 1) * d_mu(x / l);
    	}
    	return fmu[x] = cnt;
    }
    
    ll d_phi(unsigned int x)
    {
    	if (x <= lim) return prs_phi[x];
    	if (fphi[x]) return fphi[x];
    	ll cnt = 1ll * x * (x + 1) / 2;
    	for (unsigned int l = 2, r; l <= x; l = r + 1)
    	{
    		r = min(x, x / (x / l));
    		cnt = cnt - (r - l + 1) * d_phi(x / l);
    	}
    	return fphi[x] = cnt;
    }
    
    int main()
    {
    	init();
    	t = read();
    	while (t -- )
    	{
    		n = read();
    		printf("%lld %d
    ", d_phi(n), d_mu(n));
    	}
    	return 0;
    }
    
  • 相关阅读:
    创建供应商-采购模块
    定义容差组
    前台创建供应商-财务角度
    对供应商账户组分配编号范围
    拓端数据tecdat|R语言建立和可视化混合效应模型mixed effect model
    拓端数据tecdat|R语言建模收入不平等:分布函数拟合及洛伦兹曲线(Lorenz curve)
    拓端数据tecdat|R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型
    拓端数据tecdat|R语言ARIMA,SARIMA预测道路交通流量时间序列:季节性、周期性
    拓端数据tecdat|ARIMA模型预测CO2浓度时间序列
    拓端数据tecdat|R语言基于递归神经网络RNN的温度时间序列预测
  • 原文地址:https://www.cnblogs.com/andysj/p/14507575.html
Copyright © 2011-2022 走看看