Sunbelt Computer Software

PL/B Language Development and Support

Python/maths/gamma_recursive.py at master · programmin1/Python · GitHub

78 lines (70 loc) · 2.08 KB

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

"""

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

while num:

num = float(input("Gamma of: "))

print(f"gamma({num}) = {gamma(num)}")

print("\nEnter 0 to exit...")