2008-10-27 16:44:32 -06:00
|
|
|
// Copyright 2009 The Go Authors. All rights reserved.
|
|
|
|
// Use of this source code is governed by a BSD-style
|
|
|
|
// license that can be found in the LICENSE file.
|
|
|
|
|
|
|
|
package Bignum
|
|
|
|
|
|
|
|
// A package for arbitrary precision arithmethic.
|
|
|
|
// It implements the following numeric types:
|
|
|
|
//
|
|
|
|
// - Natural unsigned integer numbers
|
|
|
|
// - Integer signed integer numbers
|
|
|
|
// - Rational rational numbers
|
|
|
|
// - Number scaled rational numbers (contain exponent)
|
|
|
|
|
|
|
|
|
|
|
|
// ----------------------------------------------------------------------------
|
2008-10-28 19:42:26 -06:00
|
|
|
// Representation
|
|
|
|
|
|
|
|
type Word uint64
|
|
|
|
const LogW = 32;
|
2008-10-27 16:44:32 -06:00
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
const LogH = 4; // bits for a hex digit (= "small" number)
|
|
|
|
const H = 1 << LogH;
|
2008-10-27 16:44:32 -06:00
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
const L = LogW - LogH; // must be even (for Mul1)
|
2008-10-27 18:57:31 -06:00
|
|
|
const B = 1 << L;
|
|
|
|
const M = B - 1;
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
// For division
|
|
|
|
|
|
|
|
const (
|
|
|
|
L3 = L / 3;
|
|
|
|
B3 = 1 << L3;
|
|
|
|
M3 = B3 - 1;
|
|
|
|
)
|
|
|
|
|
|
|
|
|
|
|
|
type (
|
|
|
|
Word3 uint32;
|
|
|
|
Natural3 [] Word3;
|
|
|
|
)
|
|
|
|
|
|
|
|
|
|
|
|
// ----------------------------------------------------------------------------
|
|
|
|
// Support
|
|
|
|
|
2008-10-27 16:44:32 -06:00
|
|
|
// TODO replace this with a Go built-in assert
|
2008-10-28 19:42:26 -06:00
|
|
|
func assert(p bool) {
|
2008-10-27 16:44:32 -06:00
|
|
|
if !p {
|
2008-10-28 19:42:26 -06:00
|
|
|
panic("assert failed");
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func init() {
|
|
|
|
assert(L % 2 == 0); // L must be even
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-27 18:57:31 -06:00
|
|
|
func IsSmall(x Word) bool {
|
2008-10-28 19:42:26 -06:00
|
|
|
return x < H;
|
2008-10-27 18:57:31 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-27 16:44:32 -06:00
|
|
|
func Update(x Word) (Word, Word) {
|
|
|
|
return x & M, x >> L;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
export func Dump(x *[]Word) {
|
|
|
|
print("[", len(x), "]");
|
|
|
|
for i := len(x) - 1; i >= 0; i-- {
|
|
|
|
print(" ", x[i]);
|
|
|
|
}
|
|
|
|
println();
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-27 16:44:32 -06:00
|
|
|
// ----------------------------------------------------------------------------
|
2008-10-28 19:42:26 -06:00
|
|
|
// Natural numbers
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
export type Natural []Word;
|
|
|
|
export var NatZero *Natural = new(Natural, 0);
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
export func NewNat(x Word) *Natural {
|
|
|
|
var z *Natural;
|
|
|
|
switch {
|
|
|
|
case x == 0:
|
|
|
|
z = NatZero;
|
|
|
|
case x < B:
|
|
|
|
z = new(Natural, 1);
|
|
|
|
z[0] = x;
|
|
|
|
return z;
|
|
|
|
default:
|
|
|
|
z = new(Natural, 2);
|
|
|
|
z[0], z[1] = Update(x);
|
|
|
|
}
|
|
|
|
return z;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func Normalize(x *Natural) *Natural {
|
|
|
|
i := len(x);
|
|
|
|
for i > 0 && x[i - 1] == 0 { i-- }
|
|
|
|
if i < len(x) {
|
|
|
|
x = x[0 : i]; // trim leading 0's
|
|
|
|
}
|
|
|
|
return x;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func Normalize3(x *Natural3) *Natural3 {
|
|
|
|
i := len(x);
|
|
|
|
for i > 0 && x[i - 1] == 0 { i-- }
|
|
|
|
if i < len(x) {
|
|
|
|
x = x[0 : i]; // trim leading 0's
|
|
|
|
}
|
|
|
|
return x;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-27 16:44:32 -06:00
|
|
|
func (x *Natural) IsZero() bool {
|
|
|
|
return len(x) == 0;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Natural) Add(y *Natural) *Natural {
|
|
|
|
n := len(x);
|
|
|
|
m := len(y);
|
|
|
|
if n < m {
|
2008-10-27 16:44:32 -06:00
|
|
|
return y.Add(x);
|
|
|
|
}
|
2008-10-28 19:42:26 -06:00
|
|
|
assert(n >= m);
|
|
|
|
z := new(Natural, n + 1);
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
i := 0;
|
|
|
|
c := Word(0);
|
2008-10-28 19:42:26 -06:00
|
|
|
for i < m { z[i], c = Update(x[i] + y[i] + c); i++; }
|
|
|
|
for i < n { z[i], c = Update(x[i] + c); i++; }
|
|
|
|
z[i] = c;
|
2008-10-27 16:44:32 -06:00
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
return Normalize(z);
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Natural) Sub(y *Natural) *Natural {
|
|
|
|
n := len(x);
|
|
|
|
m := len(y);
|
|
|
|
assert(n >= m);
|
|
|
|
z := new(Natural, n);
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
i := 0;
|
|
|
|
c := Word(0);
|
2008-10-28 19:42:26 -06:00
|
|
|
for i < m { z[i], c = Update(x[i] - y[i] + c); i++; }
|
|
|
|
for i < n { z[i], c = Update(x[i] + c); i++; }
|
|
|
|
assert(c == 0); // x.Sub(y) must be called with x >= y
|
2008-10-27 16:44:32 -06:00
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
return Normalize(z);
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
// Computes x = x*a + c (in place) for "small" a's.
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x* Natural) MulAdd1(a, c Word) *Natural {
|
|
|
|
assert(IsSmall(a-1) && IsSmall(c));
|
|
|
|
if x.IsZero() || a == 0 {
|
|
|
|
return NewNat(c);
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
2008-10-28 19:42:26 -06:00
|
|
|
n := len(x);
|
2008-10-27 16:44:32 -06:00
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
z := new(Natural, n + 1);
|
|
|
|
for i := 0; i < n; i++ { z[i], c = Update(x[i] * a + c); }
|
|
|
|
z[n] = c;
|
2008-10-27 16:44:32 -06:00
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
return Normalize(z);
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
// Returns z = (x * y) div B, c = (x * y) mod B.
|
2008-10-28 19:42:26 -06:00
|
|
|
func Mul1(x, y Word) (Word, Word) {
|
|
|
|
const L2 = (L + 1) / 2; // TODO check if we can run with odd L
|
2008-10-27 16:44:32 -06:00
|
|
|
const B2 = 1 << L2;
|
|
|
|
const M2 = B2 - 1;
|
2008-10-28 19:42:26 -06:00
|
|
|
|
2008-10-27 16:44:32 -06:00
|
|
|
x0 := x & M2;
|
|
|
|
x1 := x >> L2;
|
|
|
|
|
|
|
|
y0 := y & M2;
|
|
|
|
y1 := y >> L2;
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
z0 := x0*y0;
|
|
|
|
z1 := x1*y0 + x0*y1 + z0 >> L2; z0 &= M2;
|
|
|
|
z2 := x1*y1 + z1 >> L2; z1 &= M2;
|
|
|
|
|
|
|
|
return z1 << L2 | z0, z2;
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Natural) Mul(y *Natural) *Natural {
|
2008-10-27 16:44:32 -06:00
|
|
|
if x.IsZero() || y.IsZero() {
|
|
|
|
return NatZero;
|
|
|
|
}
|
|
|
|
xl := len(x);
|
|
|
|
yl := len(y);
|
|
|
|
if xl < yl {
|
|
|
|
return y.Mul(x); // for speed
|
|
|
|
}
|
2008-10-28 19:42:26 -06:00
|
|
|
assert(xl >= yl && yl > 0);
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
// initialize z
|
|
|
|
zl := xl + yl;
|
|
|
|
z := new(Natural, zl);
|
|
|
|
|
|
|
|
for j := 0; j < yl; j++ {
|
|
|
|
d := y[j];
|
|
|
|
if d != 0 {
|
2008-10-28 19:42:26 -06:00
|
|
|
k := j;
|
2008-10-27 16:44:32 -06:00
|
|
|
c := Word(0);
|
|
|
|
for i := 0; i < xl; i++ {
|
|
|
|
// compute z[k] += x[i] * d + c;
|
|
|
|
t := z[k] + c;
|
|
|
|
var z1 Word;
|
|
|
|
z1, c = Mul1(x[i], d);
|
|
|
|
t += z1;
|
|
|
|
z[k] = t & M;
|
|
|
|
c += t >> L;
|
|
|
|
k++;
|
|
|
|
}
|
2008-10-28 19:42:26 -06:00
|
|
|
z[k] = c;
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
return Normalize(z);
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func Shl1(x Word, s int) (Word, Word) {
|
|
|
|
return 0, 0
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func Shr1(x Word, s int) (Word, Word) {
|
|
|
|
return 0, 0
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Natural) Shl(s int) *Natural {
|
|
|
|
panic("incomplete");
|
|
|
|
|
|
|
|
if s == 0 {
|
|
|
|
return x;
|
|
|
|
}
|
|
|
|
|
|
|
|
S := s/L;
|
|
|
|
s = s%L;
|
|
|
|
n := len(x) + S + 1;
|
|
|
|
z := new(Natural, n);
|
|
|
|
|
|
|
|
c := Word(0);
|
|
|
|
for i := 0; i < n; i++ {
|
|
|
|
z[i + S], c = Shl1(x[i], s);
|
|
|
|
}
|
|
|
|
z[n + S] = c;
|
|
|
|
|
|
|
|
return Normalize(z);
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func (x *Natural) Shr(s uint) *Natural {
|
|
|
|
panic("incomplete");
|
|
|
|
|
|
|
|
if s == 0 {
|
|
|
|
return x;
|
|
|
|
}
|
|
|
|
return nil
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func SplitBase(x *Natural) *Natural3 {
|
|
|
|
xl := len(x);
|
|
|
|
z := new(Natural3, xl * 3);
|
|
|
|
for i, j := 0, 0; i < xl; i, j = i + 1, j + 3 {
|
|
|
|
t := x[i];
|
|
|
|
z[j] = Word3(t & M3); t >>= L3; j++;
|
|
|
|
z[j] = Word3(t & M3); t >>= L3; j++;
|
|
|
|
z[j] = Word3(t & M3); t >>= L3; j++;
|
|
|
|
}
|
|
|
|
return Normalize3(z);
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func Scale(x *Natural, f Word) *Natural3 {
|
|
|
|
return nil;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func TrialDigit(r, d *Natural3, k, m int) Word {
|
|
|
|
km := k + m;
|
|
|
|
assert(2 <= m && m <= km);
|
|
|
|
r3 := (Word(r[km]) << L3 + Word(r[km - 1])) << L3 + Word(r[km - 2]);
|
|
|
|
d2 := Word(d[m - 1]) << L3 + Word(d[m - 2]);
|
|
|
|
qt := r3 / d2;
|
|
|
|
if qt >= B {
|
|
|
|
qt = B - 1;
|
|
|
|
}
|
|
|
|
return qt;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func DivMod(x, y *Natural) {
|
|
|
|
xl := len(x);
|
|
|
|
yl := len(y);
|
|
|
|
assert(2 <= yl && yl <= xl); // use special-case algorithm otherwise
|
|
|
|
|
|
|
|
f := B / (y[yl - 1] + 1);
|
|
|
|
r := Scale(x, f);
|
|
|
|
d := Scale(y, f);
|
|
|
|
n := len(r);
|
|
|
|
m := len(d);
|
|
|
|
|
|
|
|
for k := n - m; k >= 0; k-- {
|
|
|
|
qt := TrialDigit(r, d, k, m);
|
|
|
|
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func (x *Natural) Div(y *Natural) *Natural {
|
2008-10-27 16:44:32 -06:00
|
|
|
panic("UNIMPLEMENTED");
|
|
|
|
return nil;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Natural) Mod(y *Natural) *Natural {
|
2008-10-27 16:44:32 -06:00
|
|
|
panic("UNIMPLEMENTED");
|
|
|
|
return nil;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Natural) Cmp(y *Natural) int {
|
2008-10-27 16:44:32 -06:00
|
|
|
xl := len(x);
|
|
|
|
yl := len(y);
|
|
|
|
|
|
|
|
if xl != yl || xl == 0 {
|
|
|
|
return xl - yl;
|
|
|
|
}
|
|
|
|
|
|
|
|
i := xl - 1;
|
|
|
|
for i > 0 && x[i] == y[i] { i--; }
|
|
|
|
|
|
|
|
d := 0;
|
|
|
|
switch {
|
|
|
|
case x[i] < y[i]: d = -1;
|
|
|
|
case x[i] > y[i]: d = 1;
|
|
|
|
}
|
2008-10-28 19:42:26 -06:00
|
|
|
|
2008-10-27 16:44:32 -06:00
|
|
|
return d;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func (x *Natural) Log() int {
|
2008-10-28 19:42:26 -06:00
|
|
|
n := len(x);
|
|
|
|
if n == 0 { return 0; }
|
|
|
|
assert(n > 0);
|
|
|
|
|
|
|
|
c := (n - 1) * L;
|
|
|
|
for t := x[n - 1]; t != 0; t >>= 1 { c++ };
|
2008-10-27 16:44:32 -06:00
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
return c;
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Natural) And(y *Natural) *Natural {
|
|
|
|
n := len(x);
|
|
|
|
m := len(y);
|
|
|
|
if n < m {
|
2008-10-27 16:44:32 -06:00
|
|
|
return y.And(x);
|
|
|
|
}
|
2008-10-28 19:42:26 -06:00
|
|
|
assert(n >= m);
|
|
|
|
z := new(Natural, n);
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
i := 0;
|
2008-10-28 19:42:26 -06:00
|
|
|
for i < m { z[i] = x[i] & y[i]; i++; }
|
|
|
|
for i < n { z[i] = x[i]; i++; }
|
2008-10-27 16:44:32 -06:00
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
return Normalize(z);
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Natural) Or(y *Natural) *Natural {
|
|
|
|
n := len(x);
|
|
|
|
m := len(y);
|
|
|
|
if n < m {
|
2008-10-27 18:57:31 -06:00
|
|
|
return y.Or(x);
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
2008-10-28 19:42:26 -06:00
|
|
|
assert(n >= m);
|
|
|
|
z := new(Natural, n);
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
i := 0;
|
2008-10-28 19:42:26 -06:00
|
|
|
for i < m { z[i] = x[i] | y[i]; i++; }
|
|
|
|
for i < n { z[i] = x[i]; i++; }
|
2008-10-27 16:44:32 -06:00
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
return Normalize(z);
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Natural) Xor(y *Natural) *Natural {
|
|
|
|
n := len(x);
|
|
|
|
m := len(y);
|
|
|
|
if n < m {
|
2008-10-27 18:57:31 -06:00
|
|
|
return y.Xor(x);
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
2008-10-28 19:42:26 -06:00
|
|
|
assert(n >= m);
|
|
|
|
z := new(Natural, n);
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
i := 0;
|
2008-10-28 19:42:26 -06:00
|
|
|
for i < m { z[i] = x[i] ^ y[i]; i++; }
|
|
|
|
for i < n { z[i] = x[i]; i++; }
|
2008-10-27 16:44:32 -06:00
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
return Normalize(z);
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func Copy(x *Natural) *Natural {
|
2008-10-28 19:42:26 -06:00
|
|
|
z := new(Natural, len(x));
|
|
|
|
//*z = *x; // BUG assignment does't work yet
|
|
|
|
for i := len(x) - 1; i >= 0; i-- { z[i] = x[i]; }
|
2008-10-27 16:44:32 -06:00
|
|
|
return z;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
// Computes x = x div d (in place - the recv maybe modified) for "small" d's.
|
|
|
|
// Returns updated x and x mod d.
|
|
|
|
func (x *Natural) DivMod1(d Word) (*Natural, Word) {
|
|
|
|
assert(0 < d && IsSmall(d - 1));
|
|
|
|
|
2008-10-27 16:44:32 -06:00
|
|
|
c := Word(0);
|
2008-10-28 19:42:26 -06:00
|
|
|
for i := len(x) - 1; i >= 0; i-- {
|
|
|
|
var LL Word = L; // BUG shift broken for const L
|
|
|
|
c = c << LL + x[i];
|
|
|
|
x[i] = c / d;
|
|
|
|
c %= d;
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
return Normalize(x), c;
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Natural) String(base Word) string {
|
2008-10-27 16:44:32 -06:00
|
|
|
if x.IsZero() {
|
|
|
|
return "0";
|
|
|
|
}
|
|
|
|
|
|
|
|
// allocate string
|
2008-10-28 19:42:26 -06:00
|
|
|
// TODO n is too small for bases < 10!!!
|
|
|
|
assert(base >= 10); // for now
|
2008-10-27 16:44:32 -06:00
|
|
|
// approx. length: 1 char for 3 bits
|
|
|
|
n := x.Log()/3 + 1; // +1 (round up)
|
|
|
|
s := new([]byte, n);
|
|
|
|
|
|
|
|
// convert
|
2008-10-28 19:42:26 -06:00
|
|
|
const hex = "0123456789abcdef";
|
2008-10-27 16:44:32 -06:00
|
|
|
i := n;
|
|
|
|
x = Copy(x); // don't destroy recv
|
|
|
|
for !x.IsZero() {
|
|
|
|
i--;
|
|
|
|
var d Word;
|
2008-10-28 19:42:26 -06:00
|
|
|
x, d = x.DivMod1(base);
|
|
|
|
s[i] = hex[d];
|
2008-10-27 16:44:32 -06:00
|
|
|
};
|
|
|
|
|
|
|
|
return string(s[i : n]);
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func MulRange(a, b Word) *Natural {
|
|
|
|
switch {
|
2008-10-28 19:42:26 -06:00
|
|
|
case a > b: return NewNat(1);
|
|
|
|
case a == b: return NewNat(a);
|
|
|
|
case a + 1 == b: return NewNat(a).Mul(NewNat(b));
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
m := (a + b) >> 1;
|
2008-10-28 19:42:26 -06:00
|
|
|
assert(a <= m && m < b);
|
2008-10-27 16:44:32 -06:00
|
|
|
return MulRange(a, m).Mul(MulRange(m + 1, b));
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
export func Fact(n Word) *Natural {
|
2008-10-28 19:42:26 -06:00
|
|
|
// Using MulRange() instead of the basic for-loop
|
|
|
|
// lead to faster factorial computation.
|
|
|
|
return MulRange(2, n);
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func HexValue(ch byte) Word {
|
|
|
|
d := Word(H);
|
|
|
|
switch {
|
|
|
|
case '0' <= ch && ch <= '9': d = Word(ch - '0');
|
|
|
|
case 'a' <= ch && ch <= 'f': d = Word(ch - 'a') + 10;
|
|
|
|
case 'A' <= ch && ch <= 'F': d = Word(ch - 'A') + 10;
|
|
|
|
}
|
|
|
|
return d;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
// TODO auto-detect base if base argument is 0
|
|
|
|
export func NatFromString(s string, base Word) *Natural {
|
2008-10-27 16:44:32 -06:00
|
|
|
x := NatZero;
|
2008-10-28 19:42:26 -06:00
|
|
|
for i := 0; i < len(s); i++ {
|
|
|
|
d := HexValue(s[i]);
|
|
|
|
if d < base {
|
|
|
|
x = x.MulAdd1(base, d);
|
|
|
|
} else {
|
|
|
|
break;
|
|
|
|
}
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
return x;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
// ----------------------------------------------------------------------------
|
2008-10-28 19:42:26 -06:00
|
|
|
// Integer numbers
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
export type Integer struct {
|
|
|
|
sign bool;
|
|
|
|
mant *Natural;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Integer) Add(y *Integer) *Integer {
|
2008-10-27 16:44:32 -06:00
|
|
|
var z *Integer;
|
|
|
|
if x.sign == y.sign {
|
|
|
|
// x + y == x + y
|
|
|
|
// (-x) + (-y) == -(x + y)
|
|
|
|
z = &Integer{x.sign, x.mant.Add(y.mant)};
|
|
|
|
} else {
|
|
|
|
// x + (-y) == x - y == -(y - x)
|
|
|
|
// (-x) + y == y - x == -(x - y)
|
|
|
|
if x.mant.Cmp(y.mant) >= 0 {
|
|
|
|
z = &Integer{false, x.mant.Sub(y.mant)};
|
|
|
|
} else {
|
|
|
|
z = &Integer{true, y.mant.Sub(x.mant)};
|
|
|
|
}
|
|
|
|
}
|
|
|
|
if x.sign {
|
|
|
|
z.sign = !z.sign;
|
|
|
|
}
|
|
|
|
return z;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Integer) Sub(y *Integer) *Integer {
|
2008-10-27 16:44:32 -06:00
|
|
|
var z *Integer;
|
|
|
|
if x.sign != y.sign {
|
|
|
|
// x - (-y) == x + y
|
|
|
|
// (-x) - y == -(x + y)
|
|
|
|
z = &Integer{x.sign, x.mant.Add(y.mant)};
|
|
|
|
} else {
|
|
|
|
// x - y == x - y == -(y - x)
|
|
|
|
// (-x) - (-y) == y - x == -(x - y)
|
|
|
|
if x.mant.Cmp(y.mant) >= 0 {
|
|
|
|
z = &Integer{false, x.mant.Sub(y.mant)};
|
|
|
|
} else {
|
|
|
|
z = &Integer{true, y.mant.Sub(x.mant)};
|
|
|
|
}
|
|
|
|
}
|
|
|
|
if x.sign {
|
|
|
|
z.sign = !z.sign;
|
|
|
|
}
|
|
|
|
return z;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Integer) Mul(y *Integer) *Integer {
|
2008-10-27 16:44:32 -06:00
|
|
|
// x * y == x * y
|
|
|
|
// x * (-y) == -(x * y)
|
|
|
|
// (-x) * y == -(x * y)
|
|
|
|
// (-x) * (-y) == x * y
|
|
|
|
return &Integer{x.sign != y.sign, x.mant.Mul(y.mant)};
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Integer) Div(y *Integer) *Integer {
|
2008-10-27 16:44:32 -06:00
|
|
|
panic("UNIMPLEMENTED");
|
|
|
|
return nil;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Integer) Mod(y *Integer) *Integer {
|
2008-10-27 16:44:32 -06:00
|
|
|
panic("UNIMPLEMENTED");
|
|
|
|
return nil;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Integer) Cmp(y *Integer) int {
|
2008-10-27 16:44:32 -06:00
|
|
|
panic("UNIMPLEMENTED");
|
|
|
|
return 0;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Integer) String(base Word) string {
|
2008-10-27 18:57:31 -06:00
|
|
|
if x.mant.IsZero() {
|
|
|
|
return "0";
|
|
|
|
}
|
|
|
|
var s string;
|
|
|
|
if x.sign {
|
2008-10-28 19:42:26 -06:00
|
|
|
s = "-";
|
2008-10-27 18:57:31 -06:00
|
|
|
}
|
2008-10-28 19:42:26 -06:00
|
|
|
return s + x.mant.String(base);
|
2008-10-27 18:57:31 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
export func IntFromString(s string, base Word) *Integer {
|
2008-10-27 16:44:32 -06:00
|
|
|
// get sign, if any
|
|
|
|
sign := false;
|
|
|
|
if len(s) > 0 && (s[0] == '-' || s[0] == '+') {
|
|
|
|
sign = s[0] == '-';
|
|
|
|
}
|
2008-10-28 19:42:26 -06:00
|
|
|
return &Integer{sign, NatFromString(s[1 : len(s)], base)};
|
2008-10-27 16:44:32 -06:00
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
// ----------------------------------------------------------------------------
|
2008-10-28 19:42:26 -06:00
|
|
|
// Rational numbers
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
export type Rational struct {
|
|
|
|
a, b *Integer; // a = numerator, b = denominator
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
func NewRat(a, b *Integer) *Rational {
|
|
|
|
// TODO normalize the rational
|
|
|
|
return &Rational{a, b};
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Rational) Add(y *Rational) *Rational {
|
2008-10-27 16:44:32 -06:00
|
|
|
return NewRat((x.a.Mul(y.b)).Add(x.b.Mul(y.a)), x.b.Mul(y.b));
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Rational) Sub(y *Rational) *Rational {
|
2008-10-27 16:44:32 -06:00
|
|
|
return NewRat((x.a.Mul(y.b)).Sub(x.b.Mul(y.a)), x.b.Mul(y.b));
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Rational) Mul(y *Rational) *Rational {
|
2008-10-27 16:44:32 -06:00
|
|
|
return NewRat(x.a.Mul(y.a), x.b.Mul(y.b));
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Rational) Div(y *Rational) *Rational {
|
2008-10-27 16:44:32 -06:00
|
|
|
return NewRat(x.a.Mul(y.b), x.b.Mul(y.a));
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Rational) Mod(y *Rational) *Rational {
|
2008-10-27 16:44:32 -06:00
|
|
|
panic("UNIMPLEMENTED");
|
|
|
|
return nil;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
2008-10-28 19:42:26 -06:00
|
|
|
func (x *Rational) Cmp(y *Rational) int {
|
2008-10-27 16:44:32 -06:00
|
|
|
panic("UNIMPLEMENTED");
|
|
|
|
return 0;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
export func RatFromString(s string) *Rational {
|
|
|
|
panic("UNIMPLEMENTED");
|
|
|
|
return nil;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
// ----------------------------------------------------------------------------
|
2008-10-28 19:42:26 -06:00
|
|
|
// Scaled numbers
|
2008-10-27 16:44:32 -06:00
|
|
|
|
|
|
|
export type Number struct {
|
|
|
|
mant *Rational;
|
|
|
|
exp Integer;
|
|
|
|
}
|