<!-- canonical: https://efficientnewlanguage.org/eml-p/projects/P066-polynomials/ | updated: 2026-10-14 -->

# P066 Polynomials

Two polynomials typed the usual way, with exact fraction coefficients: add, subtract, multiply, and divide with quotient and remainder (checked back), evaluate term by term and by Horner's rule with the multiplications counted, take the derivative, and find every rational root by the rational root theorem.

EML-P project `projects/polynomials` in the EML language repo: 3 module(s), entry `main.eml`, terminal UI. There, `eml project run projects/polynomials` runs it and `eml project verify projects/polynomials` replays every session under CPython (two hash seeds) and in the interpreter; the site build replays every session in the interpreter again and publishes a session only if its screen matches.

Built on verified corpus cases: polynomial-evaluator (https://efficientnewlanguage.org/cases/164-polynomial-evaluator/).

## Sessions

### bad-input - interpreter: byte-equal to the golden

Input:

```text
0
x
1

1
3 4
x^
2y
x^21
1/0
x--1

2
0
4
1
5
7
1
0
7
5
abc
-1/2
1
1000001x + 1000003
7
9
```

Screen:

```text
== Polynomials ==
Write a polynomial the usual way, like 3x^2 - 2x + 1/2; every coefficient is an exact fraction.
P = x^3 - 6x^2 + 11x - 6
Q = x - 1

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 0
Pick a number from 1 to 9.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> x
Pick a number from 1 to 9.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 1
P = 
Cancelled.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 1
P = 3 4
Cannot read that: two numbers side by side. Write terms like 3x^2, -x, 1/2 or 0.25x^3.
P = x^
Cannot read that: a power after ^ that is not a whole number. Write terms like 3x^2, -x, 1/2 or 0.25x^3.
P = 2y
Cannot read that: the coefficient 2y. Write terms like 3x^2, -x, 1/2 or 0.25x^3.
P = x^21
Cannot read that: a power above 20. Write terms like 3x^2, -x, 1/2 or 0.25x^3.
P = 1/0
Cannot read that: the coefficient 1/0. Write terms like 3x^2, -x, 1/2 or 0.25x^3.
P = x--1
Cannot read that: a sign with no term after it. Write terms like 3x^2, -x, 1/2 or 0.25x^3.
P = 
Cancelled.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 2
Q = 0
Q = 0, degree none (the zero polynomial).

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 4
Q is 0, and nothing divides by 0.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 1
P = 5
P = 5, degree 0.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 7
P is a nonzero constant: it has no roots.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 1
P = 0
P = 0, degree none (the zero polynomial).

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 7
P is the zero polynomial: every number is a root.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 5
x = abc
Type a number like 2, -3, 1/2 or 0.5.
x = -1/2
P is the zero polynomial: it is 0 everywhere.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 1
P = 1000001x + 1000003
P = 1000001x + 1000003, degree 1.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 7
Not tried: the coefficients are too large to try every candidate.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 9
Bye.
```

### basic - interpreter: byte-equal to the golden

Input:

```text
3
4
5
4
6
7
1
2x^4 - 3x^3 - 2x^2 + 3x
7
2
x^2 + 1
4
5
1/2
1
x^4 - 4x^3 + 6x^2 - 4x + 1
7
1
x^2 - 2
7
1
x^3 - 2x
7
8
9
```

Screen:

```text
== Polynomials ==
Write a polynomial the usual way, like 3x^2 - 2x + 1/2; every coefficient is an exact fraction.
P = x^3 - 6x^2 + 11x - 6
Q = x - 1

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 3
P + Q = x^3 - 6x^2 + 12x - 7
P - Q = x^3 - 6x^2 + 10x - 5
P * Q = x^4 - 7x^3 + 17x^2 - 17x + 6

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 4
P / Q: quotient x^2 - 5x + 6, remainder 0.
Check: quotient * Q + remainder = x^3 - 6x^2 + 11x - 6.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 5
x = 4
Term by term:
  -6 -> -6
  11x -> 44
  -6x^2 -> -96
  x^3 -> 64
  sum = 6, using 9 multiplications
Horner's rule, from the top coefficient down:
  1, -2, 3, 6
  P(4) = 6, using 3 multiplications

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 6
P' = 3x^2 - 12x + 11

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 7
Rational roots: 1, 2, 3 - 8 candidates tried.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 1
P = 2x^4 - 3x^3 - 2x^2 + 3x
P = 2x^4 - 3x^3 - 2x^2 + 3x, degree 4.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 7
Rational roots: -1, 0, 1, 3/2 - 8 candidates tried.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 2
Q = x^2 + 1
Q = x^2 + 1, degree 2.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 4
P / Q: quotient 2x^2 - 3x - 4, remainder 6x + 4.
Check: quotient * Q + remainder = 2x^4 - 3x^3 - 2x^2 + 3x.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 5
x = 1/2
Term by term:
  3x -> 3/2
  -2x^2 -> -1/2
  -3x^3 -> -3/8
  2x^4 -> 1/8
  sum = 3/4, using 14 multiplications
Horner's rule, from the top coefficient down:
  2, -2, -3, 3/2, 3/4
  P(1/2) = 3/4, using 4 multiplications

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 1
P = x^4 - 4x^3 + 6x^2 - 4x + 1
P = x^4 - 4x^3 + 6x^2 - 4x + 1, degree 4.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 7
Rational roots: 1 (4 times) - 2 candidates tried.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 1
P = x^2 - 2
P = x^2 - 2, degree 2.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 7
P has no rational roots (4 candidates tried).

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 1
P = x^3 - 2x
P = x^3 - 2x, degree 3.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 7
Rational roots: 0 - 4 candidates tried.
The other 2 root(s) of this degree 3 polynomial are not rational.

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 8
P = x^2 + 1
Q = x^3 - 2x

1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit
choice> 9
Bye.
```

## Modules

### main.eml

```eml
# P066 polynomials: type two polynomials, P and Q, the usual way - like
# 3x^2 - 2x + 1/2 - and add, subtract, multiply and divide them, evaluate
# P at a number both term by term and by Horner's rule, take its
# derivative, and find its rational roots. Every coefficient is an exact
# fraction.
import frac
import poly

def trim(s):
    0 => i
    len(s) => j
    while i < j and s[i] == " ":
        i + 1 => i
    while j > i and s[j - 1] == " ":
        j - 1 => j
    return s[i:j]

def ask_poly(name):
    while True:
        trim(input(name + " = ")) => s
        if s == "":
            return [False, []]
        poly.parse(s) => r
        if r[1] == "":
            return [True, r[0]]
        ("Cannot read that: " + r[1] + ". Write terms like 3x^2, -x, 1/2 or 0.25x^3.") ^0

def show(state):
    ("P = " + poly.text(state[0])) ^0
    ("Q = " + poly.text(state[1])) ^0

def set_poly(state, k):
    "PQ"[k] => name
    ask_poly(name) => r
    if not r[0]:
        "Cancelled." ^0
        return
    r[1] => state[k]
    (name + " = " + poly.text(r[1]) + ", degree " + degree_text(r[1]) + ".") ^0

def degree_text(p):
    if len(p) == 0:
        return "none (the zero polynomial)"
    return str(poly.degree(p))

def arithmetic(state):
    state[0] => p
    state[1] => q
    ("P + Q = " + poly.text(poly.plus(p, q))) ^0
    ("P - Q = " + poly.text(poly.minus(p, q))) ^0
    ("P * Q = " + poly.text(poly.times(p, q))) ^0

def divide(state):
    if len(state[1]) == 0:
        "Q is 0, and nothing divides by 0." ^0
        return
    poly.divided(state[0], state[1]) => d
    ("P / Q: quotient " + poly.text(d[0]) + ", remainder " + poly.text(d[1]) + ".") ^0
    poly.plus(poly.times(d[0], state[1]), d[1]) => back
    ("Check: quotient * Q + remainder = " + poly.text(back) + ".") ^0

def evaluate(state):
    while True:
        trim(input("x = ")) => s
        if s == "":
            "Cancelled." ^0
            return
        "" => body
        1 => sign
        s => body
        if len(s) > 1 and s[0] == "-":
            s[1:len(s)] => body
            -1 => sign
        frac.parse(body) => x
        if len(x) > 0:
            break
        "Type a number like 2, -3, 1/2 or 0.5." ^0
    [sign * x[0], x[1]] => x
    state[0] => p
    if len(p) == 0:
        "P is the zero polynomial: it is 0 everywhere." ^0
        return
    poly.direct(p, x) => d
    poly.horner(p, x) => h
    "Term by term:" ^0
    for k in [0:len(p) - 1]:
        if not frac.is_zero(p[k]):
            ("  " + poly.text(poly_term(p[k], k)) + " -> " + frac.text(d[1][k])) ^0
    ("  sum = " + frac.text(d[0]) + ", using " + str(d[2]) + " multiplications") ^0
    "Horner's rule, from the top coefficient down:" ^0
    "" => line
    for st in h[1]:
        if line != "":
            line + ", " => line
        line + frac.text(st[1]) => line
    ("  " + line) ^0
    ("  P(" + frac.text(x) + ") = " + frac.text(h[0]) + ", using " + str(poly.degree(p)) + " multiplications") ^0

def poly_term(c, k):
    [[0, 1]] * (k + 1) => t
    c => t[k]
    return t

def derive(state):
    poly.derivative(state[0]) => d
    ("P' = " + poly.text(d)) ^0

def roots(state):
    state[0] => p
    if len(p) == 0:
        "P is the zero polynomial: every number is a root." ^0
        return
    if len(p) == 1:
        "P is a nonzero constant: it has no roots." ^0
        return
    poly.rational_roots(p) => r
    if r[2] != "":
        ("Not tried: " + r[2] + ".") ^0
        return
    if len(r[0]) == 0:
        ("P has no rational roots (" + str(r[1]) + " candidates tried).") ^0
        return
    "" => s
    0 => count
    for x in r[0]:
        if s != "":
            s + ", " => s
        s + frac.text(x[0]) => s
        if x[1] == 2:
            s + " (twice)" => s
        elif x[1] > 2:
            s + " (" + str(x[1]) + " times)" => s
        count + x[1] => count
    ("Rational roots: " + s + " - " + str(r[1]) + " candidates tried.") ^0
    if count < poly.degree(p):
        ("The other " + str(poly.degree(p) - count) + " root(s) of this degree " + str(poly.degree(p)) + " polynomial are not rational.") ^0

"== Polynomials ==" ^0
"Write a polynomial the usual way, like 3x^2 - 2x + 1/2; every coefficient is an exact fraction." ^0
# [P, Q]
[[[-6, 1], [11, 1], [-6, 1], [1, 1]], [[-1, 1], [1, 1]]] => state
show(state)
True => running
while running:
    "" ^0
    "1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit" ^0
    trim(input("choice> ")) => choice
    if choice == "1":
        set_poly(state, 0)
    elif choice == "2":
        set_poly(state, 1)
    elif choice == "3":
        arithmetic(state)
    elif choice == "4":
        divide(state)
    elif choice == "5":
        evaluate(state)
    elif choice == "6":
        derive(state)
    elif choice == "7":
        roots(state)
    elif choice == "8":
        state[0] => t
        state[1] => state[0]
        t => state[1]
        show(state)
    elif choice == "9":
        False => running
    else:
        "Pick a number from 1 to 9." ^0
"Bye." ^0
```

Python projection of main.eml:

```python
import frac
import poly

def trim(s):
    i = 0
    j = len(s)
    while i < j and s[i] == " ":
        i = i + 1
    while j > i and s[j - 1] == " ":
        j = j - 1
    return s[i:j]

def ask_poly(name):
    while True:
        s = trim(input(name + " = "))
        if s == "":
            return [False, []]
        r = poly.parse(s)
        if r[1] == "":
            return [True, r[0]]
        print("Cannot read that: " + r[1] + ". Write terms like 3x^2, -x, 1/2 or 0.25x^3.")

def show(state):
    print("P = " + poly.text(state[0]))
    print("Q = " + poly.text(state[1]))

def set_poly(state, k):
    name = "PQ"[k]
    r = ask_poly(name)
    if not r[0]:
        print("Cancelled.")
        return
    state[k] = r[1]
    print(name + " = " + poly.text(r[1]) + ", degree " + degree_text(r[1]) + ".")

def degree_text(p):
    if len(p) == 0:
        return "none (the zero polynomial)"
    return str(poly.degree(p))

def arithmetic(state):
    p = state[0]
    q = state[1]
    print("P + Q = " + poly.text(poly.plus(p, q)))
    print("P - Q = " + poly.text(poly.minus(p, q)))
    print("P * Q = " + poly.text(poly.times(p, q)))

def divide(state):
    if len(state[1]) == 0:
        print("Q is 0, and nothing divides by 0.")
        return
    d = poly.divided(state[0], state[1])
    print("P / Q: quotient " + poly.text(d[0]) + ", remainder " + poly.text(d[1]) + ".")
    back = poly.plus(poly.times(d[0], state[1]), d[1])
    print("Check: quotient * Q + remainder = " + poly.text(back) + ".")

def evaluate(state):
    while True:
        s = trim(input("x = "))
        if s == "":
            print("Cancelled.")
            return
        body = ""
        sign = 1
        body = s
        if len(s) > 1 and s[0] == "-":
            body = s[1:len(s)]
            sign = -1
        x = frac.parse(body)
        if len(x) > 0:
            break
        print("Type a number like 2, -3, 1/2 or 0.5.")
    x = [sign * x[0], x[1]]
    p = state[0]
    if len(p) == 0:
        print("P is the zero polynomial: it is 0 everywhere.")
        return
    d = poly.direct(p, x)
    h = poly.horner(p, x)
    print("Term by term:")
    for k in range(0, len(p)):
        if not frac.is_zero(p[k]):
            print("  " + poly.text(poly_term(p[k], k)) + " -> " + frac.text(d[1][k]))
    print("  sum = " + frac.text(d[0]) + ", using " + str(d[2]) + " multiplications")
    print("Horner's rule, from the top coefficient down:")
    line = ""
    for st in h[1]:
        if line != "":
            line = line + ", "
        line = line + frac.text(st[1])
    print("  " + line)
    print("  P(" + frac.text(x) + ") = " + frac.text(h[0]) + ", using " + str(poly.degree(p)) + " multiplications")

def poly_term(c, k):
    t = [[0, 1]] * (k + 1)
    t[k] = c
    return t

def derive(state):
    d = poly.derivative(state[0])
    print("P' = " + poly.text(d))

def roots(state):
    p = state[0]
    if len(p) == 0:
        print("P is the zero polynomial: every number is a root.")
        return
    if len(p) == 1:
        print("P is a nonzero constant: it has no roots.")
        return
    r = poly.rational_roots(p)
    if r[2] != "":
        print("Not tried: " + r[2] + ".")
        return
    if len(r[0]) == 0:
        print("P has no rational roots (" + str(r[1]) + " candidates tried).")
        return
    s = ""
    count = 0
    for x in r[0]:
        if s != "":
            s = s + ", "
        s = s + frac.text(x[0])
        if x[1] == 2:
            s = s + " (twice)"
        elif x[1] > 2:
            s = s + " (" + str(x[1]) + " times)"
        count = count + x[1]
    print("Rational roots: " + s + " - " + str(r[1]) + " candidates tried.")
    if count < poly.degree(p):
        print("The other " + str(poly.degree(p) - count) + " root(s) of this degree " + str(poly.degree(p)) + " polynomial are not rational.")

print("== Polynomials ==")
print("Write a polynomial the usual way, like 3x^2 - 2x + 1/2; every coefficient is an exact fraction.")
state = [[[-6, 1], [11, 1], [-6, 1], [1, 1]], [[-1, 1], [1, 1]]]
show(state)
running = True
while running:
    print("")
    print("1) set P  2) set Q  3) P + Q, P - Q, P * Q  4) P / Q  5) P at x  6) P'  7) roots of P  8) swap P and Q  9) quit")
    choice = trim(input("choice> "))
    if choice == "1":
        set_poly(state, 0)
    elif choice == "2":
        set_poly(state, 1)
    elif choice == "3":
        arithmetic(state)
    elif choice == "4":
        divide(state)
    elif choice == "5":
        evaluate(state)
    elif choice == "6":
        derive(state)
    elif choice == "7":
        roots(state)
    elif choice == "8":
        t = state[0]
        state[0] = state[1]
        state[1] = t
        show(state)
    elif choice == "9":
        running = False
    else:
        print("Pick a number from 1 to 9.")
print("Bye.")
```

### poly.eml

```eml
# P066 polynomials - a polynomial is the list of its coefficients, exact
# fractions, from the constant up: [c0, c1, ..., cn] is c0 + c1 x + ... +
# cn x^n, with no zero at the end (the zero polynomial is []).
import frac

def trimmed(p):
    len(p) => n
    while n > 0 and frac.is_zero(p[n - 1]):
        n - 1 => n
    return p[0:n]

def degree(p):
    return len(p) - 1

def coefficient(p, k):
    if k < len(p):
        return p[k]
    return [0, 1]

def plus(p, q):
    [] => out
    len(p) => n
    if len(q) > n:
        len(q) => n
    for k in [0:n - 1]:
        out + [frac.plus(coefficient(p, k), coefficient(q, k))] => out
    return trimmed(out)

def minus(p, q):
    [] => out
    len(p) => n
    if len(q) > n:
        len(q) => n
    for k in [0:n - 1]:
        out + [frac.minus(coefficient(p, k), coefficient(q, k))] => out
    return trimmed(out)

def times(p, q):
    if len(p) == 0 or len(q) == 0:
        return []
    [[0, 1]] * (len(p) + len(q) - 1) => out
    for i in [0:len(p) - 1]:
        for j in [0:len(q) - 1]:
            frac.plus(out[i + j], frac.times(p[i], q[j])) => out[i + j]
    return trimmed(out)

def divided(p, q):
    # Long division: [quotient, remainder] with p = quotient * q + remainder
    # and the remainder of lower degree than q (q not zero).
    [] => quo
    p[0:len(p)] => r
    if len(r) >= len(q):
        [[0, 1]] * (len(r) - len(q) + 1) => quo
    while len(r) >= len(q) and len(r) > 0:
        len(r) - len(q) => shift
        frac.over(r[len(r) - 1], q[len(q) - 1]) => c
        c => quo[shift]
        for k in [0:len(q) - 1]:
            frac.minus(r[shift + k], frac.times(c, q[k])) => r[shift + k]
        trimmed(r[0:len(r) - 1]) => r
    return [trimmed(quo), trimmed(r)]

def derivative(p):
    [] => out
    for k in [1:len(p) - 1]:
        out + [frac.times([k, 1], p[k])] => out
    return trimmed(out)

def horner(p, x):
    # [value, steps]: start from the top coefficient, then again and again
    # multiply by x and add the next one down - n multiplications for degree
    # n. Each step is [coefficient added, value so far].
    [0, 1] => v
    [] => steps
    for k in [0:len(p) - 1]:
        len(p) - 1 - k => i
        frac.plus(frac.times(v, x), p[i]) => v
        steps + [[p[i], v]] => steps
    return [v, steps]

def direct(p, x):
    # The corpus case polynomial-evaluator's way: each term c x^k on its
    # own, x^k by k multiplications, then the terms added up. Returns
    # [value, each term's value, multiplications used].
    [0, 1] => total
    [] => parts
    0 => mults
    for k in [0:len(p) - 1]:
        [1, 1] => xk
        for i in [1:k]:
            frac.times(xk, x) => xk
            mults + 1 => mults
        frac.times(p[k], xk) => term
        if k > 0:
            mults + 1 => mults
        parts + [term] => parts
        frac.plus(total, term) => total
    return [total, parts, mults]

def text(p):
    # 3x^2 - 2x + 1/2; "0" for the zero polynomial.
    if len(p) == 0:
        return "0"
    "" => s
    for k in [0:len(p) - 1]:
        len(p) - 1 - k => i
        p[i] => c
        if not frac.is_zero(c):
            c[0] < 0 => negative
            frac.text([abs(c[0]), c[1]]) => mag
            if i > 0 and mag == "1":
                "" => mag
            elif i > 0 and c[1] != 1:
                "(" + mag + ")" => mag
            "" => xpart
            if i == 1:
                "x" => xpart
            elif i > 1:
                "x^" + str(i) => xpart
            if s == "":
                if negative:
                    "-" => s
            elif negative:
                s + " - " => s
            else:
                s + " + " => s
            s + mag + xpart => s
    return s

def digits_value(s):
    if s == "" or len(s) > 9:
        return -1
    0 => n
    for c in s:
        if not (c >= "0" and c <= "9"):
            return -1
        n * 10 + int(c) => n
    return n

def digit_like(ch):
    return (ch >= "0" and ch <= "9") or ch == "." or ch == "/"

def parse(s):
    # "3x^2 - 2x + 1/2", "x^3 - x", "0.5x + 2" as [polynomial, ""], or
    # [[], what is wrong]. Terms with the same power are added; powers up to
    # 20.
    # Spaces are dropped - but two numbers with only spaces between them are
    # an error, not one number: "3 4" is not 34.
    "" => t
    "" => before
    False => gap
    for ch in s:
        if ch == " ":
            if before != "":
                True => gap
        else:
            if gap and digit_like(before) and digit_like(ch):
                return [[], "two numbers side by side"]
            False => gap
            ch => before
            t + ch => t
    if t == "":
        return [[], "it is empty"]
    [[0, 1]] * 21 => coeffs
    0 => i
    while i < len(t):
        [1, 1] => sign
        if t[i] == "+" or t[i] == "-":
            if t[i] == "-":
                [-1, 1] => sign
            i + 1 => i
        elif i > 0:
            return [[], "a term without + or - before it"]
        i => j
        while j < len(t) and t[j] != "x" and t[j] != "+" and t[j] != "-":
            j + 1 => j
        t[i:j] => num
        # a fraction written in brackets, as text() writes one: (9/2)x
        if len(num) > 2 and num[0] == "(" and num[len(num) - 1] == ")":
            num[1:len(num) - 1] => num
        [1, 1] => c
        if num != "":
            frac.parse(num) => c
            if len(c) == 0:
                return [[], "the coefficient " + num]
        j => i
        0 => power
        if i < len(t) and t[i] == "x":
            1 => power
            i + 1 => i
            if i < len(t) and t[i] == "^":
                i + 1 => i
                i => j
                while j < len(t) and t[j] >= "0" and t[j] <= "9":
                    j + 1 => j
                digits_value(t[i:j]) => power
                if power < 0:
                    return [[], "a power after ^ that is not a whole number"]
                j => i
        elif num == "":
            return [[], "a sign with no term after it"]
        if power > 20:
            return [[], "a power above 20"]
        frac.plus(coeffs[power], frac.times(sign, c)) => coeffs[power]
    return [trimmed(coeffs), ""]

def divisors(n):
    # The positive divisors of n > 0, by trial up to its square root.
    [] => small
    [] => large
    1 => d
    while d * d <= n:
        if n % d == 0:
            small + [d] => small
            if d * d != n:
                [frac.quotient(n, d)] + large => large
        d + 1 => d
    return small + large

def rational_roots(p):
    # The rational roots of p with their multiplicities, by the rational
    # root theorem: after scaling p to whole coefficients, a root a/b in
    # lowest terms has a dividing the lowest nonzero coefficient and b the
    # leading one. Every such candidate is tried. Returns [[root, times],
    # ...] in increasing order, the candidates tried, and "" or why not.
    # The smallest whole scale that clears every denominator: their least
    # common multiple, by exact whole-number division (a float would lose
    # digits once the denominators multiply up).
    1 => scale
    for c in p:
        frac.quotient(scale * c[1], frac.gcd(scale, c[1])) => scale
    [] => a
    for c in p:
        a + [c[0] * frac.quotient(scale, c[1])] => a
    0 => zeros
    while zeros < len(a) and a[zeros] == 0:
        zeros + 1 => zeros
    [] => roots
    if zeros > 0:
        roots + [[[0, 1], zeros]] => roots
    a[zeros:len(a)] => a
    if abs(a[0]) > 1000000 or abs(a[len(a) - 1]) > 1000000:
        return [roots, 0, "the coefficients are too large to try every candidate"]
    [] => cands
    for num in divisors(abs(a[0])):
        for den in divisors(abs(a[len(a) - 1])):
            for sg in [1, -1]:
                frac.make(sg * num, den) => r
                False => seen
                for x in cands:
                    if x[0] == r[0] and x[1] == r[1]:
                        True => seen
                if not seen:
                    cands + [r] => cands
    p[zeros:len(p)] => rest
    for r in sorted_fracs(cands):
        0 => times
        True => going
        while going and len(rest) > 1:
            divided(rest, [frac.minus([0, 1], r), [1, 1]]) => d
            if len(d[1]) == 0:
                d[0] => rest
                times + 1 => times
            else:
                False => going
        if times > 0:
            roots + [[r, times]] => roots
    return [sorted_roots(roots), len(cands), ""]

def sorted_fracs(xs):
    xs[0:len(xs)] => out
    for i in [1:len(out) - 1]:
        out[i] => x
        i => j
        while j > 0 and frac.less(x, out[j - 1]):
            out[j - 1] => out[j]
            j - 1 => j
        x => out[j]
    return out

def sorted_roots(rs):
    rs[0:len(rs)] => out
    for i in [1:len(out) - 1]:
        out[i] => x
        i => j
        while j > 0 and frac.less(x[0], out[j - 1][0]):
            out[j - 1] => out[j]
            j - 1 => j
        x => out[j]
    return out
```

Python projection of poly.eml:

```python
import frac

def trimmed(p):
    n = len(p)
    while n > 0 and frac.is_zero(p[n - 1]):
        n = n - 1
    return p[0:n]

def degree(p):
    return len(p) - 1

def coefficient(p, k):
    if k < len(p):
        return p[k]
    return [0, 1]

def plus(p, q):
    out = []
    n = len(p)
    if len(q) > n:
        n = len(q)
    for k in range(0, n):
        out = out + [frac.plus(coefficient(p, k), coefficient(q, k))]
    return trimmed(out)

def minus(p, q):
    out = []
    n = len(p)
    if len(q) > n:
        n = len(q)
    for k in range(0, n):
        out = out + [frac.minus(coefficient(p, k), coefficient(q, k))]
    return trimmed(out)

def times(p, q):
    if len(p) == 0 or len(q) == 0:
        return []
    out = [[0, 1]] * (len(p) + len(q) - 1)
    for i in range(0, len(p)):
        for j in range(0, len(q)):
            out[i + j] = frac.plus(out[i + j], frac.times(p[i], q[j]))
    return trimmed(out)

def divided(p, q):
    quo = []
    r = p[0:len(p)]
    if len(r) >= len(q):
        quo = [[0, 1]] * (len(r) - len(q) + 1)
    while len(r) >= len(q) and len(r) > 0:
        shift = len(r) - len(q)
        c = frac.over(r[len(r) - 1], q[len(q) - 1])
        quo[shift] = c
        for k in range(0, len(q)):
            r[shift + k] = frac.minus(r[shift + k], frac.times(c, q[k]))
        r = trimmed(r[0:len(r) - 1])
    return [trimmed(quo), trimmed(r)]

def derivative(p):
    out = []
    for k in range(1, len(p)):
        out = out + [frac.times([k, 1], p[k])]
    return trimmed(out)

def horner(p, x):
    v = [0, 1]
    steps = []
    for k in range(0, len(p)):
        i = len(p) - 1 - k
        v = frac.plus(frac.times(v, x), p[i])
        steps = steps + [[p[i], v]]
    return [v, steps]

def direct(p, x):
    total = [0, 1]
    parts = []
    mults = 0
    for k in range(0, len(p)):
        xk = [1, 1]
        for i in range(1, k+1):
            xk = frac.times(xk, x)
            mults = mults + 1
        term = frac.times(p[k], xk)
        if k > 0:
            mults = mults + 1
        parts = parts + [term]
        total = frac.plus(total, term)
    return [total, parts, mults]

def text(p):
    if len(p) == 0:
        return "0"
    s = ""
    for k in range(0, len(p)):
        i = len(p) - 1 - k
        c = p[i]
        if not frac.is_zero(c):
            negative = c[0] < 0
            mag = frac.text([abs(c[0]), c[1]])
            if i > 0 and mag == "1":
                mag = ""
            elif i > 0 and c[1] != 1:
                mag = "(" + mag + ")"
            xpart = ""
            if i == 1:
                xpart = "x"
            elif i > 1:
                xpart = "x^" + str(i)
            if s == "":
                if negative:
                    s = "-"
            elif negative:
                s = s + " - "
            else:
                s = s + " + "
            s = s + mag + xpart
    return s

def digits_value(s):
    if s == "" or len(s) > 9:
        return -1
    n = 0
    for c in s:
        if not (c >= "0" and c <= "9"):
            return -1
        n = n * 10 + int(c)
    return n

def digit_like(ch):
    return ch >= "0" and ch <= "9" or ch == "." or ch == "/"

def parse(s):
    t = ""
    before = ""
    gap = False
    for ch in s:
        if ch == " ":
            if before != "":
                gap = True
        else:
            if gap and digit_like(before) and digit_like(ch):
                return [[], "two numbers side by side"]
            gap = False
            before = ch
            t = t + ch
    if t == "":
        return [[], "it is empty"]
    coeffs = [[0, 1]] * 21
    i = 0
    while i < len(t):
        sign = [1, 1]
        if t[i] == "+" or t[i] == "-":
            if t[i] == "-":
                sign = [-1, 1]
            i = i + 1
        elif i > 0:
            return [[], "a term without + or - before it"]
        j = i
        while j < len(t) and t[j] != "x" and t[j] != "+" and t[j] != "-":
            j = j + 1
        num = t[i:j]
        if len(num) > 2 and num[0] == "(" and num[len(num) - 1] == ")":
            num = num[1:len(num) - 1]
        c = [1, 1]
        if num != "":
            c = frac.parse(num)
            if len(c) == 0:
                return [[], "the coefficient " + num]
        i = j
        power = 0
        if i < len(t) and t[i] == "x":
            power = 1
            i = i + 1
            if i < len(t) and t[i] == "^":
                i = i + 1
                j = i
                while j < len(t) and t[j] >= "0" and t[j] <= "9":
                    j = j + 1
                power = digits_value(t[i:j])
                if power < 0:
                    return [[], "a power after ^ that is not a whole number"]
                i = j
        elif num == "":
            return [[], "a sign with no term after it"]
        if power > 20:
            return [[], "a power above 20"]
        coeffs[power] = frac.plus(coeffs[power], frac.times(sign, c))
    return [trimmed(coeffs), ""]

def divisors(n):
    small = []
    large = []
    d = 1
    while d * d <= n:
        if n % d == 0:
            small = small + [d]
            if d * d != n:
                large = [frac.quotient(n, d)] + large
        d = d + 1
    return small + large

def rational_roots(p):
    scale = 1
    for c in p:
        scale = frac.quotient(scale * c[1], frac.gcd(scale, c[1]))
    a = []
    for c in p:
        a = a + [c[0] * frac.quotient(scale, c[1])]
    zeros = 0
    while zeros < len(a) and a[zeros] == 0:
        zeros = zeros + 1
    roots = []
    if zeros > 0:
        roots = roots + [[[0, 1], zeros]]
    a = a[zeros:len(a)]
    if abs(a[0]) > 1000000 or abs(a[len(a) - 1]) > 1000000:
        return [roots, 0, "the coefficients are too large to try every candidate"]
    cands = []
    for num in divisors(abs(a[0])):
        for den in divisors(abs(a[len(a) - 1])):
            for sg in [1, -1]:
                r = frac.make(sg * num, den)
                seen = False
                for x in cands:
                    if x[0] == r[0] and x[1] == r[1]:
                        seen = True
                if not seen:
                    cands = cands + [r]
    rest = p[zeros:len(p)]
    for r in sorted_fracs(cands):
        times = 0
        going = True
        while going and len(rest) > 1:
            d = divided(rest, [frac.minus([0, 1], r), [1, 1]])
            if len(d[1]) == 0:
                rest = d[0]
                times = times + 1
            else:
                going = False
        if times > 0:
            roots = roots + [[r, times]]
    return [sorted_roots(roots), len(cands), ""]

def sorted_fracs(xs):
    out = xs[0:len(xs)]
    for i in range(1, len(out)):
        x = out[i]
        j = i
        while j > 0 and frac.less(x, out[j - 1]):
            out[j] = out[j - 1]
            j = j - 1
        out[j] = x
    return out

def sorted_roots(rs):
    out = rs[0:len(rs)]
    for i in range(1, len(out)):
        x = out[i]
        j = i
        while j > 0 and frac.less(x[0], out[j - 1][0]):
            out[j] = out[j - 1]
            j = j - 1
        out[j] = x
    return out
```

### frac.eml

```eml
# P066 polynomials - exact fractions, from P007 (calculator) by way of
# P038, P044, P056, P060 and P062.
# A number is a list [n, d]: n / d in lowest terms with d > 0. Integers have
# no size limit, but EML has no //, and a / b goes through a float that cannot
# hold a large quotient exactly, so whole-number division is written out.

def quotient(a, b):
    # a // b for whole numbers a >= 0 and b > 0, exact at any size. Long
    # division by doubling: take away the largest b * 2^k that still fits.
    0 => q
    while a >= b:
        b => m
        1 => k
        while m + m <= a:
            m + m => m
            k + k => k
        a - m => a
        q + k => q
    return q

def gcd(a, b):
    while b != 0:
        a % b => r
        b => a
        r => b
    return a

def make(n, d):
    # n / d in lowest terms with a positive denominator (d != 0).
    if d < 0:
        0 - n => n
        0 - d => d
    abs(n) => a
    gcd(a, d) => g
    if n < 0:
        return [0 - quotient(a, g), quotient(d, g)]
    return [quotient(a, g), quotient(d, g)]

def plus(x, y):
    return make(x[0] * y[1] + y[0] * x[1], x[1] * y[1])

def minus(x, y):
    return make(x[0] * y[1] - y[0] * x[1], x[1] * y[1])

def times(x, y):
    return make(x[0] * y[0], x[1] * y[1])

def over(x, y):
    # x / y; the caller has checked that y is not zero.
    return make(x[0] * y[1], x[1] * y[0])

def is_zero(x):
    return x[0] == 0

def digits_value(s):
    # The value of a string of 1 to 9 digits, otherwise -1.
    if s == "" or len(s) > 9:
        return -1
    0 => n
    for c in s:
        if not (c in "0123456789"):
            return -1
        n * 10 + int(c) => n
    return n

def parse(s):
    # "3", "-2", "3/4", "-0.25" as a fraction; [] if it is not a number.
    1 => sign
    if len(s) > 0 and s[0] == "-":
        -1 => sign
        s[1:len(s)] => s
    0 => k
    while k < len(s) and s[k] != "/" and s[k] != ".":
        k + 1 => k
    if k == len(s):
        digits_value(s) => n
        if n < 0:
            return []
        return make(sign * n, 1)
    digits_value(s[0:k]) => whole
    s[k + 1:len(s)] => rest
    digits_value(rest) => part
    if whole < 0 or part < 0:
        return []
    if s[k] == "/":
        if part == 0:
            return []
        return make(sign * whole, part)
    # a decimal point: 0.25 is 25 / 100
    1 => scale
    for c in rest:
        scale * 10 => scale
    return make(sign * (whole * scale + part), scale)

def less(x, y):
    return x[0] * y[1] < y[0] * x[1]

def decimal(x, places):
    # x to the given number of decimal places, halves rounded away from zero,
    # with "about " in front when the value does not end there exactly.
    "" => sign
    abs(x[0]) => n
    if x[0] < 0:
        "-" => sign
    1 => scale
    for k in [1:places]:
        scale * 10 => scale
    quotient(2 * n * scale + x[1], 2 * x[1]) => t
    "" => about
    if (n * scale) % x[1] != 0:
        "about " => about
    str(quotient(t, scale)) => whole
    str(t % scale) => part
    while len(part) < places:
        "0" + part => part
    # drop trailing zeros of an exact value
    if about == "":
        while len(part) > 0 and part[len(part) - 1] == "0":
            part[0:len(part) - 1] => part
    if t == 0:
        "" => sign
    if part == "":
        return about + sign + whole
    return about + sign + whole + "." + part

def text(x):
    # "3", "-3", "3/4", "-3/4".
    if x[1] == 1:
        return str(x[0])
    return str(x[0]) + "/" + str(x[1])
```

Python projection of frac.eml:

```python
def quotient(a, b):
    q = 0
    while a >= b:
        m = b
        k = 1
        while m + m <= a:
            m = m + m
            k = k + k
        a = a - m
        q = q + k
    return q

def gcd(a, b):
    while b != 0:
        r = a % b
        a = b
        b = r
    return a

def make(n, d):
    if d < 0:
        n = 0 - n
        d = 0 - d
    a = abs(n)
    g = gcd(a, d)
    if n < 0:
        return [0 - quotient(a, g), quotient(d, g)]
    return [quotient(a, g), quotient(d, g)]

def plus(x, y):
    return make(x[0] * y[1] + y[0] * x[1], x[1] * y[1])

def minus(x, y):
    return make(x[0] * y[1] - y[0] * x[1], x[1] * y[1])

def times(x, y):
    return make(x[0] * y[0], x[1] * y[1])

def over(x, y):
    return make(x[0] * y[1], x[1] * y[0])

def is_zero(x):
    return x[0] == 0

def digits_value(s):
    if s == "" or len(s) > 9:
        return -1
    n = 0
    for c in s:
        if not c in "0123456789":
            return -1
        n = n * 10 + int(c)
    return n

def parse(s):
    sign = 1
    if len(s) > 0 and s[0] == "-":
        sign = -1
        s = s[1:len(s)]
    k = 0
    while k < len(s) and s[k] != "/" and s[k] != ".":
        k = k + 1
    if k == len(s):
        n = digits_value(s)
        if n < 0:
            return []
        return make(sign * n, 1)
    whole = digits_value(s[0:k])
    rest = s[k + 1:len(s)]
    part = digits_value(rest)
    if whole < 0 or part < 0:
        return []
    if s[k] == "/":
        if part == 0:
            return []
        return make(sign * whole, part)
    scale = 1
    for c in rest:
        scale = scale * 10
    return make(sign * (whole * scale + part), scale)

def less(x, y):
    return x[0] * y[1] < y[0] * x[1]

def decimal(x, places):
    sign = ""
    n = abs(x[0])
    if x[0] < 0:
        sign = "-"
    scale = 1
    for k in range(1, places+1):
        scale = scale * 10
    t = quotient(2 * n * scale + x[1], 2 * x[1])
    about = ""
    if n * scale % x[1] != 0:
        about = "about "
    whole = str(quotient(t, scale))
    part = str(t % scale)
    while len(part) < places:
        part = "0" + part
    if about == "":
        while len(part) > 0 and part[len(part) - 1] == "0":
            part = part[0:len(part) - 1]
    if t == 0:
        sign = ""
    if part == "":
        return about + sign + whole
    return about + sign + whole + "." + part

def text(x):
    if x[1] == 1:
        return str(x[0])
    return str(x[0]) + "/" + str(x[1])
```

## README

# P066 - Polynomials

Two polynomials, P and Q, typed the usual way, like 3x^2 - 2x + 1/2, with
exact fraction coefficients. You can add, subtract and multiply them, and
divide P by Q, getting a quotient and remainder that are checked back. P
can be evaluated at a number two ways: term by term, and by Horner's rule,
with the multiplications counted. It can also be differentiated, and its
rational roots are all found.

- `main.eml` - the menu, reading polynomials and numbers, and the messages
- `poly.eml` - the arithmetic, division, both ways of evaluating, the
  derivative, reading and writing polynomials, and the rational roots
- `frac.eml` - exact fractions, from P007 by way of P038, P044, P056, P060
  and P062

How each part works:

- A polynomial is the list of its coefficients from the constant up, so
  3x^2 - 2x + 1/2 is [1/2, -2, 3]. Adding and subtracting go coefficient by
  coefficient, and multiplying adds each product a x^i times b x^j into
  the place for x^(i+j).
- Dividing is long division: take the leading term of what is left over
  the leading term of Q, take that many Qs away, and repeat until what is
  left has a lower degree than Q.
- Term by term, as in the corpus case `polynomial-evaluator`, each term is
  worked out on its own and then they are added, which shows what every
  term gives. Horner's rule rewrites 3x^2 - 2x + 1/2 as (3x - 2)x + 1/2: it
  starts from the top coefficient and then multiplies by x and adds the
  next one, n multiplications for degree n where term by term needs about
  n^2 / 2. Both give the same value.
- Rational roots: scale P to whole coefficients and take out any factor x.
  By the rational root theorem, a root a/b in lowest terms has a dividing
  the lowest coefficient and b dividing the highest. Every such candidate
  is tried, and each root found is divided out as many times as it goes, so
  its multiplicity is known. Whatever degree is left over belongs to roots
  that are not rational.

What is checked: menu choices 1 to 9. Polynomials are written as terms in
x with whole, decimal or fraction coefficients (a fraction may be in
brackets, as the program writes it), powers up to 20, and no two numbers
side by side. A number x can be whole, decimal, a fraction, or negative.
Root search needs the scaled end coefficients to be at most a million. An
empty answer cancels.

Sessions: `sessions/basic.in` starts from P = x^3 - 6x^2 + 11x - 6, which
is (x-1)(x-2)(x-3), and Q = x - 1.
- P + Q, P - Q and P * Q.
- P / Q = x^2 - 5x + 6 with remainder 0.
- P at 4: 6, from 9 multiplications term by term and 3 by Horner.
- P' and the roots 1, 2, 3.
- 2x^4 - 3x^3 - 2x^2 + 3x: roots -1, 0, 1 and 3/2. Divided by x^2 + 1, and
  P(1/2) = 3/4.
- (x - 1)^4 written out: 1 four times.
- x^2 - 2: no rational root.
- x^3 - 2x: 0, with two roots that are not rational.
- P and Q swapped.

`sessions/bad-input.in` gives:
- menu choices 0 and x;
- the inputs "3 4", "x^", 2y, x^21, 1/0 and x--1;
- Q = 0, then a division by it;
- P = 5, which has no roots, and P = 0, which has every number as a root;
- x = abc, then -1/2 for the zero polynomial;
- 1000001x + 1000003, too large to search.

Built on the verified corpus case `polynomial-evaluator` (a polynomial
worked out term by term, each term's value shown).