|
| 1 | +""" |
| 2 | +Gamma function is a very useful tool in math and physics. |
| 3 | +It helps calculating complex integral in a convenient way. |
| 4 | +for more info: https://en.wikipedia.org/wiki/Gamma_function |
| 5 | +
|
| 6 | +Python's Standard Library math.gamma() function overflows around gamma(171.624). |
| 7 | +""" |
| 8 | +from math import pi, sqrt |
| 9 | + |
| 10 | + |
| 11 | +def gamma(num: float) -> float: |
| 12 | + """ |
| 13 | + Calculates the value of Gamma function of num |
| 14 | + where num is either an integer (1, 2, 3..) or a half-integer (0.5, 1.5, 2.5 ...). |
| 15 | + Implemented using recursion |
| 16 | + Examples: |
| 17 | + >>> from math import isclose, gamma as math_gamma |
| 18 | + >>> gamma(0.5) |
| 19 | + 1.7724538509055159 |
| 20 | + >>> gamma(2) |
| 21 | + 1.0 |
| 22 | + >>> gamma(3.5) |
| 23 | + 3.3233509704478426 |
| 24 | + >>> gamma(171.5) |
| 25 | + 9.483367566824795e+307 |
| 26 | + >>> all(isclose(gamma(num), math_gamma(num)) for num in (0.5, 2, 3.5, 171.5)) |
| 27 | + True |
| 28 | + >>> gamma(0) |
| 29 | + Traceback (most recent call last): |
| 30 | + ... |
| 31 | + ValueError: math domain error |
| 32 | + >>> gamma(-1.1) |
| 33 | + Traceback (most recent call last): |
| 34 | + ... |
| 35 | + ValueError: math domain error |
| 36 | + >>> gamma(-4) |
| 37 | + Traceback (most recent call last): |
| 38 | + ... |
| 39 | + ValueError: math domain error |
| 40 | + >>> gamma(172) |
| 41 | + Traceback (most recent call last): |
| 42 | + ... |
| 43 | + OverflowError: math range error |
| 44 | + >>> gamma(1.1) |
| 45 | + Traceback (most recent call last): |
| 46 | + ... |
| 47 | + NotImplementedError: num must be an integer or a half-integer |
| 48 | + """ |
| 49 | + if num <= 0: |
| 50 | + raise ValueError("math domain error") |
| 51 | + if num > 171.5: |
| 52 | + raise OverflowError("math range error") |
| 53 | + elif num - int(num) not in (0, 0.5): |
| 54 | + raise NotImplementedError("num must be an integer or a half-integer") |
| 55 | + elif num == 0.5: |
| 56 | + return sqrt(pi) |
| 57 | + else: |
| 58 | + return 1.0 if num == 1 else (num - 1) * gamma(num - 1) |
| 59 | + |
| 60 | + |
| 61 | +def test_gamma() -> None: |
| 62 | + """ |
| 63 | + >>> test_gamma() |
| 64 | + """ |
| 65 | + assert gamma(0.5) == sqrt(pi) |
| 66 | + assert gamma(1) == 1.0 |
| 67 | + assert gamma(2) == 1.0 |
| 68 | + |
| 69 | + |
| 70 | +if __name__ == "__main__": |
| 71 | + from doctest import testmod |
| 72 | + |
| 73 | + testmod() |
| 74 | + num = 1 |
| 75 | + while num: |
| 76 | + num = float(input("Gamma of: ")) |
| 77 | + print(f"gamma({num}) = {gamma(num)}") |
| 78 | + print("\nEnter 0 to exit...") |
0 commit comments