公式

A - Two Arithmetic Progressions 解説 by toam


解法 1

互除法の要領で,\(C=0\) に帰着することができます.

while C != 0:
    int k = floor(A / C)
    A = A - k * C
    B = B - k * D
    swap(A, C)
    swap(B, D)

この処理後の \(A,B,C,D\)\(A',B',C'(=0),D'\) と置くと,答えは \(\displaystyle \sum_{i=1}^N \gcd(A'i+B',|D'|)\) です.ここで,\(|D'|\leq 10^8\) です(\(|AD-BC|\) は常に一定であるため).

\(D'=0\) のときは簡単です.

\(D'\neq 0\) のときは,\(D'\) の各正の約数 \(m\) について,\(A'i+B'\)\(m\) の倍数になるような \(i(1\leq i\leq N)\) の個数を不定方程式を解くことで求めることができます.あとは約数包除することで,\(\gcd(Ai+B,|D'|)=m\) となる \(m\) の個数を求めることができます.

\(10^8\) 以下の整数の約数の個数は 高々 \(768\) 個なので,約数包除に \(2\) 乗時間かけても間に合います.


解法 2

\(C(Ai+B)-A(Ci+D)=(CB-AD)\) なので,最大公約数は \(|CB-AD|\) の約数です.\(CB-AD=0\) のときは簡単です.

\(|CB-AD|\) の各約数 \(m\) について,\(Ai+B,Ci+D\) がともに \(m\) の倍数になるような \(i\) の条件は,不定方程式を繰り返し解くことで求めることができます.あとは解法 \(1\) と同様です.


解法 1 の実装例

from atcoder.math import inv_mod
from math import gcd


def divisors(n):
    lower, upper = [], []
    i = 1
    while i * i <= n:
        if n % i == 0:
            lower.append(i)
            if i * i != n:
                upper.append(n // i)
        i += 1
    return lower + upper[::-1]


mod = 998244353


def solve():
    n, A, B, C, D = map(int, input().split())
    B += A
    D += C
    while C != 0:
        k = A // C
        A -= k * C
        B -= k * D
        A, B, C, D = C, D, A, B
    if D == 0:
        print((n * (n - 1) // 2 * A + B * n) % mod)
        return
    divs = divisors(abs(D))
    sz = len(divs)
    cnt = [0] * sz
    ans = 0
    for i in range(sz - 1, -1, -1):
        m = divs[i]
        g = gcd(A, m)
        if B % g == 0:
            p, q, r = A // g, B // g, m // g
            mn = (-q * inv_mod(p, r)) % r
            cnt[i] = max(0, (n - 1 - mn) // r + 1)
        for j in range(i + 1, sz):
            if divs[j] % m == 0:
                cnt[i] -= cnt[j]
        ans += m * cnt[i]
    print(ans % mod)


for _ in range(int(input())):
    solve()

投稿日時:
最終更新: