56 lines
1.5 KiB
Plaintext
56 lines
1.5 KiB
Plaintext
# -*- ObjectIcon -*-
|
|
|
|
import exception
|
|
import io
|
|
|
|
procedure main ()
|
|
test_euclid_div ()
|
|
io.write (inverse (42, 2017))
|
|
end
|
|
|
|
procedure inverse (a, n) # FAILS if there is no inverse.
|
|
local t, newt, r, newr, quotient, tmp
|
|
|
|
if n <= 0 then throw ("non-positive modulus")
|
|
t := 0; newt := 1
|
|
r := n; newr := a
|
|
while newr ~= 0 do
|
|
{
|
|
quotient := euclid_div (r, newr)
|
|
tmp := newt; newt := t - (quotient * newt); t := tmp
|
|
tmp := newr; newr := r - (quotient * newr); r := tmp
|
|
}
|
|
r <= 1 | fail
|
|
return (if t < 0 then t + n else t)
|
|
end
|
|
|
|
procedure euclid_div (x, y)
|
|
# This kind of integer division always gives a remainder between 0
|
|
# and abs(y)-1, inclusive. Thus the remainder is always a LEAST
|
|
# RESIDUE modulo abs(y). (If y is a positive modulus, then only the
|
|
# floor division branch is used.)
|
|
return \
|
|
if 0 <= y then # Do floor division.
|
|
(if 0 <= x then x / y
|
|
else if (-x) % y = 0 then -((-x) / y)
|
|
else -((-x) / y) - 1)
|
|
else # Do ceiling division.
|
|
(if 0 <= x then -(x / (-y))
|
|
else if (-x) % (-y) = 0 then ((-x) / (-y))
|
|
else ((-x) / (-y)) + 1)
|
|
end
|
|
|
|
procedure test_euclid_div ()
|
|
local x, y, q, r
|
|
|
|
every x := -100 to 100 do
|
|
every y := -100 to 100 & y ~= 0 do
|
|
{
|
|
q := euclid_div (x, y)
|
|
r := x - (q * y)
|
|
if r < 0 | abs (y) <= r then
|
|
# A remainder was outside the expected range.
|
|
throw ("Test of euclid_div failed.")
|
|
}
|
|
end
|