From 4e99bb593594c623c18aaee1de4a3736f99c7573 Mon Sep 17 00:00:00 2001
From: Niels Dossche <7771979+nielsdos@users.noreply.github.com>
Date: Mon, 20 May 2024 15:10:58 +0200
Subject: [PATCH] Faster BCD into integer parsing

---
 ext/bcmath/libbcmath/src/recmul.c | 72 +++++++++++++++++++++++++++----
 1 file changed, 64 insertions(+), 8 deletions(-)

diff --git a/ext/bcmath/libbcmath/src/recmul.c b/ext/bcmath/libbcmath/src/recmul.c
index 869bc9b0fcc..6fb2776170e 100644
--- a/ext/bcmath/libbcmath/src/recmul.c
+++ b/ext/bcmath/libbcmath/src/recmul.c
@@ -58,12 +58,68 @@ static inline void bc_digits_adjustment(BC_UINT_T *prod_uint, size_t prod_arr_si
 	}
 }
 
+/* This is based on the technique described in https://kholdstare.github.io/technical/2020/05/26/faster-integer-parsing.html.
+ * This function transforms AABBCCDD into 1000 * AA + 100 * BB + 10 * CC + DD,
+ * with the caveat that all components must be in the interval [0, 25] to prevent overflow
+ * due to the multiplication by power of 10 (10 * 25 = 250 is the largest number that fits in a byte).
+ * The advantage of this method instead of using shifts + 3 multiplications is that this is cheaper
+ * due to its divide-and-conquer nature.
+ */
+#if SIZEOF_SIZE_T == 4
+static uint32_t bc_parse_chunk_chars(const char *str)
+{
+	uint32_t tmp;
+	memcpy(&tmp, str, sizeof(tmp));
+#if !BC_LITTLE_ENDIAN
+	tmp = BC_BSWAP(tmp);
+#endif
+
+	uint32_t lower_digits = (tmp & 0x0f000f00) >> 8;
+	uint32_t upper_digits = (tmp & 0x000f000f) * 10;
+
+	tmp = lower_digits + upper_digits;
+
+	lower_digits = (tmp & 0x00ff0000) >> 16;
+	upper_digits = (tmp & 0x000000ff) * 100;
+
+	return lower_digits + upper_digits;
+}
+#elif SIZEOF_SIZE_T == 8
+static uint64_t bc_parse_chunk_chars(const char *str)
+{
+	uint64_t tmp;
+	memcpy(&tmp, str, sizeof(tmp));
+#if !BC_LITTLE_ENDIAN
+	tmp = BC_BSWAP(tmp);
+#endif
+
+	uint64_t lower_digits = (tmp & 0x0f000f000f000f00) >> 8;
+	uint64_t upper_digits = (tmp & 0x000f000f000f000f) * 10;
+
+	tmp = lower_digits + upper_digits;
+
+	lower_digits = (tmp & 0x00ff000000ff0000) >> 16;
+	upper_digits = (tmp & 0x000000ff000000ff) * 100;
+
+	tmp = lower_digits + upper_digits;
+
+	lower_digits = (tmp & 0x0000ffff00000000) >> 32;
+	upper_digits = (tmp & 0x000000000000ffff) * 10000;
+
+	return lower_digits + upper_digits;
+}
+#endif
+
 /*
  * Converts BCD to uint, going backwards from pointer n by the number of
  * characters specified by len.
  */
 static inline BC_UINT_T bc_partial_convert_to_uint(const char *n, size_t len)
 {
+	if (len == BC_MUL_UINT_DIGITS) {
+		return bc_parse_chunk_chars(n - BC_MUL_UINT_DIGITS + 1);
+	}
+
 	BC_UINT_T num = 0;
 	BC_UINT_T base = 1;
 
@@ -226,14 +282,14 @@ static void bc_standard_mul(bc_num n1, size_t n1len, bc_num n2, size_t n2len, bc
 	char *pend = pptr + prodlen - 1;
 	i = 0;
 	while (i < prod_arr_size - 1) {
-		if (BC_MUL_UINT_DIGITS == 4) {
-			bc_write_bcd_representation(prod_uint[i], pend - 3);
-			pend -= 4;
-		} else {
-			bc_write_bcd_representation(prod_uint[i] / 10000, pend - 7);
-			bc_write_bcd_representation(prod_uint[i] % 10000, pend - 3);
-			pend -= 8;
-		}
+#if BC_MUL_UINT_DIGITS == 4
+		bc_write_bcd_representation(prod_uint[i], pend - 3);
+		pend -= 4;
+#else
+		bc_write_bcd_representation(prod_uint[i] / 10000, pend - 7);
+		bc_write_bcd_representation(prod_uint[i] % 10000, pend - 3);
+		pend -= 8;
+#endif
 		i++;
 	}