Loj#6682. 梦中的数论(min25筛)

Loj#6682. 梦中的数论(min25筛),第1张

Loj#6682. 梦中的数论(min25筛)

传送门:https://loj.ac/p/6682
题目大意:给定 n n n,求 ∑ i = 1 n ∑ j = 1 n ∑ k = 1 n [ j ∣ i ∧ ( j + k ) ∣ i ] sumlimits_{i=1}^nsumlimits_{j=1}^nsumlimits_{k=1}^{n}[jmid iland (j+k)mid i] i=1∑n​j=1∑n​k=1∑n​[j∣i∧(j+k)∣i],其中 1 ≤ n ≤ 1 0 10 1leq nleq 10^{10} 1≤n≤1010。
题解:不难发现对于 i i i的一个因数 j j j,上式相等于计数还有多少个比 j j j大的 i i i的因数,那么可以转化为:
∑ i = 1 n ∑ j ∣ i ∑ j < k , k ∣ i 1 = ∑ i = 1 n ( τ ( i ) − 1 + τ ( i ) − 2 + ⋯ + 2 + 1 + 0 ) = ∑ i = 1 n τ ( i ) ( τ ( i ) − 1 ) 2 = 1 2 ( ∑ i = 1 n τ 2 ( i ) − τ ( i ) ) = 1 2 ( S ⁡ τ 2 ( n ) − S ⁡ τ ( n ) ) begin{aligned} &sumlimits_{i=1}^nsumlimits_{jmid i}sumlimits_{j 其中 τ ( n ) tau(n) τ(n)是约数个数函数, S ⁡ f ( n ) = ∑ i = 1 n f ( i ) operatorname{S}_{f}(n)=sumlimits_{i=1}^nf(i) Sf​(n)=i=1∑n​f(i)。现在的问题就是要如何快速求出 S ⁡ τ 2 ( n ) operatorname{S}_{tau^2}(n) Sτ2​(n)了。方法有很多,这里我只用min25筛去做,如果有空的话就补上杜教筛和PN筛的写法(
考虑 τ 2 ( n ) tau^2(n) τ2(n)在质数次方时的取值,显然 τ 2 ( p k ) = ( k + 1 ) 2 tau^2(p^k)=(k+1)^2 τ2(pk)=(k+1)2且对于所有质数都恒有 τ 2 ( p ) = 4 = 4 × 1 tau^2(p)=4=4times1 τ2(p)=4=4×1。那么对于min25中的 g ( i , j ) g(i,j) g(i,j)函数,可以作如下转移: g ( i , j ) = g ( i , j − 1 ) − ( g ( ⌊ i p j ⌋ , j − 1 ) − s u m ( p j − 1 ) ) ,     p j 2 ≤ i g(i,j)=g(i,j-1)-(g(lfloorfrac{i}{p_j}rfloor,j-1)-sum(p_{j-1})),spacespacespace p_j^2leq i g(i,j)=g(i,j−1)−(g(⌊pj​i​⌋,j−1)−sum(pj−1​)),   pj2​≤i且初始状态 g ( i , 0 ) = 4 ( i − 1 ) g(i,0)=4(i-1) g(i,0)=4(i−1)。对于 s ( i , j ) s(i,j) s(i,j),有 s ( i , j ) = g ( i ) − s u m ( p j − 1 ) + ∑ k = j + 1 p k 2 ≤ i ∑ r = 1 p k r ≤ i ( r + 1 ) 2 ( s ( ⌊ i p k r ⌋ , k ) + [ r > 1 ] ) s(i,j)=g(i)-sum(p_{j-1})+sumlimits_{k=j+1}^{p_k^2leq i}sumlimits_{r=1}^{p_k^rleq i}(r+1)^2(s(lfloorfrac{i}{p_k^r}rfloor,k)+[r>1]) s(i,j)=g(i)−sum(pj−1​)+k=j+1∑pk2​≤i​r=1∑pkr​≤i​(r+1)2(s(⌊pkr​i​⌋,k)+[r>1])。
剩下的 S ⁡ τ ( n )   O ⁡ ( n ) operatorname{S}_{tau}(n)space operatorname{O}(sqrt n) Sτ​(n) O(n ​)快速求即可。
(不知道为啥数组改int会tle,奇奇怪怪)

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
using namespace std;
#define debug(x) cerr<<#x<<":"< vis;
unordered_map g;
void pre(ll n) {
    ll sqn = sqrt(n);

    for (int i = 2; i <= sqn; ++i) {
        if (!vis[i]) {
            sum[i] = add(sum[i - 1], 4);
            prime[++cnt] = i;
        } else
            sum[i] = sum[i - 1];

        for (int j = 1; i * prime[j] <= sqn; ++j) {
            vis[i * prime[j]] = 1;

            if (i % prime[j] == 0)
                break;
        }
    }

    for (ll l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        g[n / l] = mul(4, sub(n / l, 1));
    }

    for (int i = 1; i <= cnt; ++i) {
        for (ll l = 1, r; l <= n && prime[i]*prime[i] <= n / l; l = r + 1) {
            ll t = n / l;
            r = n / t;
            g[t] = sub(g[t], sub(g[t / prime[i]], sum[prime[i - 1]]));
        }
    }
}
ll s(ll i, ll j) {
    if (prime[j] > i)
        return 0;

    ll res = sub(g[i], sum[prime[j]]);

    for (ll k = j + 1; k <= cnt && prime[k]*prime[k] <= i; ++k) {
        for (ll r = 1, t = prime[k]; t <= i; ++r, t *= prime[k]) {
            res = add(res, mul(mul(r + 1, r + 1), add(s(i / t, k), r > 1)));
        }
    }

    return res;
}
int main() {
    ll n;
    cin >> n;
    pre(n);
    ll res = s(n, 0) + 1;

    for (ll l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        res = sub(res, mul(sub(r, sub(l, 1)), n / l));
    }

    res = mul(res, M / 2 + 1);
    cout << res << endl;
}

欢迎分享,转载请注明来源:内存溢出

原文地址: http://outofmemory.cn/zaji/5670054.html

(0)
打赏 微信扫一扫 微信扫一扫 支付宝扫一扫 支付宝扫一扫
上一篇 2022-12-16
下一篇 2022-12-16

发表评论

登录后才能评论

评论列表(0条)

保存