Skip to main content
1 of 3
Peter Mueller
  • 22.5k
  • 1
  • 75
  • 107

Below is the Python code which for instance computes $\nu(6)$. (Initially, I had written a SageMath script using Sturm's Theorem, however the polsturm function from gp/pari is somewhat faster.)

The script, if it is named mo.py, is best called in a Linux console via nohup python3 mo.py nr cpus start &, where nr is a lower bound for the number of real roots we want, cpus is the number of kernels we want to use, and start should be $i$ for $0\le i\le cpus-1$.

For instance, in order to compute $\nu(6)$ without parallelization, just call python3 mo.py 6 1 0.

The script identifies the binary expansion $p=\sum a_i2^i$ with the polynomial $\sum a_ix^i$. Note that we need not look at polynomials which are divisible by $x^2$. Another reduction, by about the factor $2$, uses that the number of nonzero roots of a polynomial equals the number of its reciprocal.

In case that someone wants to look for a pattern, here is a list of all $0$-$1$-polynomials of degree $28$ with $6$ distinct real roots. To save space, we give the list as the sequence (see above) of numbers: 437545878, 437594646, 437594838, 437619606, 437643798, 440494998, 440875158, 443812758, 443836950, 443837334, 443910678, 443911062, 443916438, 443923062, 443923350, 444033942, 450177558, 462908310. (In SageMath, if K is a polynomial ring, and p is one of these numbers, K(p.bits()) yields the corresponding polynomial.)

And here is the Python code (which requires the library cypari2):

from sys import argv
import cypari2
pari = cypari2.Pari()

nr, cpus, start = [int(z) for z in argv[1:]]

def int2pol(p):
    if p&3 == 0:
        return None
    i = 0
    l = []
    c = p
    rec = 0
    while c != 0:
        if c&1 != 0:
            l.append(f'x^{i}')
        rec = (rec << 1) + (c & 1)
        i += 1
        c >>= 1
    if p&1 == 0:
        rec <<= 1
    return pari('+'.join(l)) if p <= rec else None

n = -1
file = open(f'mo_461631_nr={nr}_cpus={cpus}_start={start}.log', 'w')
while True:
    n += 1
    print(f'degree = {n}')
    file.write(f'degree = {n}\n')
    file.flush()
    for p in range(2**n + start, 2**(n+1), cpus):
        f = int2pol(p)
        if f != None:
            m = f.polsturm()
            if m >= nr:
                print(m, p, f)
                file.write(str((m, p, f)) + '\n')
                file.close()
                break
    else:
        continue
    break
Peter Mueller
  • 22.5k
  • 1
  • 75
  • 107