Upgrade bc to 6.6.1 am: 306341a225 am: cdb14c26d9 am: 9e6368af49

Original change: https://android-review.googlesource.com/c/platform/external/bc/+/2768152

Change-Id: Iae388bc78f6759b2c75b600ee00e7e1a4e00ea4d
Signed-off-by: Automerger Merge Worker <android-build-automerger-merge-worker@system.gserviceaccount.com>
diff --git a/.gitignore b/.gitignore
index ea0bf73..d2b5712 100644
--- a/.gitignore
+++ b/.gitignore
@@ -70,6 +70,7 @@
 *.gcno
 *.gcov
 *.html
+*.css
 *.profraw
 
 core.*
diff --git a/METADATA b/METADATA
index 3c6ce98..d56b217 100644
--- a/METADATA
+++ b/METADATA
@@ -1,6 +1,6 @@
 # This project was upgraded with external_updater.
 # Usage: tools/external_updater/updater.sh update bc
-# For more info, check https://cs.android.com/android/platform/superproject/+/master:tools/external_updater/README.md
+# For more info, check https://cs.android.com/android/platform/superproject/+/main:tools/external_updater/README.md
 
 name: "bc"
 description: "An implementation of the POSIX bc calculator with GNU extensions and dc."
@@ -9,11 +9,11 @@
     type: GIT
     value: "https://github.com/gavinhoward/bc"
   }
-  version: "6.6.0"
+  version: "6.6.1"
   license_type: NOTICE
   last_upgrade_date {
     year: 2023
-    month: 5
-    day: 31
+    month: 9
+    day: 29
   }
 }
diff --git a/Makefile.in b/Makefile.in
index 55e2e4a..e1309cd 100644
--- a/Makefile.in
+++ b/Makefile.in
@@ -554,7 +554,7 @@
 clean_coverage:
 	@printf 'Cleaning coverage files...\n'
 	@$(RM) -f *.gcov
-	@$(RM) -f *.html
+	@$(RM) -f *.html *.css
 	@$(RM) -f *.gcda *.gcno
 	@$(RM) -f *.profraw
 	@$(RM) -f $(GCDA) $(GCNO)
diff --git a/NEWS.md b/NEWS.md
index de3b350..28639ea 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,12 @@
 # News
 
+## 6.6.1
+
+This is a production release with an improved `p()` function in the [extended
+math library][16].
+
+Users who don't care do not need to upgrade.
+
 ## 6.6.0
 
 This is a production release with two bug fixes and one change.
