Ý tưởng

Gọi `nu(t)` là số các ước của số tự nhiên `t`

Xét phương trình: `n^2+x^2=y^2 => n^2=y^2-x^2=(y-x)(y+x)`

Mà `y-x` và `y+x` là các số nguyên

`=>` `y-x` và `y+x` là ước của `n^2`

Thêm vào đó, ta nhận thấy `y + x > y - x > 0` và `y > x > 0`

Do `y+x>y-x>0`, `y+x,y-x` là hai ước của `n^2` và cặp nghiệm `x,y` sao cho `y+x=y-x` không thỏa bài toán nên số lượng cặp `(y+x;y-x)` thỏa mãn là `floor((nu(n^2))/2)`

Lần lượt gọi `y+x=a` và `y-x=b` `(a > b)`

`=>` `{(2y=a+b),(x=y-b):}`

`=>` `{(y=(a+b)/2),(x=(a+b)/2-b=(a-b)/2):}`

Dễ thấy cả `x` và `y` dều dương

Lại xét: `y=(a+b)/2`; mà `y` là số nguyên dương `=>` `(a+b) \ vdots \ 2`

(Ta không cần xét cho `x` do `(a+b) \ vdots \ 2` thì `(a-b) \ vdots \ 2`)

`+)` Với `n` là số lẻ thì `n^2` là số lẻ nên mọi ước của `n^2` là số lẻ

`=>` Hai số `a,b` là các số lẻ

`=>` `(a+b) \ vdots 2`

`=>` Kết quả: `floor((nu(n^2))/2)`

`+)` Với `n` là số chẵn thì `n^2` là số chẵn

Khi đó ta có thể viết `n^2=2^r*d` với `r` nguyên dương và `d` là số lẻ

Ta sẽ phải trừ đi các trường hợp khi `a` hoặc `b` là ước của `d` vì khi đó trong `a` và `b` sẽ có một số chẵn và một số lẻ dẫn đến `a+b` không chia hết cho `2`

`=>` Kết quả: `floor((nu(n^2))/2)-nu(d)`

Do giới hạn của `n` khá to (lên đến `10^12`) lên ta sẽ vận dụng thuật toán đếm ước trong `O(n^(1/3))` (Bạn có thể tham khảo trên mạng nha)

Code

from math import sqrt
def sieve(maxn):
    isPrime = [True]*(maxn + 1)
    isPrime[0] = False
    isPrime[1] = False
    i = 2
    while i*i <= maxn:
        if isPrime[i]:
            j = i*i
            while j <= maxn:
                isPrime[j] = False
                j += i
        i += 1
    primes = []
    for i in range(2, maxn + 1):
        if isPrime[i]:
            primes.append(i)
    return primes
def test(a, n, k, m):
    mod = pow(a, m, n)
    if mod == 1 or mod == n - 1:
        return True
    for l in range(1, k):
        mod = (mod*mod) % n
        if mod == n - 1:
            return True
    return False
def isPrime(n):
    if n in [2,3,5,7,11,13,17,19,23,29,31,37]:
        return True
    elif n < 37:
        return False
    k = 0
    m = n - 1
    while m % 2 == 0:
        m //= 2
        k += 1
    checkset = [2,3,5,7,11,13,17,19,23,29,31,37]
    for a in checkset:
        if not test(a, n, k, m):
            return False
    return True
n = int(input())
primes = sieve(10**4)
m = 1 if n % 2 == 0 else 0
d = 1
res = 1
for p in primes:
    if p**3 > n:
        break
    cnt = 0
    while n % p == 0:
        n //= p
        cnt += 1
    if p % 2 != 0:
        d *= 2*cnt + 1
    res *= 2*cnt + 1
if isPrime(n):
    res *= 3
    if n % 2 != 0:
        d *= 3
else:
    squareRoot = int(sqrt(n))
    if squareRoot*squareRoot == n and isPrime(squareRoot):
        if n % 2 != 0:
            d *= 5
        res *= 5
    elif n != 1:
        if n % 2 != 0:
            d *= 9
        else:
            d *= 3
        res *= 9
print(res // 2 - d*m)