estimate area under a curve defined by non-negative real-valued continuous function within a continuous interval using monte-carlo (#1785)

* estimate area under a curve defined by non-negative real-valued continuous function within a continuous interval using monte-carlo

* run black; update comments

* Use f”strings” and drop unnecessary returns

Co-authored-by: Christian Clauss <cclauss@me.com>

maths/monte_carlo.py

@@ -1,9 +1,10 @@
 """
 @author: MatteoRaso
 """
-from numpy import pi, sqrt
+from math import pi, sqrt
 from random import uniform
 from statistics import mean
+from typing import Callable
 
 
 def pi_estimator(iterations: int):
@@ -18,7 +19,7 @@ def pi_estimator(iterations: int):
     6. Print the estimated and numpy value of pi
     """
     # A local function to see if a dot lands in the circle.
-    def in_circle(x: float, y: float) -> bool:
+    def is_in_circle(x: float, y: float) -> bool:
         distance_from_centre = sqrt((x ** 2) + (y ** 2))
         # Our circle has a radius of 1, so a distance
         # greater than 1 would land outside the circle.
@@ -26,50 +27,99 @@ def pi_estimator(iterations: int):
 
     # The proportion of guesses that landed in the circle
     proportion = mean(
-        int(in_circle(uniform(-1.0, 1.0), uniform(-1.0, 1.0))) for _ in range(iterations)
+        int(is_in_circle(uniform(-1.0, 1.0), uniform(-1.0, 1.0)))
+        for _ in range(iterations)
     )
     # The ratio of the area for circle to square is pi/4.
     pi_estimate = proportion * 4
-    print("The estimated value of pi is ", pi_estimate)
+    print(f"The estimated value of pi is {pi_estimate}")
-    print("The numpy value of pi is ", pi)
+    print(f"The numpy value of pi is {pi}")
-    print("The total error is ", abs(pi - pi_estimate))
+    print(f"The total error is {abs(pi - pi_estimate)}")
 
 
-def area_under_line_estimator(iterations: int,
+def area_under_curve_estimator(
-                              min_value: float=0.0,
+    iterations: int,
-                              max_value: float=1.0) -> float:
+    function_to_integrate: Callable[[float], float],
+    min_value: float = 0.0,
+    max_value: float = 1.0,
+) -> float:
     """
     An implementation of the Monte Carlo method to find area under
-       y = x where x lies between min_value to max_value
+       a single variable non-negative real-valued continuous function, 
-    1. Let x be a uniformly distributed random variable between min_value to max_value
+       say f(x), where x lies within a continuous bounded interval, 
-    2. Expected value of x = (integration of x from min_value to max_value) / (max_value - min_value)
+       say [min_value, max_value], where min_value and max_value are 
-    3. Finding expected value of x:
+       finite numbers
+    1. Let x be a uniformly distributed random variable between min_value to 
+       max_value
+    2. Expected value of f(x) = 
+       (integrate f(x) from min_value to max_value)/(max_value - min_value)
+    3. Finding expected value of f(x):
         a. Repeatedly draw x from uniform distribution
-        b. Expected value = average of those values
+        b. Evaluate f(x) at each of the drawn x values
-    4. Actual value = (max_value^2 - min_value^2) / 2
+        c. Expected value = average of the function evaluations
+    4. Estimated value of integral = Expected value * (max_value - min_value)
     5. Returns estimated value
     """
-    return mean(uniform(min_value, max_value) for _ in range(iterations)) * (max_value - min_value)
+
+    return mean(
+        function_to_integrate(uniform(min_value, max_value)) for _ in range(iterations)
+    ) * (max_value - min_value)
 
 
-def area_under_line_estimator_check(iterations: int,
+def area_under_line_estimator_check(
-                                    min_value: float=0.0,
+    iterations: int, min_value: float = 0.0, max_value: float = 1.0
-                                    max_value: float=1.0) -> None:
+) -> None:
     """
-    Checks estimation error for area_under_line_estimator func
+    Checks estimation error for area_under_curve_estimator function
-    1. Calls "area_under_line_estimator" function
+    for f(x) = x where x lies within min_value to max_value
+    1. Calls "area_under_curve_estimator" function
     2. Compares with the expected value
     3. Prints estimated, expected and error value
     """
 
-    estimated_value = area_under_line_estimator(iterations, min_value, max_value)
+    def identity_function(x: float) -> float:
-    expected_value = (max_value*max_value - min_value*min_value) / 2
+        """
+        Represents identity function
+        >>> [function_to_integrate(x) for x in [-2.0, -1.0, 0.0, 1.0, 2.0]]
+        [-2.0, -1.0, 0.0, 1.0, 2.0]
+        """
+        return x
+
+    estimated_value = area_under_curve_estimator(
+        iterations, identity_function, min_value, max_value
+    )
+    expected_value = (max_value * max_value - min_value * min_value) / 2
 
     print("******************")
-    print("Estimating area under y=x where x varies from ",min_value, " to ",max_value)
+    print(f"Estimating area under y=x where x varies from {min_value} to {max_value}")
-    print("Estimated value is ", estimated_value)
+    print(f"Estimated value is {estimated_value}")
-    print("Expected value is ", expected_value)
+    print(f"Expected value is {expected_value}")
-    print("Total error is ", abs(estimated_value - expected_value))
+    print(f"Total error is {abs(estimated_value - expected_value)}")
+    print("******************")
+
+
+def pi_estimator_using_area_under_curve(iterations: int) -> None:
+    """
+    Area under curve y = sqrt(4 - x^2) where x lies in 0 to 2 is equal to pi
+    """
+
+    def function_to_integrate(x: float) -> float:
+        """
+        Represents semi-circle with radius 2
+        >>> [function_to_integrate(x) for x in [-2.0, 0.0, 2.0]]
+        [0.0, 2.0, 0.0]
+        """
+        return sqrt(4.0 - x * x)
+
+    estimated_value = area_under_curve_estimator(
+        iterations, function_to_integrate, 0.0, 2.0
+    )
+
+    print("******************")
+    print("Estimating pi using area_under_curve_estimator")
+    print(f"Estimated value is {estimated_value}")
+    print(f"Expected value is {pi}")
+    print(f"Total error is {abs(estimated_value - pi)}")
     print("******************")