diff --git a/gen/lib2.bc b/gen/lib2.bc
index ba3f76b..1d7d48c 100644
--- a/gen/lib2.bc
+++ b/gen/lib2.bc
@@ -34,10 +34,34 @@
  */
 
 define p(x,y){
-	auto a
+	auto a,i,s,z
+	if(y==0)return 1@scale
+	if(x==0){
+		if(y>0)return 0
+		return 1/0
+	}
 	a=y$
 	if(y==a)return(x^a)@scale
-	return e(y*l(x))
+	z=0
+	if(x<1){
+		y=-y
+		a=-a
+		z=x
+		x=1/x
+	}
+	if(y<0){
+		return e(y*l(x))
+	}
+	i=x^a
+	s=scale
+	scale+=length(i)+5
+	if(z){
+		x=1/z
+		i=x^a
+	}
+	i*=e((y-a)*l(x))
+	scale=s
+	return i@scale
 }
 define r(x,p){
 	auto t,n
diff --git a/include/bcl.h b/include/bcl.h
index 0908e21..d3a9f42 100644
--- a/include/bcl.h
+++ b/include/bcl.h
@@ -36,9 +36,6 @@
 #ifndef BC_BCL_H
 #define BC_BCL_H
 
-// TODO: Add a generation index when building with Valgrind to check for
-// use-after-free's or double frees.
-
 #include <stdbool.h>
 #include <stdlib.h>
 #include <limits.h>
diff --git a/include/version.h b/include/version.h
index a4df383..c7d2449 100644
--- a/include/version.h
+++ b/include/version.h
@@ -37,6 +37,6 @@
 #define BC_VERSION_H
 
 /// The current version.
-#define VERSION 6.6.0
+#define VERSION 6.6.1
 
 #endif // BC_VERSION_H
diff --git a/manuals/algorithms.md b/manuals/algorithms.md
index 4d7a0ed..e0a5e8a 100644
--- a/manuals/algorithms.md
+++ b/manuals/algorithms.md
@@ -193,6 +193,74 @@
 
 It has a complexity of `O(n^3)` because both `e()` and `l()` do.
 
+However, there are details to this algorithm, described by the author,
+TediusTimmy, in GitHub issue #69.
+
+First, check if the exponent is 0. If it is, return 1 at the appropriate
+`scale`.
+
+Next, check if the number is 0. If so, check if the exponent is greater than
+zero; if it is, return 0. If the exponent is less than 0, error (with a divide
+by 0) because that is undefined.
+
+Next, check if the exponent is actually an integer, and if it is, use the
+exponentiation operator.
+
+At the `z=0` line is the start of the meat of the new code.
+
+`z` is set to zero as a flag and as a value. What I mean by that will be clear
+later.
+
+Then we check if the number is less than 0. If it is, we negate the exponent
+(and the integer version of the exponent, which we calculated earlier to check
+if it was an integer). We also save the number in `z`; being non-zero is a flag
+for later and a value to be used. Then we store the reciprocal of the number in
+itself.
+
+All of the above paragraph will not make sense unless you remember the
+relationship `l(x) == -l(1/x)`; we negated the exponent, which is equivalent to
+the negative sign in that relationship, and we took the reciprocal of the
+number, which is equivalent to the reciprocal in the relationship.
+
+But what if the number is negative? We ignore that for now because we eventually
+call `l(x)`, which will raise an error if `x` is negative.
+
+Now, we can keep going.
+
+If at this point, the exponent is negative, we need to use the original formula
+(`e(y * l(x))`) and return that result because the result will go to zero
+anyway.
+
+But if we did *not* return, we know the exponent is *not* negative, so we can
+get clever.
+
+We then compute the integral portion of the power by computing the number to
+power of the integral portion of the exponent.
+
+Then we have the most clever trick: we add the length of that integer power (and
+a little extra) to the `scale`. Why? Because this will ensure that the next part
+is calculated to at least as many digits as should be in the integer *plus* any
+extra `scale` that was wanted.
+
+Then we check `z`, which, if it is not zero, is the original value of the
+number. If it is not zero, we need to take the take the reciprocal *again*
+because now we have the correct `scale`. And we *also* have to calculate the
+integer portion of the power again.
+
+Then we need to calculate the fractional portion of the number. We do this by
+using the original formula, but we instead of calculating `e(y * l(x))`, we
+calculate `e((y - a) * l(x))`, where `a` is the integer portion of `y`. It's
+easy to see that `y - a` will be just the fractional portion of `y` (the
+exponent), so this makes sense.
+
+But then we *multiply* it into the integer portion of the power. Why? Because
+remember: we're dealing with an exponent and a power; the relationship is
+`x^(y+z) == (x^y)*(x^z)`.
+
+So we multiply it into the integer portion of the power.
+
+Finally, we set the result to the `scale`.
+
 ### Rounding (`bc` Math Library 2 Only)
 
 This is implemented in the function `r(x,p)`.
diff --git a/scripts/fuzz_prep.sh b/scripts/fuzz_prep.sh
index 2c57d94..1198b7f 100755
--- a/scripts/fuzz_prep.sh
+++ b/scripts/fuzz_prep.sh
@@ -83,6 +83,11 @@
 	CFLAGS="$CFLAGS -fsanitize=address"
 fi
 
+# These are to get better instrumentation.
+export AFL_LLVM_LAF_SPLIT_SWITCHES=1
+export AFL_LLVM_LAF_TRANSFORM_COMPARES=1
+export AFL_LLVM_LAF_SPLIT_COMPARES=1
+
 # We want a debug build because asserts are counted as crashes too.
 CC="$CC" CFLAGS="$CFLAGS" ./configure.sh -gO3 -z
 
diff --git a/src/num.c b/src/num.c
index 0a59707..e45aa62 100644
--- a/src/num.c
+++ b/src/num.c
@@ -2929,14 +2929,14 @@
 #endif // BC_ENABLE_EXTRA_MATH
 
 /**
- * Converts a number from limbs with base BC_BASE_POW to base @a pow, where
- * @a pow is obase^N.
+ * Takes a number with limbs with base BC_BASE_POW and converts the limb at the
+ * given index to base @a pow, where @a pow is obase^N.
  * @param n    The number to convert.
  * @param rem  BC_BASE_POW - @a pow.
  * @param pow  The power of obase we will convert the number to.
  * @param idx  The index of the number to start converting at. Doing the
  *             conversion is O(n^2); we have to sweep through starting at the
- *             least significant limb
+ *             least significant limb.
  */
 static void
 bc_num_printFixup(BcNum* restrict n, BcBigDig rem, BcBigDig pow, size_t idx)
@@ -2998,8 +2998,8 @@
 }
 
 /**
- * Prepares a number for printing in a base that is not a divisor of
- * BC_BASE_POW. This basically converts the number from having limbs of base
+ * Prepares a number for printing in a base that does not have BC_BASE_POW as a
+ * power. This basically converts the number from having limbs of base
  * BC_BASE_POW to limbs of pow, where pow is obase^N.
  * @param n    The number to prepare for printing.
  * @param rem  The remainder of BC_BASE_POW when divided by a power of the base.
diff --git a/tests/bc/lib2.txt b/tests/bc/lib2.txt
index 0032da1..74e1256 100644
--- a/tests/bc/lib2.txt
+++ b/tests/bc/lib2.txt
@@ -1,6 +1,7 @@
 p(2, 8.0000)
 p(2, 8.0001)
 p(2, -8.0001)
+p(1024,32.1)
 r(0, 0)
 r(0, 1)
 r(0, 100)
diff --git a/tests/bc/lib2_results.txt b/tests/bc/lib2_results.txt
index f0753af..e5ddb51 100644
--- a/tests/bc/lib2_results.txt
+++ b/tests/bc/lib2_results.txt
@@ -1,6 +1,8 @@
 256.00000000000000000000
-256.01774518281640169821
+256.01774518281640171326
 .00390597924876622489
+42719740718418201647900434123391042292054090447133055398940832156444\
+39451561281100045924173873151.99999999999999999999
 0
 0
 0