def main():

    mod = 998244353

    max_A = 10**6

    N = int(eval(input()))

    A = list(map(int, input().split()))



    is_prime = [False, True, False, False, False, True] * (max_A//6+1)

    del is_prime[max_A+1:]

    is_prime[1:4] = False, True, True

    for i in range(5, int(max_A**0.5)+1):

        if is_prime[i]:

            is_prime[i*i::i] = [False] * (max_A//i-i+1)

    primes = [p for p, is_p in enumerate(is_prime) if is_p]  # max_A 以下の素数のリスト



    g = [0] * (max_A+1)

    for a in A:

        g[a] += a

    for p in primes:

        # 倍数集合の高速ゼータ変換みたいなやつ  O((max_A)loglog(max_A))

        # 参考: http://noshi91.hatenablog.com/entry/2018/12/27/121649

        # 大量に約数列挙みたいなことをするときはこれで高速化できる場合が多そう？（みんぷろ 2018 本戦 A - Uncommon など）

        for k in range(max_A//p, 0, -1):

            g[k] += g[k*p]

    g = [v * v % mod for v in g]

    for p in primes:

        for k, g_kp in enumerate(g[p::p], 1):

            g[k] -= g_kp



    modinv = [0, 1]

    for i in range(2, max_A+1):

        modinv.append(mod - mod//i * modinv[mod%i] % mod)

    ans = sum((gg * minv % mod for gg, minv in zip(g, modinv)))

    ans %= mod

    ans -= sum(A)

    ans *= modinv[2]

    ans %= mod

    print(ans)

main()

# 解説 AC



# 参考1: https://twitter.com/rickytheta/status/1175412019006074880

# 参考2: https://maspypy.com/atcoder-%E5%8F%82%E5%8A%A0%E6%84%9F%E6%83%B3-2019-09-21agc-038



def make_prime_checker(n):

    is_prime = [False, True, False, False, False, True] * (n//6+1)

    del is_prime[n+1:]

    is_prime[1:4] = False, True, True

    for i in range(5, int(n**0.5)+1):

        if is_prime[i]:

            is_prime[i*i::i] = [False] * (n//i-i+1)

    return is_prime



def make_modinv_list(n, mod=10**9+7):

    # 0 から n までの mod 逆元のリストを返す O(n)

    modinv = [0, 1]

    for i in range(2, n+1):

        modinv.append(mod - mod//i * modinv[mod%i] % mod)

    return modinv



def main():

    mod = 998244353

    max_A = 10**6

    N = int(eval(input()))

    A = list(map(int, input().split()))

    

    # nまでの自然数が素数かどうかを表すリストを返す  O(nloglogn)

    

    primes = [p for p, is_prime in enumerate(make_prime_checker(max_A)) if is_prime]  # max_A 以下の素数のリスト

    

    g = [0] * (max_A+1)

    for a in A:

        g[a] += a

    for p in primes:

        # 倍数集合の高速ゼータ変換みたいなやつ  O((max_A)loglog(max_A))

        # 参考: http://noshi91.hatenablog.com/entry/2018/12/27/121649

        # 大量に約数列挙みたいなことをするときはこれで高速化できる場合が多そう？（みんぷろ 2018 本戦 A - Uncommon など）

        for k in range(max_A//p, 0, -1):

            g[k] += g[k*p]

    

    # この時点で g[d] = Σ_{d|a} a  (A の要素のうち d の倍数であるものの総和)

    

    g = [v * v % mod for v in g]

    

    # この時点で

    # g[d] = (Σ_{d|a} a)(Σ_{d|b} b)

    #      = Σ_{(d|a)∧(d|b)} ab

    #      = Σ_{d | gcd(a,b)} ab

    

    # この式変形天才すぎないか？（これを使った一連の操作を gcd 畳み込み というらしい？）

    

    for p in primes:

        # 倍数集合の高速メビウス変換みたいなやつ  O((max_A)loglog(max_A))

    

        # for k in range(1, max_A//p+1):

        #     g[k] -= g[k*p]  # 包除原理ヤバい

    

        for k, g_kp in enumerate(g[p::p], 1):  # 高速化（ゼータ変換の方は途中で g が変化するので高速化できない）

            g[k] -= g_kp  # 包除原理ヤバい

    

    # この時点で g[d] = Σ_{gcd(a,b)=d} ab

    

    modinv_list = make_modinv_list(max_A, mod)

    ans = sum((gg * minv % mod for gg, minv in zip(g, modinv_list)))

    ans %= mod

    

    # この時点で

    # ans = Σ_d ( Σ_{gcd(a,b)=d} ab/d )

    #     = Σ_d ( Σ_{gcd(a,b)=d} lcm(a,b) )

    #     = Σ_a Σ_b lcm(a,b)

    

    ans -= sum(A)

    ans *= modinv_list[2]

    ans %= mod

    

    # この時点で ans = Σ_{i<j} lcm(a,b)

    

    print(ans)



main()

def main():
    mod = 998244353
    max_A = 10**6
    N = int(eval(input()))
    A = list(map(int, input().split()))
    is_prime = [False, True, False, False, False, True] * (max_A // 6 + 1)
    del is_prime[max_A + 1 :]
    is_prime[1:4] = False, True, True
    for i in range(5, int(max_A**0.5) + 1):
        if is_prime[i]:
            is_prime[i * i :: i] = [False] * (max_A // i - i + 1)
    primes = [p for p, is_p in enumerate(is_prime) if is_p]  # max_A 以下の素数のリスト
    g = [0] * (max_A + 1)
    for a in A:
        g[a] += a
    for p in primes:
        # 倍数集合の高速ゼータ変換みたいなやつ  O((max_A)loglog(max_A))
        # 参考: http://noshi91.hatenablog.com/entry/2018/12/27/121649
        # 大量に約数列挙みたいなことをするときはこれで高速化できる場合が多そう？（みんぷろ 2018 本戦 A - Uncommon など）
        for k in range(max_A // p, 0, -1):
            g[k] += g[k * p]
    g = [v * v % mod for v in g]
    for p in primes:
        for k, g_kp in enumerate(g[p::p], 1):
            g[k] -= g_kp
    modinv = [0, 1]
    for i in range(2, max_A + 1):
        modinv.append(mod - mod // i * modinv[mod % i] % mod)
    ans = sum((gg * minv % mod for gg, minv in zip(g, modinv)))
    ans %= mod
    ans -= sum(A)
    ans *= modinv[2]
    ans %= mod
    print(ans)


main()

# 解説 AC
# 参考1: https://twitter.com/rickytheta/status/1175412019006074880
# 参考2: https://maspypy.com/atcoder-%E5%8F%82%E5%8A%A0%E6%84%9F%E6%83%B3-2019-09-21agc-038
def make_prime_checker(n):
    is_prime = [False, True, False, False, False, True] * (n // 6 + 1)
    del is_prime[n + 1 :]
    is_prime[1:4] = False, True, True
    for i in range(5, int(n**0.5) + 1):
        if is_prime[i]:
            is_prime[i * i :: i] = [False] * (n // i - i + 1)
    return is_prime


def make_modinv_list(n, mod=10**9 + 7):
    # 0 から n までの mod 逆元のリストを返す O(n)
    modinv = [0, 1]
    for i in range(2, n + 1):
        modinv.append(mod - mod // i * modinv[mod % i] % mod)
    return modinv


def main():
    mod = 998244353
    max_A = 10**6
    N = int(eval(input()))
    A = list(map(int, input().split()))
    # nまでの自然数が素数かどうかを表すリストを返す  O(nloglogn)
    primes = [
        p for p, is_prime in enumerate(make_prime_checker(max_A)) if is_prime
    ]  # max_A 以下の素数のリスト
    g = [0] * (max_A + 1)
    for a in A:
        g[a] += a
    for p in primes:
        # 倍数集合の高速ゼータ変換みたいなやつ  O((max_A)loglog(max_A))
        # 参考: http://noshi91.hatenablog.com/entry/2018/12/27/121649
        # 大量に約数列挙みたいなことをするときはこれで高速化できる場合が多そう？（みんぷろ 2018 本戦 A - Uncommon など）
        for k in range(max_A // p, 0, -1):
            g[k] += g[k * p]
    # この時点で g[d] = Σ_{d|a} a  (A の要素のうち d の倍数であるものの総和)
    g = [v * v % mod for v in g]
    # この時点で
    # g[d] = (Σ_{d|a} a)(Σ_{d|b} b)
    #      = Σ_{(d|a)∧(d|b)} ab
    #      = Σ_{d | gcd(a,b)} ab
    # この式変形天才すぎないか？（これを使った一連の操作を gcd 畳み込み というらしい？）
    for p in primes:
        # 倍数集合の高速メビウス変換みたいなやつ  O((max_A)loglog(max_A))
        # for k in range(1, max_A//p+1):
        #     g[k] -= g[k*p]  # 包除原理ヤバい
        for k, g_kp in enumerate(g[p::p], 1):  # 高速化（ゼータ変換の方は途中で g が変化するので高速化できない）
            g[k] -= g_kp  # 包除原理ヤバい
    # この時点で g[d] = Σ_{gcd(a,b)=d} ab
    modinv_list = make_modinv_list(max_A, mod)
    ans = sum((gg * minv % mod for gg, minv in zip(g, modinv_list)))
    ans %= mod
    # この時点で
    # ans = Σ_d ( Σ_{gcd(a,b)=d} ab/d )
    #     = Σ_d ( Σ_{gcd(a,b)=d} lcm(a,b) )
    #     = Σ_a Σ_b lcm(a,b)
    ans -= sum(A)
    ans *= modinv_list[2]
    ans %= mod
    # この時点で ans = Σ_{i<j} lcm(a,b)
    print(ans)


main()

+# 解説 AC +# 参考1: https://twitter.com/rickytheta/status/1175412019006074880 +# 参考2: https://maspypy.com/atcoder-%E5%8F%82%E5%8A%A0%E6%84%9F%E6%83%B3-2019-09-21agc-038 +def make_prime_checker(n): +    is_prime = [False, True, False, False, False, True] * (n // 6 + 1) +    del is_prime[n + 1 :] +    is_prime[1:4] = False, True, True +    for i in range(5, int(n**0.5) + 1): +        if is_prime[i]: +            is_prime[i * i :: i] = [False] * (n // i - i + 1) +    return is_prime + + +def make_modinv_list(n, mod=10**9 + 7): +    # 0 から n までの mod 逆元のリストを返す O(n) +    modinv = [0, 1] +    for i in range(2, n + 1): +        modinv.append(mod - mod // i * modinv[mod % i] % mod) +    return modinv + + -    is_prime = [False, True, False, False, False, True] * (max_A // 6 + 1) -    del is_prime[max_A + 1 :] -    is_prime[1:4] = False, True, True -    for i in range(5, int(max_A**0.5) + 1): -        if is_prime[i]: -            is_prime[i * i :: i] = [False] * (max_A // i - i + 1) -    primes = [p for p, is_p in enumerate(is_prime) if is_p]  # max_A 以下の素数のリスト +    # nまでの自然数が素数かどうかを表すリストを返す  O(nloglogn) +    primes = [ +        p for p, is_prime in enumerate(make_prime_checker(max_A)) if is_prime +    ]  # max_A 以下の素数のリスト +    # この時点で g[d] = Σ_{d|a} a  (A の要素のうち d の倍数であるものの総和) +    # この時点で +    # g[d] = (Σ_{d|a} a)(Σ_{d|b} b) +    #      = Σ_{(d|a)∧(d|b)} ab +    #      = Σ_{d | gcd(a,b)} ab +    # この式変形天才すぎないか？（これを使った一連の操作を gcd 畳み込み というらしい？） -        for k, g_kp in enumerate(g[p::p], 1): -            g[k] -= g_kp -    modinv = [0, 1] -    for i in range(2, max_A + 1): -        modinv.append(mod - mod // i * modinv[mod % i] % mod) -    ans = sum((gg * minv % mod for gg, minv in zip(g, modinv))) +        # 倍数集合の高速メビウス変換みたいなやつ  O((max_A)loglog(max_A)) +        # for k in range(1, max_A//p+1): +        #     g[k] -= g[k*p]  # 包除原理ヤバい +        for k, g_kp in enumerate(g[p::p], 1):  # 高速化（ゼータ変換の方は途中で g が変化するので高速化できない） +            g[k] -= g_kp  # 包除原理ヤバい +    # この時点で g[d] = Σ_{gcd(a,b)=d} ab +    modinv_list = make_modinv_list(max_A, mod) +    ans = sum((gg * minv % mod for gg, minv in zip(g, modinv_list))) +    # この時点で +    # ans = Σ_d ( Σ_{gcd(a,b)=d} ab/d ) +    #     = Σ_d ( Σ_{gcd(a,b)=d} lcm(a,b) ) +    #     = Σ_a Σ_b lcm(a,b) -    ans *= modinv[2] +    ans *= modinv_list[2] +    # この時点で ans = Σ_{i<j} lcm(a,b)