Project 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.
3 modules · 2 recorded sessionstext-menu UI in the terminalupdated 2026-10-14
Every screen below was recorded under CPython. When this page was built, the EML interpreter replayed each session from the same input and printed the same bytes.
About
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 messagespoly.eml- the arithmetic, division, both ways of evaluating, the derivative, reading and writing polynomials, and the rational rootsfrac.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).
Recorded sessions
What the screen shows while someone uses the program. Each typed line appears after its prompt, the way a terminal shows it.
bad-input
interpreter: byte-equal== 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.
What was typed (28 lines)
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
basic
interpreter: byte-equal== 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.
What was typed (25 lines)
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
Modules
The program as written, entry module first. Each module transpiles to its own Python file, which is what eml project run executes.
main.eml(entry)
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 (main.py)
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 (poly.py)
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 (frac.py)
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])