project-euler/Sidef/745 Sum of Squares -- v2.sf

#!/usr/bin/ruby

# Daniel "Trizen" Șuteu
# Date: 31 January 2021
# https://github.com/trizen

# Sum of Squares
# https://projecteuler.net/problem=745

# Formula:
#   S(n) = Sum_{k=1..floor(sqrt(n))} k^2 * R(floor(n/k^2))

# Where R(n) is the number of squarefree numbers <= n:
#   R(n) = Sum_{k=1..floor(sqrt(n))} moebius(k) * floor(n/k^2)

# Faster formula:
#   S(n) = Sum_{k=1..floor(sqrt(n))} J_2(k) * floor(n / k^2)
# where J_n(x) is the Jordan totient function.

# S(10^1)  = 24
# S(10^2)  = 767
# S(10^3)  = 22606
# S(10^4)  = 722592
# S(10^5)  = 22910120
# S(10^6)  = 725086120
# S(10^7)  = 22910324448
# S(10^8)  = 724475280152
# S(10^9)  = 22907428923832
# S(10^10) = 724420596049320
# S(10^11) = 22908061437420776
# S(10^12) = 724418227020757048
# S(10^13) = 22908104289912800016

# Runtime: 1 minute, 53 seconds (previously: ~3 minutes).

define MOD = 1_000_000_007

func S(n) {
    sum(1..n.isqrt, {|k|
        mulmod(idiv(n, k.sqr), k.jordan_totient(2), MOD)
    })
}

say (S(10**14) % MOD)