math: faster factorial function

vlib/math/math.v

@@ -7,7 +7,7 @@ module math
 
 // NOTE
 // When adding a new function, please make sure it's in the right place.
-// All functions are sorted alphabetically. 
+// All functions are sorted alphabetically.
 
 const (
 	E   = 2.71828182845904523536028747135266249775724709369995957496696763
@@ -34,13 +34,13 @@ const (
         MinI16  = -1 << 15
         MaxI32  = (1<<31) - 1
         MinI32  = -1 << 31
-//        MaxI64  = ((1<<63) - 1) 
+//        MaxI64  = ((1<<63) - 1)
-//        MinI64  = (-(1 << 63) ) 
+//        MinI64  = (-(1 << 63) )
-        MaxU8  = (1<<8) - 1 
+        MaxU8  = (1<<8) - 1
         MaxU16 = (1<<16) - 1
         MaxU32 = (1<<32) - 1
         MaxU64 = (1<<64) - 1
-) 
+)
 
 // Returns the absolute value.
 pub fn abs(a f64) f64 {
@@ -76,7 +76,7 @@ pub fn cbrt(a f64) f64 {
 }
 
 // ceil returns the nearest integer greater or equal to the provided value.
-pub fn ceil(a f64) f64 {
+pub fn ceil(a f64) int {
 	return C.ceil(a)
 }
 
@@ -131,15 +131,45 @@ pub fn exp2(a f64) f64 {
 }
 
 // factorial calculates the factorial of the provided value.
-pub fn factorial(a int) i64 {
+fn recursive_product( n int, current_number_ptr &int) int{
-	if a < 0 {
+    mut m := n / 2
-		panic('factorial: Cannot find factorial of negative number')
+    if (m == 0){
-	}
+        return *current_number_ptr += 2
-	mut prod := 1
+    }
-	for i:= 0; i < a; i++ {
+    if (n == 2){
-		prod *= (i+1)
+        return (*current_number_ptr += 2) * (*current_number_ptr += 2)
-	}
+    }
-	return prod
+    return recursive_product((n - m), *current_number_ptr) * recursive_product(m, *current_number_ptr)
+}
+
+pub fn factorial(n int) i64 {
+    if n < 0 {
+        panic('factorial: Cannot find factorial of negative number')
+    }
+    if n < 2 {
+        return i64(1)
+    }
+    mut r := 1
+    mut p := 1
+    mut current_number := 1
+    mut h := 0
+    mut shift := 0
+    mut high := 1
+    mut len := high
+    mut log2n := int(floor(log2(n)))
+    for ;h != n; {
+        shift += h
+        h = n >> log2n
+        log2n -= 1
+        len = high
+        high = (h - 1) | 1
+        len = (high - len)/2
+        if (len > 0){
+            p *= recursive_product(len, &current_number)
+            r *= p
+        }
+    }
+    return i64((r << shift))
 }
 
 // floor returns the nearest integer lower or equal of the provided value.

vlib/math/math_test.v

@@ -28,6 +28,7 @@ fn test_digits() {
 }
 
 fn test_factorial() {
+	assert math.factorial(12) == 479001600
 	assert math.factorial(5) == 120
 	assert math.factorial(0) == 1
 }