## Trial Division

### October 25, 2016

Here is our solution:

```(define (factors n)
(let loop ((n n) (f 2) (fs (list)))
(if (< n (* f f))
(reverse (cons n fs))
(if (zero? (modulo n f))
(loop (/ n f) f (cons f fs))
(loop n (+ f 1) fs)))))```

There’s no simpler way to factor integers. Here are some examples:

```> (factors 16)
(2 2 2 2)
> (factors 17)
(17)
> (factors 13290059)
(3119 4261)
> (factors 600851475143)
(71 839 1471 6857)```

You can run the program at http://ideone.com/r2xHK0, where you will also see the primality checker from the previous exercise, which has exactly the same form and structure.

Pages: 1 2

### 2 Responses to “Trial Division”

1. Alex B said

Very basic Python 3 solution that divides out by 2 and then by each odd number until the reduced n becomes less than the square of the factor to be tested.
I used Python’s builtin divmod instead of separate modulo and division.

```def factors_td(n):
'''Return a list of the prime factors of n'''
factors = []

while n % 2 == 0:
factors.append(2)
n //= 2

f = 3
while n > f*f:
q, r = divmod(n, f)
if r == 0:
factors.append(f)
n = q
else:
f += 2
if n != 1:
factors.append(n)
return factors
```
2. matthew said

Here’s another Python solution, to make things more interesting, it factorizes Gaussian integers (ie. integral complex numbers):

```import heapq

# Usual complex division but only for integers
def div(z,w):
(a,b),(c,d) = z,w
e = a*c + b*d
f = b*c - a*d
k = c*c + d*d
if e%k == 0 and f%k == 0: return (e//k,f//k)

# norm(a*b) = norm(a)*norm(b)
def norm(z):
(a,b) = z
return a*a + b*b

def factorize(z):
h, factors = [],[]
# Heap of all pairs, ordered by norm
def extend(z):
(a,b) = z
extend((1,0))
while True:
(n,w) = heapq.heappop(h)
if n*n > norm(z): break
extend(w)
while True:
quot = div(z,w)
if quot is None: break
factors.append(w)
z = quot
if z != (1,0): factors.append(z)
return factors

for a in range(1,5):
for b in range(0,5):
z = (a,b)
print("%s: %s"%(z,factorize(z)))
```
```(1, 0): []
(1, 1): [(1, 1)]
(1, 2): [(1, 2)]
(1, 3): [(1, 1), (2, 1)]
(1, 4): [(1, 4)]
(2, 0): [(1, 1), (1, 1), (0, -1)]
(2, 1): [(2, 1)]
(2, 2): [(1, 1), (1, 1), (1, 1), (0, -1)]
(2, 3): [(2, 3)]
(2, 4): [(1, 1), (1, 1), (2, -1)]
(3, 0): [(3, 0)]
(3, 1): [(1, 1), (2, -1)]
(3, 2): [(3, 2)]
(3, 3): [(1, 1), (3, 0)]
(3, 4): [(2, 1), (2, 1)]
(4, 0): [(1, 1), (1, 1), (1, 1), (1, 1), (-1, 0)]
(4, 1): [(4, 1)]
(4, 2): [(1, 1), (1, 1), (1, -2)]
(4, 3): [(1, 2), (1, 2), (0, -1)]
(4, 4): [(1, 1), (1, 1), (1, 1), (1, 1), (1, 1), (-1, 0)]
```