Revision history [back]

There is a command for exactly that in Sage, it is called trial_division.

It will find the smallest factor of a number, optionally limited to some bound.

For example,

sage: a = 207411784165069788403685789342707568622276103

Find the smallest factor:

sage: trial_division(a)
19

Find the smallest factor if it at most some bound, otherwise return the unfactored number:

sage: trial_division(a, bound=19)
19
sage: trial_division(a, bound=18)
207411784165069788403685789342707568622276103

One could use that function repeatedly to "factor as you go" until the user interrupts the factoring attempt.

For that, one way would be to define the following functions:

def inshort(m, l=5):
    numbertype = "primelike" if m.is_prime(proof=False) else "composite"
    return ("({:" + "{}".format(l) + "}-digit " +
            "{})").format(m.ndigits(), numbertype)

@cached_function
def trial_div(a, b):
    d = a.trial_division(start=b)
    k = a.valuation(d)
    return (d, k)

def trial_factor(n, start=None):
    print("\nTrying primelike factorisation of\n# = {}\n".format(n))
    if n == 1:
        print("# = 1 has no prime factors")
    m = Integer(n)
    d = Integer(1 if start is None else start)
    factors = []
    l = m.ndigits().ndigits()
    format_exponent = lambda k: '^{}'.format(k) if k > 1 else ''
    format_factor = lambda d, k: '{}{}'.format(d, format_exponent(k))
    sfirst = " * found new prime divisor: "
    snext = " * (prime divisors above) * "
    print("# = {}".format(inshort(m, l=l)))
    try:
        while m != 1:
            d, k = trial_div(m, d + 1)
            m = m // d^k
            if m > 1:
                a = inshort(m, l=l)
                b = (snext if factors else sfirst) + format_factor(d, k)
                print("# = {}{}".format(a, b))
            factors.append((d, k))
            if m == 1:
                print('\nComplete factorisation:\n# = ' +
                      ' * '.join(format_factor(*f) for f in factors))
                return(factors)
    except KeyboardInterrupt:
        print('\nFactoring interrupted. So far, obtained:\n# = ' +
              ' * '.join(format_factor(*f) for f in factors) +
              ' * {}'.format(inshort(m, l=l)))
        return factors + [(m, 1)]

Sometimes we get something that is probably a prime factorisation:

sage: a = 207411784165069788403685789342707568622276103
sage: f = trial_factor(a)

Trying primelike factorisation of
# = 207411784165069788403685789342707568622276103

# = (45-digit composite)
# = (44-digit composite) * found new prime divisor: 19
# = (42-digit composite) * (prime divisors above) * 29
# = (41-digit composite) * (prime divisors above) * 37
# = (38-digit composite) * (prime divisors above) * 929
# = (32-digit primelike) * (prime divisors above) * 518717
^C
Factoring interrupted. So far, obtained:
# = 19 * 29 * 37 * 929 * 518717 * (32-digit primelike)
where the final (32-digit primelike) is
21112220300468085766530282304433
sage: f
[(19, 1),
 (29, 1),
 (37, 1),
 (929, 1),
 (518717, 1),
 (21112220300468085766530282304433, 1)]

Sometimes the remaining factor is still known to be composite.

sage: a
872277058918929207644191411339288743538056112242
sage: f = trial_factor(a)

Trying primelike factorisation of
# = 872277058918929207644191411339288743538056112242

# = (48-digit composite)
# = (48-digit composite) * found new prime divisor: 2
# = (48-digit composite) * (prime divisors above) * 3
^C
Factoring interrupted. So far, obtained:
# = 2 * 3 * (48-digit composite)
where the final (48-digit composite) is
145379509819821534607365235223214790589676018707
sage: f
[(2, 1), (3, 1), (145379509819821534607365235223214790589676018707, 1)]