Project P114
Quadratic solver
ax^2 + bx + c = 0 with exact rational, irrational and complex roots, square roots and divisions worked digit by digit on whole numbers; and the small root of x^2 + bx + 1 in floats, by the textbook formula and a stable one, against its exact digits.
2 modules · 2 recorded sessionstext-menu UI in the terminalupdated 2026-10-16
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
ax² + bx + c = 0, solved exactly: rational, irrational and complex roots. The program also solves x² + bx + 1 in floats, by the textbook formula and by a stable one, against the root's exact digits.
main.eml- the menu, the equation and its roots, and the float comparisonquad.eml- reading coefficients, square roots and divisions digit by digit, the solver, a float square root, and the small root two ways and exactly
How each part works:
- Coefficients have at most three decimals. They are scaled by 1000 to whole numbers, which leaves the roots as they are, so the discriminant b² - 4ac is exact.
- The kind of root, as in the corpus case
quadratic-solver, here with complex roots too: - a discriminant of 0 gives one double root; - the square of a fraction gives two rational roots; - any other positive one gives two irrational roots; - a negative one gives complex roots.
A leading coefficient below zero is turned positive by changing every sign, which leaves the roots as they are. - Square roots use the digit-pair method and divisions use long division, both on whole numbers. Only adding, subtracting, multiplying and comparing are used, so the 12 decimal places are exact: 18 are worked out, and then rounded half away from zero. - Floats, as in the corpus case catastrophic-cancellation: for x² + bx + 1 with b large, the textbook formula takes -b + √(b² - 4) when the two are nearly equal, and the digits they share cancel. The stable formula works out the big root, where nothing cancels, and divides: the two roots multiply to 1. The float square root is Newton's method from above, so it is the same in every Python. The program measures each printed float's relative error against 20 exact digits.
What is checked: menu choices 1 to 3; coefficients with at most 9 digits before the point and 3 after; b from 3 to 10¹² for the float comparison. An empty answer cancels.
Sessions: sessions/basic.in gives: - x² - 3x + 2 and 2x² - 7x + 3, with rational roots; - x² - 2, with irrational roots, and x² + 2x + 5, with complex roots; - x² - 2x + 1, a double root; - 2x - 4, 0 = 0 and 5 = 0, where a, or everything, is 0; - -x² + 2, and 0.5x² - 0.125x + 0.007; - the float comparison for b = 3, 10⁸ (the textbook root off by 49%) and 10¹² (off by 100%: it comes out as 0).
sessions/bad-input.in gives: - menu choices 0 and x; - coefficients: a cancel, abc, four decimals, two later cancels, a comma and ten digits; - the float comparison with a cancel, 2, 10¹² + 1 and x.
Built on the verified corpus cases quadratic-solver (the quadratic formula and its discriminant) and catastrophic-cancellation (the same algebra, and one arrangement keeps its digits).
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== Quadratic solver ==
ax^2 + bx + c = 0, solved exactly - and the same algebra in floats losing its small root.
1) solve 2) floats and cancellation 3) quit
choice> 0
Pick a number from 1 to 3.
1) solve 2) floats and cancellation 3) quit
choice> x
Pick a number from 1 to 3.
1) solve 2) floats and cancellation 3) quit
choice> 1
a>
Cancelled.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> abc
Type a number with at most 9 digits before the point and 3 after, such as -2 or 0.125.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 1.2345
Type a number with at most 9 digits before the point and 3 after, such as -2 or 0.125.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 1
b>
Cancelled.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 1
b> 2
c>
Cancelled.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 1
b> 2
c> 1,5
Type a number with at most 9 digits before the point and 3 after, such as -2 or 0.125.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 9999999999
Type a number with at most 9 digits before the point and 3 after, such as -2 or 0.125.
1) solve 2) floats and cancellation 3) quit
choice> 2
b, a whole number from 3 to 1000000000000>
Cancelled.
1) solve 2) floats and cancellation 3) quit
choice> 2
b, a whole number from 3 to 1000000000000> 2
Type a whole number from 3 to 1000000000000.
1) solve 2) floats and cancellation 3) quit
choice> 2
b, a whole number from 3 to 1000000000000> 1000000000001
Type a whole number from 3 to 1000000000000.
1) solve 2) floats and cancellation 3) quit
choice> 2
b, a whole number from 3 to 1000000000000> x
Type a whole number from 3 to 1000000000000.
1) solve 2) floats and cancellation 3) quit
choice> 3
Bye.
What was typed (30 lines)
0
x
1
1
abc
1
1.2345
1
1
1
1
2
1
1
2
1,5
1
9999999999
2
2
2
2
1000000000001
2
x
3
basic
interpreter: byte-equal== Quadratic solver ==
ax^2 + bx + c = 0, solved exactly - and the same algebra in floats losing its small root.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 1
b> -3
c> 2
x^2 - 3x + 2 = 0
Discriminant b^2 - 4ac = 1.
It is the square of a fraction: two rational roots, x = 1 and x = 2.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 2
b> -7
c> 3
2x^2 - 7x + 3 = 0
Discriminant b^2 - 4ac = 25.
It is the square of a fraction: two rational roots, x = 1/2 and x = 3.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 1
b> 0
c> -2
x^2 - 2 = 0
Discriminant b^2 - 4ac = 8.
It is positive and not the square of a fraction: two irrational roots, x = (-b - sqrt(D)) / 2a and (-b + sqrt(D)) / 2a:
x = -1.414213562373
x = 1.414213562373
to 12 places, each digit worked out exactly.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 1
b> 2
c> 5
x^2 + 2x + 5 = 0
Discriminant b^2 - 4ac = -16.
It is negative: two complex roots, x = (-b / 2a) - i (sqrt(-D) / 2a) and its conjugate:
x = -1 - 2.000000000000i
x = -1 + 2.000000000000i
the imaginary part to 12 places.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 1
b> -2
c> 1
x^2 - 2x + 1 = 0
Discriminant b^2 - 4ac = 0.
It is zero: one double root, x = 1.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 0
b> 2
c> -4
2x - 4 = 0
With a = 0 it is not quadratic: one root, x = 2.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 0
b> 0
c> 0
0 = 0
Every x solves it: all three coefficients are zero.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 0
b> 0
c> 5
5 = 0
No x solves it: it says a number other than zero equals zero.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> -1
b> 0
c> 2
-x^2 + 2 = 0
Discriminant b^2 - 4ac = 8.
It is positive and not the square of a fraction: two irrational roots, x = (-b - sqrt(D)) / 2a and (-b + sqrt(D)) / 2a:
x = -1.414213562373
x = 1.414213562373
to 12 places, each digit worked out exactly.
1) solve 2) floats and cancellation 3) quit
choice> 1
a> 0.5
b> -0.125
c> 0.007
0.5x^2 - 0.125x + 0.007 = 0
Discriminant b^2 - 4ac = 0.001625.
It is positive and not the square of a fraction: two irrational roots, x = (-b - sqrt(D)) / 2a and (-b + sqrt(D)) / 2a:
x = 0.084688711259
x = 0.165311288741
to 12 places, each digit worked out exactly.
1) solve 2) floats and cancellation 3) quit
choice> 2
b, a whole number from 3 to 1000000000000> 3
x^2 + 3x + 1 = 0, its root nearest zero:
exact, 20 digits: -3.8196601125010515179e-1
textbook formula: -0.3819660112501051 (right to about 15 significant digits)
stable formula: -0.38196601125010515 (right to 17 significant digits or more)
The textbook formula takes -b + sqrt(b^2 - 4) when the two are nearly equal, and the digits they share cancel. The stable one works out the big root, where nothing cancels, and divides: the product of the roots is c / a = 1.
1) solve 2) floats and cancellation 3) quit
choice> 2
b, a whole number from 3 to 1000000000000> 100000000
x^2 + 100000000x + 1 = 0, its root nearest zero:
exact, 20 digits: -1.0000000000000001000e-8
textbook formula: -1.4901161193847656e-08 (off by 49%)
stable formula: -1.0000000000000002e-08 (right to about 16 significant digits)
The textbook formula takes -b + sqrt(b^2 - 4) when the two are nearly equal, and the digits they share cancel. The stable one works out the big root, where nothing cancels, and divides: the product of the roots is c / a = 1.
1) solve 2) floats and cancellation 3) quit
choice> 2
b, a whole number from 3 to 1000000000000> 1000000000000
x^2 + 1000000000000x + 1 = 0, its root nearest zero:
exact, 20 digits: -1.0000000000000000000e-12
textbook formula: 0.0 (off by 100%)
stable formula: -1e-12 (right to 17 significant digits or more)
The textbook formula takes -b + sqrt(b^2 - 4) when the two are nearly equal, and the digits they share cancel. The stable one works out the big root, where nothing cancels, and divides: the product of the roots is c / a = 1.
1) solve 2) floats and cancellation 3) quit
choice> 3
Bye.
What was typed (47 lines)
1
1
-3
2
1
2
-7
3
1
1
0
-2
1
1
2
5
1
1
-2
1
1
0
2
-4
1
0
0
0
1
0
0
5
1
-1
0
2
1
0.5
-0.125
0.007
2
3
2
100000000
2
1000000000000
3
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# P114 quadratic solver: exact roots of ax^2 + bx + c = 0, real or complex,
# and the small root of x^2 + bx + 1 in floats, by the textbook formula and
# the stable one, against its exact digits.
import quad
def plain(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 number(s, most):
if len(s) == 0 or len(s) > most:
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 ask_coef(name):
# [coefficient x 1000, text], [] to cancel, ["!"] after a complaint
plain(input(name + "> ")) => t
if t == "":
return []
quad.parse_coef(t) => c
if len(c) == 0:
"Type a number with at most 9 digits before the point and 3 after, such as -2 or 0.125." ^0
return ["!"]
return [c[0], quad.trimmed(quad.scaled_text(c[0], 3))]
def equation(a, b, c):
"" => out
for part in [[a, "x^2"], [b, "x"], [c, ""]]:
part[0] => v
if v[0] != 0:
v[1] => t
"+" => sign
if t[0] == "-":
"-" => sign
t[1:len(t)] => t
if (t == "1") and part[1] != "":
"" => t
if out == "":
if sign == "-":
"-" => out
out + t + part[1] => out
else:
out + " " + sign + " " + t + part[1] => out
if out == "":
"0" => out
return out + " = 0"
def solve_menu():
ask_coef("a") => a
if len(a) == 0:
"Cancelled." ^0
return
if a[0] == "!":
return
ask_coef("b") => b
if len(b) == 0:
"Cancelled." ^0
return
if b[0] == "!":
return
ask_coef("c") => c
if len(c) == 0:
"Cancelled." ^0
return
if c[0] == "!":
return
equation(a, b, c) ^0
quad.solve(a[0], b[0], c[0]) => r
if r[0] == "all":
"Every x solves it: all three coefficients are zero." ^0
return
if r[0] == "none":
"No x solves it: it says a number other than zero equals zero." ^0
return
if r[0] == "linear":
("With a = 0 it is not quadratic: one root, x = " + r[1] + ".") ^0
return
b[0] * b[0] - 4 * a[0] * c[0] => D
("Discriminant b^2 - 4ac = " + quad.trimmed(quad.scaled_text(D, 6)) + ".") ^0
if r[0] == "double":
("It is zero: one double root, x = " + r[1] + ".") ^0
elif r[0] == "two rational":
("It is the square of a fraction: two rational roots, x = " + r[1] + " and x = " + r[2] + ".") ^0
elif r[0] == "two real":
"It is positive and not the square of a fraction: two irrational roots, x = (-b - sqrt(D)) / 2a and (-b + sqrt(D)) / 2a:" ^0
(" x = " + r[2]) ^0
(" x = " + r[3]) ^0
"to 12 places, each digit worked out exactly." ^0
else:
"It is negative: two complex roots, x = (-b / 2a) - i (sqrt(-D) / 2a) and its conjugate:" ^0
(" x = " + r[1] + " - " + r[2] + "i") ^0
(" x = " + r[1] + " + " + r[2] + "i") ^0
"the imaginary part to 12 places." ^0
def printed_value(x):
# [M, e] with the printed float equal to M x 10^e exactly, M >= 0
str(x) => t
if t[0] == "-":
t[1:len(t)] => t
0 => e
-1 => at
for i in [0:len(t) - 1]:
if at == -1 and t[i] == "e":
i => at
if at >= 0:
t[at + 1:len(t)] => es
1 => sign
if es[0] == "-":
0 - 1 => sign
es[1:len(es)] => es
elif es[0] == "+":
es[1:len(es)] => es
sign * int(es) => e
t[0:at] => t
"" => ds
0 => after
False => dot
for c in t:
if c == ".":
True => dot
else:
ds + c => ds
if dot:
after + 1 => after
return [int(ds), e - after]
def ten_to(k):
"1" => t
for i in [1:k]:
t + "0" => t
return int(t)
def quality(x, ex):
# how close the printed float x is to the exact root, -digits x
# 10^(lead - 19): [k, -1] with k the largest up to 17 whose relative
# error is at most 10^-k, or [0, whole percent] when it is over 10%
printed_value(x) => p
int(ex[0]) => m
ex[1] - 19 => e
p[1] => lo
if e < lo:
e => lo
p[0] * ten_to(p[1] - lo) => xv
if x < 0.0:
0 - xv => xv
0 - m * ten_to(e - lo) => ev
xv - ev => diff
if diff < 0:
0 - diff => diff
0 - ev => size
0 => k
while k < 17 and diff * ten_to(k + 1) <= size:
k + 1 => k
if k >= 1:
return [k, 0 - 1]
return [0, quad.quotient_digits(diff * 100, size, 0)]
def verdict(q):
if q[1] >= 0:
return "off by " + str(q[1]) + "%"
if q[0] >= 17:
return "right to 17 significant digits or more"
return "right to about " + str(q[0]) + " significant digits"
def cancellation():
plain(input("b, a whole number from 3 to 1000000000000> ")) => t
if t == "":
"Cancelled." ^0
return
number(t, 13) => b
if b < 3 or b > 1000000000000:
"Type a whole number from 3 to 1000000000000." ^0
return
quad.small_root(b) => r
quad.exact_small(b, 20) => ex
("x^2 + " + str(b) + "x + 1 = 0, its root nearest zero:") ^0
(" exact, 20 digits: -" + ex[0][0:1] + "." + ex[0][1:len(ex[0])] + "e" + str(ex[1])) ^0
(" textbook formula: " + str(r[0]) + " (" + verdict(quality(r[0], ex)) + ")") ^0
(" stable formula: " + str(r[1]) + " (" + verdict(quality(r[1], ex)) + ")") ^0
"The textbook formula takes -b + sqrt(b^2 - 4) when the two are nearly equal, and the digits they share cancel. The stable one works out the big root, where nothing cancels, and divides: the product of the roots is c / a = 1." ^0
"== Quadratic solver ==" ^0
"ax^2 + bx + c = 0, solved exactly - and the same algebra in floats losing its small root." ^0
True => running
while running:
"" ^0
"1) solve 2) floats and cancellation 3) quit" ^0
plain(input("choice> ")) => choice
if choice == "1":
solve_menu()
elif choice == "2":
cancellation()
elif choice == "3":
False => running
else:
"Pick a number from 1 to 3." ^0
"Bye." ^0
Python projection (main.py)
import quad
def plain(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 number(s, most):
if len(s) == 0 or len(s) > most:
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 ask_coef(name):
t = plain(input(name + "> "))
if t == "":
return []
c = quad.parse_coef(t)
if len(c) == 0:
print("Type a number with at most 9 digits before the point and 3 after, such as -2 or 0.125.")
return ["!"]
return [c[0], quad.trimmed(quad.scaled_text(c[0], 3))]
def equation(a, b, c):
out = ""
for part in [[a, "x^2"], [b, "x"], [c, ""]]:
v = part[0]
if v[0] != 0:
t = v[1]
sign = "+"
if t[0] == "-":
sign = "-"
t = t[1:len(t)]
if t == "1" and part[1] != "":
t = ""
if out == "":
if sign == "-":
out = "-"
out = out + t + part[1]
else:
out = out + " " + sign + " " + t + part[1]
if out == "":
out = "0"
return out + " = 0"
def solve_menu():
a = ask_coef("a")
if len(a) == 0:
print("Cancelled.")
return
if a[0] == "!":
return
b = ask_coef("b")
if len(b) == 0:
print("Cancelled.")
return
if b[0] == "!":
return
c = ask_coef("c")
if len(c) == 0:
print("Cancelled.")
return
if c[0] == "!":
return
print(equation(a, b, c))
r = quad.solve(a[0], b[0], c[0])
if r[0] == "all":
print("Every x solves it: all three coefficients are zero.")
return
if r[0] == "none":
print("No x solves it: it says a number other than zero equals zero.")
return
if r[0] == "linear":
print("With a = 0 it is not quadratic: one root, x = " + r[1] + ".")
return
D = b[0] * b[0] - 4 * a[0] * c[0]
print("Discriminant b^2 - 4ac = " + quad.trimmed(quad.scaled_text(D, 6)) + ".")
if r[0] == "double":
print("It is zero: one double root, x = " + r[1] + ".")
elif r[0] == "two rational":
print("It is the square of a fraction: two rational roots, x = " + r[1] + " and x = " + r[2] + ".")
elif r[0] == "two real":
print("It is positive and not the square of a fraction: two irrational roots, x = (-b - sqrt(D)) / 2a and (-b + sqrt(D)) / 2a:")
print(" x = " + r[2])
print(" x = " + r[3])
print("to 12 places, each digit worked out exactly.")
else:
print("It is negative: two complex roots, x = (-b / 2a) - i (sqrt(-D) / 2a) and its conjugate:")
print(" x = " + r[1] + " - " + r[2] + "i")
print(" x = " + r[1] + " + " + r[2] + "i")
print("the imaginary part to 12 places.")
def printed_value(x):
t = str(x)
if t[0] == "-":
t = t[1:len(t)]
e = 0
at = -1
for i in range(0, len(t)):
if at == -1 and t[i] == "e":
at = i
if at >= 0:
es = t[at + 1:len(t)]
sign = 1
if es[0] == "-":
sign = 0 - 1
es = es[1:len(es)]
elif es[0] == "+":
es = es[1:len(es)]
e = sign * int(es)
t = t[0:at]
ds = ""
after = 0
dot = False
for c in t:
if c == ".":
dot = True
else:
ds = ds + c
if dot:
after = after + 1
return [int(ds), e - after]
def ten_to(k):
t = "1"
for i in range(1, k+1):
t = t + "0"
return int(t)
def quality(x, ex):
p = printed_value(x)
m = int(ex[0])
e = ex[1] - 19
lo = p[1]
if e < lo:
lo = e
xv = p[0] * ten_to(p[1] - lo)
if x < 0.0:
xv = 0 - xv
ev = 0 - m * ten_to(e - lo)
diff = xv - ev
if diff < 0:
diff = 0 - diff
size = 0 - ev
k = 0
while k < 17 and diff * ten_to(k + 1) <= size:
k = k + 1
if k >= 1:
return [k, 0 - 1]
return [0, quad.quotient_digits(diff * 100, size, 0)]
def verdict(q):
if q[1] >= 0:
return "off by " + str(q[1]) + "%"
if q[0] >= 17:
return "right to 17 significant digits or more"
return "right to about " + str(q[0]) + " significant digits"
def cancellation():
t = plain(input("b, a whole number from 3 to 1000000000000> "))
if t == "":
print("Cancelled.")
return
b = number(t, 13)
if b < 3 or b > 1000000000000:
print("Type a whole number from 3 to 1000000000000.")
return
r = quad.small_root(b)
ex = quad.exact_small(b, 20)
print("x^2 + " + str(b) + "x + 1 = 0, its root nearest zero:")
print(" exact, 20 digits: -" + ex[0][0:1] + "." + ex[0][1:len(ex[0])] + "e" + str(ex[1]))
print(" textbook formula: " + str(r[0]) + " (" + verdict(quality(r[0], ex)) + ")")
print(" stable formula: " + str(r[1]) + " (" + verdict(quality(r[1], ex)) + ")")
print("The textbook formula takes -b + sqrt(b^2 - 4) when the two are nearly equal, and the digits they share cancel. The stable one works out the big root, where nothing cancels, and divides: the product of the roots is c / a = 1.")
print("== Quadratic solver ==")
print("ax^2 + bx + c = 0, solved exactly - and the same algebra in floats losing its small root.")
running = True
while running:
print("")
print("1) solve 2) floats and cancellation 3) quit")
choice = plain(input("choice> "))
if choice == "1":
solve_menu()
elif choice == "2":
cancellation()
elif choice == "3":
running = False
else:
print("Pick a number from 1 to 3.")
print("Bye.")
quad.eml
eml# P114 quadratic solver - ax^2 + bx + c = 0 with coefficients of at most
# three decimals, scaled by 1000 to whole numbers so the discriminant
# b^2 - 4ac is exact (the corpus case quadratic-solver, here with complex
# roots too). Square roots and divisions are done digit by digit on whole
# numbers - adding, subtracting, multiplying and comparing only - so every
# printed digit is exact. And as in the corpus case
# catastrophic-cancellation, the textbook formula in floats loses the small
# root of x^2 + bx + 1 when b is large, while the stable one keeps it.
def gcd(a, b):
if a < 0:
0 - a => a
if b < 0:
0 - b => b
while b != 0:
a % b => r
b => a
r => b
return a
def parse_coef(t):
# the coefficient times 1000, or [] (at most 9 digits before the point
# and 3 after)
if len(t) == 0:
return []
1 => sign
t => s
if s[0] == "-":
0 - 1 => sign
s[1:len(s)] => s
"" => whole
"" => frac
False => dot
for c in s:
if c == ".":
if dot:
return []
True => dot
elif c >= "0" and c <= "9":
if dot:
frac + c => frac
else:
whole + c => whole
else:
return []
if whole == "" or len(frac) > 3 or (dot and frac == "") or len(whole) > 9:
return []
while len(frac) < 3:
frac + "0" => frac
return [sign * (int(whole) * 1000 + int(frac))]
def root_digits(n, k):
# floor(sqrt(n) x 10^k) for n >= 0, by the digit-pair method
str(n) => s
for i in [1:k]:
s + "00" => s
if len(s) % 2 == 1:
"0" + s => s
0 => p
0 => r
0 => i
while i < len(s):
r * 100 + int(s[i:i + 2]) => r
0 => d
while (20 * p + d + 1) * (d + 1) <= r:
d + 1 => d
r - (20 * p + d) * d => r
p * 10 + d => p
i + 2 => i
return p
def is_square(n):
root_digits(n, 0) => s
return s * s == n
def quotient_digits(num, den, k):
# floor(num / den x 10^k) for num >= 0, den > 0, by long division
str(num) => s
for i in [1:k]:
s + "0" => s
0 => q
0 => r
for c in s:
r * 10 + int(c) => r
0 => d
while r >= den:
r - den => r
d + 1 => d
q * 10 + d => q
return q
def scaled_text(n, k):
# n / 10^k as a decimal with k places
"" => sign
if n < 0:
"-" => sign
0 - n => n
str(n) => s
while len(s) <= k:
"0" + s => s
if k == 0:
return sign + s
return sign + s[0:len(s) - k] + "." + s[len(s) - k:len(s)]
def trimmed(t):
# drop trailing zeros after the point, and the point if nothing follows
t => s
False => dot
for c in s:
if c == ".":
True => dot
if not dot:
return s
while s[len(s) - 1] == "0":
s[0:len(s) - 1] => s
if s[len(s) - 1] == ".":
s[0:len(s) - 1] => s
return s
def abs_(n):
if n < 0:
return 0 - n
return n
def fraction_text(num, den):
# num / den reduced, the sign in front; den is not 0
if den < 0:
0 - num => num
0 - den => den
gcd(num, den) => g
if g == 0:
1 => g
quotient_digits(abs_(num), g, 0) => n
quotient_digits(den, g, 0) => d
"" => sign
if num < 0:
"-" => sign
if d == 1:
return sign + str(n)
return sign + str(n) + "/" + str(d)
def decimal12(top, root_sign, d, den):
# (top + root_sign x sqrt(d)) / den to 12 places, rounded half away from
# zero; den > 0 and d >= 0. Worked to 18 places first: the floor of the
# root is off by less than one unit there, far below the 12th place.
root_digits(d, 18) => s
"1" => t
for i in [1:18]:
t + "0" => t
top * int(t) + root_sign * s => num
"" => sign
if num < 0:
"-" => sign
0 - num => num
str(quotient_digits(num, den, 0)) => q
while len(q) < 7:
"0" + q => q
int(q[0:len(q) - 6]) => whole
if int(q[len(q) - 6:len(q)]) >= 500000:
whole + 1 => whole
if whole == 0:
"" => sign
return sign + scaled_text(whole, 12)
def solve(A, B, C):
# [kind, ...] for whole coefficients (already scaled):
# ["all"] or ["none"] when A = B = 0; ["linear", root]
# ["double", root]; ["two rational", root, root]
# ["two real", D, root, root] with D not a square, roots to 12 places
# ["complex", real part, imaginary part]
if A == 0:
if B == 0:
if C == 0:
return ["all"]
return ["none"]
return ["linear", fraction_text(0 - C, B)]
if A < 0:
# the same roots, with a positive leading coefficient
0 - A => A
0 - B => B
0 - C => C
B * B - 4 * A * C => D
2 * A => den
if D == 0:
return ["double", fraction_text(0 - B, den)]
if D > 0 and is_square(D):
root_digits(D, 0) => s
return ["two rational", fraction_text(0 - B - s, den), fraction_text(0 - B + s, den)]
if D > 0:
return ["two real", D, decimal12(0 - B, 0 - 1, D, den), decimal12(0 - B, 1, D, den)]
return ["complex", fraction_text(0 - B, den), decimal12(0, 1, 0 - D, den)]
def fsqrt(x):
# a float square root by Newton's method from above, using only
# adding, dividing and comparing - the same in every Python
if x <= 0.0:
return 0.0
x => r
if r < 1.0:
1.0 => r
for i in [1:200]:
(r + x / r) / 2.0 => y
if y >= r:
return r
y => r
return r
def small_root(b):
# x^2 + bx + 1 = 0 for whole b >= 3: [textbook, stable] in floats,
# the root nearer zero
b * 1.0 => fb
fsqrt(fb * fb - 4.0) => s
(0.0 - fb + s) / 2.0 => naive
(0.0 - fb - s) / 2.0 => big
1.0 / big => stable
return [naive, stable]
def exact_small(b, n):
# the root of x^2 + bx + 1 nearest zero is -2 / (b + sqrt(b^2 - 4));
# [its first n significant digits, the power of ten of the first one]
n + len(str(b)) + 5 => k
"1" => t
for i in [1:k]:
t + "0" => t
int(t) => ten
root_digits(b * b - 4, k) => s
str(quotient_digits(2 * ten * ten, b * ten + s, 0)) => q
return [q[0:n], len(q) - k - 1]
Python projection (quad.py)
def gcd(a, b):
if a < 0:
a = 0 - a
if b < 0:
b = 0 - b
while b != 0:
r = a % b
a = b
b = r
return a
def parse_coef(t):
if len(t) == 0:
return []
sign = 1
s = t
if s[0] == "-":
sign = 0 - 1
s = s[1:len(s)]
whole = ""
frac = ""
dot = False
for c in s:
if c == ".":
if dot:
return []
dot = True
elif c >= "0" and c <= "9":
if dot:
frac = frac + c
else:
whole = whole + c
else:
return []
if whole == "" or len(frac) > 3 or dot and frac == "" or len(whole) > 9:
return []
while len(frac) < 3:
frac = frac + "0"
return [sign * (int(whole) * 1000 + int(frac))]
def root_digits(n, k):
s = str(n)
for i in range(1, k+1):
s = s + "00"
if len(s) % 2 == 1:
s = "0" + s
p = 0
r = 0
i = 0
while i < len(s):
r = r * 100 + int(s[i:i + 2])
d = 0
while (20 * p + d + 1) * (d + 1) <= r:
d = d + 1
r = r - (20 * p + d) * d
p = p * 10 + d
i = i + 2
return p
def is_square(n):
s = root_digits(n, 0)
return s * s == n
def quotient_digits(num, den, k):
s = str(num)
for i in range(1, k+1):
s = s + "0"
q = 0
r = 0
for c in s:
r = r * 10 + int(c)
d = 0
while r >= den:
r = r - den
d = d + 1
q = q * 10 + d
return q
def scaled_text(n, k):
sign = ""
if n < 0:
sign = "-"
n = 0 - n
s = str(n)
while len(s) <= k:
s = "0" + s
if k == 0:
return sign + s
return sign + s[0:len(s) - k] + "." + s[len(s) - k:len(s)]
def trimmed(t):
s = t
dot = False
for c in s:
if c == ".":
dot = True
if not dot:
return s
while s[len(s) - 1] == "0":
s = s[0:len(s) - 1]
if s[len(s) - 1] == ".":
s = s[0:len(s) - 1]
return s
def abs_(n):
if n < 0:
return 0 - n
return n
def fraction_text(num, den):
if den < 0:
num = 0 - num
den = 0 - den
g = gcd(num, den)
if g == 0:
g = 1
n = quotient_digits(abs_(num), g, 0)
d = quotient_digits(den, g, 0)
sign = ""
if num < 0:
sign = "-"
if d == 1:
return sign + str(n)
return sign + str(n) + "/" + str(d)
def decimal12(top, root_sign, d, den):
s = root_digits(d, 18)
t = "1"
for i in range(1, 19):
t = t + "0"
num = top * int(t) + root_sign * s
sign = ""
if num < 0:
sign = "-"
num = 0 - num
q = str(quotient_digits(num, den, 0))
while len(q) < 7:
q = "0" + q
whole = int(q[0:len(q) - 6])
if int(q[len(q) - 6:len(q)]) >= 500000:
whole = whole + 1
if whole == 0:
sign = ""
return sign + scaled_text(whole, 12)
def solve(A, B, C):
if A == 0:
if B == 0:
if C == 0:
return ["all"]
return ["none"]
return ["linear", fraction_text(0 - C, B)]
if A < 0:
A = 0 - A
B = 0 - B
C = 0 - C
D = B * B - 4 * A * C
den = 2 * A
if D == 0:
return ["double", fraction_text(0 - B, den)]
if D > 0 and is_square(D):
s = root_digits(D, 0)
return ["two rational", fraction_text(0 - B - s, den), fraction_text(0 - B + s, den)]
if D > 0:
return ["two real", D, decimal12(0 - B, 0 - 1, D, den), decimal12(0 - B, 1, D, den)]
return ["complex", fraction_text(0 - B, den), decimal12(0, 1, 0 - D, den)]
def fsqrt(x):
if x <= 0.0:
return 0.0
r = x
if r < 1.0:
r = 1.0
for i in range(1, 201):
y = (r + x / r) / 2.0
if y >= r:
return r
r = y
return r
def small_root(b):
fb = b * 1.0
s = fsqrt(fb * fb - 4.0)
naive = (0.0 - fb + s) / 2.0
big = (0.0 - fb - s) / 2.0
stable = 1.0 / big
return [naive, stable]
def exact_small(b, n):
k = n + len(str(b)) + 5
t = "1"
for i in range(1, k+1):
t = t + "0"
ten = int(t)
s = root_digits(b * b - 4, k)
q = str(quotient_digits(2 * ten * ten, b * ten + s, 0))
return [q[0:n], len(q) - k - 1]