aboutsummaryrefslogtreecommitdiff
path: root/v_windows/v/old/vlib/math/factorial/factorial.v
diff options
context:
space:
mode:
Diffstat (limited to 'v_windows/v/old/vlib/math/factorial/factorial.v')
-rw-r--r--v_windows/v/old/vlib/math/factorial/factorial.v80
1 files changed, 80 insertions, 0 deletions
diff --git a/v_windows/v/old/vlib/math/factorial/factorial.v b/v_windows/v/old/vlib/math/factorial/factorial.v
new file mode 100644
index 0000000..9668d5d
--- /dev/null
+++ b/v_windows/v/old/vlib/math/factorial/factorial.v
@@ -0,0 +1,80 @@
+// Copyright (c) 2019-2021 Alexander Medvednikov. All rights reserved.
+// Use of this source code is governed by an MIT license
+// that can be found in the LICENSE file.
+
+// Module created by Ulises Jeremias Cornejo Fandos based on
+// the definitions provided in https://scientificc.github.io/cmathl/
+
+module factorial
+
+import math
+
+// factorial calculates the factorial of the provided value.
+pub fn factorial(n f64) f64 {
+ // For a large postive argument (n >= FACTORIALS.len) return max_f64
+
+ if n >= factorials_table.len {
+ return math.max_f64
+ }
+
+ // Otherwise return n!.
+ if n == f64(i64(n)) && n >= 0.0 {
+ return factorials_table[i64(n)]
+ }
+
+ return math.gamma(n + 1.0)
+}
+
+// log_factorial calculates the log-factorial of the provided value.
+pub fn log_factorial(n f64) f64 {
+ // For a large postive argument (n < 0) return max_f64
+
+ if n < 0 {
+ return -math.max_f64
+ }
+
+ // If n < N then return ln(n!).
+
+ if n != f64(i64(n)) {
+ return math.log_gamma(n + 1)
+ } else if n < log_factorials_table.len {
+ return log_factorials_table[i64(n)]
+ }
+
+ // Otherwise return asymptotic expansion of ln(n!).
+
+ return log_factorial_asymptotic_expansion(int(n))
+}
+
+fn log_factorial_asymptotic_expansion(n int) f64 {
+ m := 6
+ mut term := []f64{}
+ xx := f64((n + 1) * (n + 1))
+ mut xj := f64(n + 1)
+
+ log_factorial := log_sqrt_2pi - xj + (xj - 0.5) * math.log(xj)
+
+ mut i := 0
+
+ for i = 0; i < m; i++ {
+ term << b_numbers[i] / xj
+ xj *= xx
+ }
+
+ mut sum := term[m - 1]
+
+ for i = m - 2; i >= 0; i-- {
+ if math.abs(sum) <= math.abs(term[i]) {
+ break
+ }
+
+ sum = term[i]
+ }
+
+ for i >= 0 {
+ sum += term[i]
+ i--
+ }
+
+ return log_factorial + sum
+}