From 2884fa55060ba870b217d38f2e3d3129f1108a55 Mon Sep 17 00:00:00 2001
From: Jeroen van Rijn <Kelimion@users.noreply.github.com>
Date: Sat, 24 Jul 2021 19:01:55 +0200
Subject: [PATCH] big: add div by 3.

---
 core/math/big/basic.odin   | 104 +++++++++++++++++++++++++++++++++++++
 core/math/big/example.odin |  21 +++++---
 2 files changed, 119 insertions(+), 6 deletions(-)

diff --git a/core/math/big/basic.odin b/core/math/big/basic.odin
index 53733fdb0..00411f80a 100644
--- a/core/math/big/basic.odin
+++ b/core/math/big/basic.odin
@@ -953,4 +953,108 @@ _int_sqr :: proc(dest, src: ^Int) -> (err: Error) {
 	swap(dest, t);
 	destroy(t);
 	return err;
+}
+
+/*
+	Fivide by three (based on routine from MPI and the GMP manual).
+*/
+_int_div_3 :: proc(quotient, numerator: ^Int) -> (remainder: int, err: Error) {
+	/*
+		b = 2**MP_DIGIT_BIT / 3
+	*/
+ 	b := _WORD(1) << _WORD(_DIGIT_BITS) / _WORD(3);
+
+	q := &Int{};
+	if err = grow(q, numerator.used); err != .None { return -1, err; }
+	q.used = numerator.used;
+	q.sign = numerator.sign;
+
+	w, t: _WORD;
+	for ix := numerator.used; ix >= 0; ix -= 1 {
+		w = (w << _WORD(_DIGIT_BITS)) | _WORD(numerator.digit[ix]);
+		if w >= 3 {
+			/*
+				Multiply w by [1/3].
+			*/
+			t = (w * b) >> _WORD(_DIGIT_BITS);
+
+			/*
+				Now subtract 3 * [w/3] from w, to get the remainder.
+			*/
+			w -= t+t+t;
+
+			/*
+				Fixup the remainder as required since the optimization is not exact.
+			*/
+			for w >= 3 {
+				t += 1;
+				w -= 3;
+			}
+		} else {
+			t = 0;
+		}
+		q.digit[ix] = DIGIT(t);
+	}
+
+	remainder = int(w);
+
+	/*
+		[optional] store the quotient.
+	*/
+	if quotient != nil {
+		err = clamp(q);
+ 		swap(q, quotient);
+ 	}
+	destroy(q);
+	return remainder, .None;
+}
+
+/*
+	Slower bit-bang division... also smaller.
+*/
+_int_div_small :: proc(quotient, remainder, numerator, denominator: ^Int) -> (err: Error) {
+
+// 	mp_int ta, tb, tq, q;
+// 	int n;
+// 	bool neg;
+// 	mp_err err;
+
+// 	/* init our temps */
+// 	if ((err = mp_init_multi(&ta, &tb, &tq, &q, NULL)) != MP_OKAY) {
+// 		return err;
+// 	}
+
+// 	mp_set(&tq, 1uL);
+// 	n = mp_count_bits(a) - mp_count_bits(b);
+// 	if ((err = mp_abs(a, &ta)) != MP_OKAY)                         goto LBL_ERR;
+// 	if ((err = mp_abs(b, &tb)) != MP_OKAY)                         goto LBL_ERR;
+// 	if ((err = mp_mul_2d(&tb, n, &tb)) != MP_OKAY)                 goto LBL_ERR;
+// 	if ((err = mp_mul_2d(&tq, n, &tq)) != MP_OKAY)                 goto LBL_ERR;
+
+// 	while (n-- >= 0) {
+// 		if (mp_cmp(&tb, &ta) != MP_GT) {
+// 			if ((err = mp_sub(&ta, &tb, &ta)) != MP_OKAY)            goto LBL_ERR;
+// 			if ((err = mp_add(&q, &tq, &q)) != MP_OKAY)              goto LBL_ERR;
+// 		}
+// 		if ((err = mp_div_2d(&tb, 1, &tb, NULL)) != MP_OKAY)        goto LBL_ERR;
+// 		if ((err = mp_div_2d(&tq, 1, &tq, NULL)) != MP_OKAY)        goto LBL_ERR;
+// 	}
+
+// 	/* now q == quotient and ta == remainder */
+
+// 	neg = (a->sign != b->sign);
+// 	if (c != NULL) {
+// 		mp_exch(c, &q);
+// 		c->sign = ((neg && !mp_iszero(c)) ? MP_NEG : MP_ZPOS);
+// 	}
+// 	if (d != NULL) {
+// 		mp_exch(d, &ta);
+// 		d->sign = (mp_iszero(d) ? MP_ZPOS : a->sign);
+// 	}
+// LBL_ERR:
+// 	mp_clear_multi(&ta, &tb, &tq, &q, NULL);
+// 	return err;
+
+
+	return .None;
 }
\ No newline at end of file
diff --git a/core/math/big/example.odin b/core/math/big/example.odin
index c0de1d7c3..103259b3e 100644
--- a/core/math/big/example.odin
+++ b/core/math/big/example.odin
@@ -26,6 +26,7 @@ print_configation :: proc() {
 	SQR_KARATSUBA_CUTOFF %v
 	MUL_TOOM_CUTOFF      %v
 	SQR_TOOM_CUTOFF      %v
+
 `, _DIGIT_BITS,
 _MIN_DIGIT_COUNT,
 _MAX_DIGIT_COUNT,
@@ -37,8 +38,6 @@ _SQR_KARATSUBA_CUTOFF,
 _MUL_TOOM_CUTOFF,
 _SQR_TOOM_CUTOFF,
 );
-
-	fmt.println();
 }
 
 print :: proc(name: string, a: ^Int, base := i8(16)) {
@@ -54,11 +53,21 @@ print :: proc(name: string, a: ^Int, base := i8(16)) {
 
 demo :: proc() {
 	err: Error;
-	a, b, c := &Int{}, &Int{}, &Int{};
-	defer destroy(a, b, c);
+	i:   int;
+	destination, source, quotient, remainder, numerator, denominator := &Int{}, &Int{}, &Int{}, &Int{}, &Int{}, &Int{};
+	defer destroy(destination, source, quotient, remainder, numerator, denominator);
 
-	err = set(a, 5125);
-	print("a", a);
+	err = set(numerator,   15625);
+	err = set(denominator,     3);
+
+	print("numerator  ", numerator,   10);
+	print("denominator", denominator, 10);
+
+	i, err = _int_div_3(quotient, numerator);
+
+	print("quotient   ", quotient,    10);
+	fmt.println("remainder  ", i);
+	fmt.println("error", err);
 }
 
 main :: proc() {