From 995c5533c645250c120b11f0eddc53909fc3d012 Mon Sep 17 00:00:00 2001 From: fxdup <47389903+fxdup@users.noreply.github.com> Date: Fri, 6 Oct 2023 14:46:58 -0400 Subject: [PATCH] Consolidate gamma (#9769) * refactor(gamma): Append _iterative to func name * refactor(gamma): Consolidate implementations * refactor(gamma): Redundant test function removal * Update maths/gamma.py --------- Co-authored-by: Tianyi Zheng --- maths/gamma.py | 91 ++++++++++++++++++++++++++++++++++------ maths/gamma_recursive.py | 77 ---------------------------------- 2 files changed, 79 insertions(+), 89 deletions(-) delete mode 100644 maths/gamma_recursive.py diff --git a/maths/gamma.py b/maths/gamma.py index d5debc587..822bbc744 100644 --- a/maths/gamma.py +++ b/maths/gamma.py @@ -1,35 +1,43 @@ +""" +Gamma function is a very useful tool in math and physics. +It helps calculating complex integral in a convenient way. +for more info: https://en.wikipedia.org/wiki/Gamma_function +In mathematics, the gamma function is one commonly +used extension of the factorial function to complex numbers. +The gamma function is defined for all complex numbers except +the non-positive integers +Python's Standard Library math.gamma() function overflows around gamma(171.624). +""" import math from numpy import inf from scipy.integrate import quad -def gamma(num: float) -> float: +def gamma_iterative(num: float) -> float: """ - https://en.wikipedia.org/wiki/Gamma_function - In mathematics, the gamma function is one commonly - used extension of the factorial function to complex numbers. - The gamma function is defined for all complex numbers except the non-positive - integers - >>> gamma(-1) + Calculates the value of Gamma function of num + where num is either an integer (1, 2, 3..) or a half-integer (0.5, 1.5, 2.5 ...). + + >>> gamma_iterative(-1) Traceback (most recent call last): ... ValueError: math domain error - >>> gamma(0) + >>> gamma_iterative(0) Traceback (most recent call last): ... ValueError: math domain error - >>> gamma(9) + >>> gamma_iterative(9) 40320.0 >>> from math import gamma as math_gamma - >>> all(.99999999 < gamma(i) / math_gamma(i) <= 1.000000001 + >>> all(.99999999 < gamma_iterative(i) / math_gamma(i) <= 1.000000001 ... for i in range(1, 50)) True - >>> gamma(-1)/math_gamma(-1) <= 1.000000001 + >>> gamma_iterative(-1)/math_gamma(-1) <= 1.000000001 Traceback (most recent call last): ... ValueError: math domain error - >>> gamma(3.3) - math_gamma(3.3) <= 0.00000001 + >>> gamma_iterative(3.3) - math_gamma(3.3) <= 0.00000001 True """ if num <= 0: @@ -42,7 +50,66 @@ def integrand(x: float, z: float) -> float: return math.pow(x, z - 1) * math.exp(-x) +def gamma_recursive(num: float) -> float: + """ + Calculates the value of Gamma function of num + where num is either an integer (1, 2, 3..) or a half-integer (0.5, 1.5, 2.5 ...). + Implemented using recursion + Examples: + >>> from math import isclose, gamma as math_gamma + >>> gamma_recursive(0.5) + 1.7724538509055159 + >>> gamma_recursive(1) + 1.0 + >>> gamma_recursive(2) + 1.0 + >>> gamma_recursive(3.5) + 3.3233509704478426 + >>> gamma_recursive(171.5) + 9.483367566824795e+307 + >>> all(isclose(gamma_recursive(num), math_gamma(num)) + ... for num in (0.5, 2, 3.5, 171.5)) + True + >>> gamma_recursive(0) + Traceback (most recent call last): + ... + ValueError: math domain error + >>> gamma_recursive(-1.1) + Traceback (most recent call last): + ... + ValueError: math domain error + >>> gamma_recursive(-4) + Traceback (most recent call last): + ... + ValueError: math domain error + >>> gamma_recursive(172) + Traceback (most recent call last): + ... + OverflowError: math range error + >>> gamma_recursive(1.1) + Traceback (most recent call last): + ... + NotImplementedError: num must be an integer or a half-integer + """ + if num <= 0: + raise ValueError("math domain error") + if num > 171.5: + raise OverflowError("math range error") + elif num - int(num) not in (0, 0.5): + raise NotImplementedError("num must be an integer or a half-integer") + elif num == 0.5: + return math.sqrt(math.pi) + else: + return 1.0 if num == 1 else (num - 1) * gamma_recursive(num - 1) + + if __name__ == "__main__": from doctest import testmod testmod() + num = 1.0 + while num: + num = float(input("Gamma of: ")) + print(f"gamma_iterative({num}) = {gamma_iterative(num)}") + print(f"gamma_recursive({num}) = {gamma_recursive(num)}") + print("\nEnter 0 to exit...") diff --git a/maths/gamma_recursive.py b/maths/gamma_recursive.py deleted file mode 100644 index 3d6b8c5e8..000000000 --- a/maths/gamma_recursive.py +++ /dev/null @@ -1,77 +0,0 @@ -""" -Gamma function is a very useful tool in math and physics. -It helps calculating complex integral in a convenient way. -for more info: https://en.wikipedia.org/wiki/Gamma_function -Python's Standard Library math.gamma() function overflows around gamma(171.624). -""" -from math import pi, sqrt - - -def gamma(num: float) -> float: - """ - Calculates the value of Gamma function of num - where num is either an integer (1, 2, 3..) or a half-integer (0.5, 1.5, 2.5 ...). - Implemented using recursion - Examples: - >>> from math import isclose, gamma as math_gamma - >>> gamma(0.5) - 1.7724538509055159 - >>> gamma(2) - 1.0 - >>> gamma(3.5) - 3.3233509704478426 - >>> gamma(171.5) - 9.483367566824795e+307 - >>> all(isclose(gamma(num), math_gamma(num)) for num in (0.5, 2, 3.5, 171.5)) - True - >>> gamma(0) - Traceback (most recent call last): - ... - ValueError: math domain error - >>> gamma(-1.1) - Traceback (most recent call last): - ... - ValueError: math domain error - >>> gamma(-4) - Traceback (most recent call last): - ... - ValueError: math domain error - >>> gamma(172) - Traceback (most recent call last): - ... - OverflowError: math range error - >>> gamma(1.1) - Traceback (most recent call last): - ... - NotImplementedError: num must be an integer or a half-integer - """ - if num <= 0: - raise ValueError("math domain error") - if num > 171.5: - raise OverflowError("math range error") - elif num - int(num) not in (0, 0.5): - raise NotImplementedError("num must be an integer or a half-integer") - elif num == 0.5: - return sqrt(pi) - else: - return 1.0 if num == 1 else (num - 1) * gamma(num - 1) - - -def test_gamma() -> None: - """ - >>> test_gamma() - """ - assert gamma(0.5) == sqrt(pi) - assert gamma(1) == 1.0 - assert gamma(2) == 1.0 - - -if __name__ == "__main__": - from doctest import testmod - - testmod() - num = 1.0 - while num: - num = float(input("Gamma of: ")) - print(f"gamma({num}) = {gamma(num)}") - print("\nEnter 0 to exit...")