Day 118-Q2 | 強力数法（Powerful Number Method）

問題

正の整数 $N$ が与えられる。1以上N以下の整数のうち平方因子を持たない数（squarefree number）の総和

$$S(N) = \sum_{n=1}^{N} n \cdot [n \text{ は平方因子を持たない}]$$

を 998244353 で割った余りを求めよ。

入力形式

制約

$1 \le N \le 10^{11}$

入出力例

入力例1

出力例1

入力例2

出力例2

1,2,3,5,6,7,10が平方因子を持たない数で、その和は$1+2+3+5+6+7+10=34$。4,8,9は平方因子($2^2$または$3^2$)を持つため除外。

概念図: 強力数への帰着

ヒント（段階的開示）

ヒント1: 方向性

素直にN個を全部調べるとO(N)で、$N\le10^{11}$には全く間に合わない。$f(n)=n\cdot[n\text{が平方因子を持たない}]$は乗法的関数（互いに素な$a,b$に対し$f(ab)=f(a)f(b)$）である。乗法的関数の総和を高速に求める代表的な枠組みの一つに「強力数法（Powerful Number Method）」がある。

ヒント2: アプローチ

強力数法は、計算したい乗法的関数$f$に対して「素数$p$で$f(p)$と一致し、累積和が簡単に求まる乗法的関数$g$」を選び、$h=f*g^{-1}$（ディリクレ逆元との畳み込み）を考える。素数$p$で$h(p)=0$になるように選んだので、$h$が非零になるのは「全ての素因数が2乗以上の指数を持つ数（強力数）」に限られ、強力数はN以下にO(√N)個しか存在しない。このとき $S(N)=\sum_{d:\text{強力数},d\le N}h(d)\cdot G(\lfloor N/d\rfloor)$ が成り立つ。今回は$g(n)=n$（累積和$G(x)=x(x+1)/2$）を選ぶと、非零な強力数は「平方因子を持たない$d$の2乗$d^2$」の形に限られ、$h(d^2)=\mu(d)\cdot d^2$になる。

ヒント3: 誘導（コード骨格）

# 1. B = floor(sqrt(N)) までのメビウス関数 mu[1..B] を線形篩で前計算
# 2. answer = 0
# 3. d = 1..B について、mu[d] == 0 でなければ
#      m = N // (d*d)
#      G(m) = m*(m+1)/2 mod MOD  （逆元 inv2 を使う）
#      answer += mu[d] * d^2 * G(m)   (mod MOD)
# 4. answer mod MOD を出力

模範解答 (Python)

import sys

MOD = 998244353


def solve():
    N = int(sys.stdin.readline())

    B = int(N ** 0.5)
    while (B + 1) * (B + 1) <= N:
        B += 1
    while B * B > N:
        B -= 1

    mu = [0] * (B + 1)
    if B >= 1:
        mu[1] = 1
    primes = []
    is_comp = [False] * (B + 1)
    for i in range(2, B + 1):
        if not is_comp[i]:
            primes.append(i)
            mu[i] = -1
        for p in primes:
            if i * p > B:
                break
            is_comp[i * p] = True
            if i % p == 0:
                mu[i * p] = 0
                break
            else:
                mu[i * p] = -mu[i]

    inv2 = pow(2, MOD - 2, MOD)
    ans = 0
    for d in range(1, B + 1):
        if mu[d] == 0:
            continue
        m = N // (d * d)
        Gm = (m % MOD) * ((m + 1) % MOD) % MOD * inv2 % MOD
        term = mu[d] * (d * d % MOD)
        ans = (ans + term * Gm) % MOD

    print(ans % MOD)


solve()

計算量: 線形篩 $O(\sqrt N)$、強力数の列挙が$O(\sqrt N)$、全体で$O(\sqrt N)$（一般の強力数法は$O(N^{2/3})$だが本問は単純な形に落ちるためこれより速い）。ブルートフォースと200通り以上の乱択$N$で一致を確認済み。$N=10^{11}$でも0.1秒程度で完了することを確認済み。

Step-by-Step 解説

1乗法的関数の総和問題とは
$f(n)$が乗法的で$f(p^e)$が高速に計算できるとき、N以下の総和を高速に求める一般的な枠組み(Min_25篩・杜教篩・強力数法など)の一つが強力数法である。

2gの選び方とhの導出
$g(p)=f(p)$となる$g$を選ぶと、ディリクレ畳み込みの逆元$h=f*g^{-1}$は$h(p)=0$を満たし、強力数にのみ非零値を持つ。

3本問での具体化
$g(n)=n$とすると$h(p)=0$。$f(p^e)=0(e\ge2)$から逆算すると$h(p^2)=-p^2$、$h(p^e)=0(e\ne2)$となり、非零な強力数は「平方因子を持たない$d$の2乗$d^2$」のみで$h(d^2)=\mu(d)d^2$。

4累積和G(x)の利用
$g(n)=n$の累積和は$G(x)=x(x+1)/2$で定数時間。$d=1\dots\sqrt N$について$h(d^2)G(\lfloor N/d^2\rfloor)$を足し合わせるだけで$S(N)$が求まる。

5計算量の確認
全体で$O(\sqrt N)$。$N\le10^{11}$でも高速に完了する。

よくあるミス

ミス	原因	正しい書き方
Bの計算を`int(N**0.5)`のみに頼り境界で1ずれる	浮動小数点誤差でsqrtが正確な整数境界からずれることがある	計算後に(B+1)²≤NやB²>Nをwhileループで補正する
G(m)の計算でmをそのままMODで割ろうとする	m自体は10¹¹規模でMOD未満とは限らないと錯覚する	m%MOD, (m+1)%MODにしてから掛け算し、最後にinv2を掛ける
μ(d)=0のdもd²として計算に含めてしまう	強力数の条件を見落とす	`if mu[d]==0: continue`で明示的にスキップする
強力数法を「N^(2/3)が常に最速」と思い込む	一般論の計算量だけを覚えてしまう	単純な形に帰着できる場合はO(√N)まで落ちることもあると理解する

次のステップ

発展: $f(n) = n^2\cdot[n\text{が平方因子を持たない}]$ や、より一般の乗法的関数に強力数法を適用してみる。
次回予告: 縮約階層法（Contraction Hierarchies）

自己評価

完全理解おおむね理解要復習

Day 118-Q2 — 強力数法（Powerful Number Method）

問題