Refactor Monte Carlo estimation functions and add doctest

- Refactored all functions to return values instead of printing output.

- Added doctests for each fucntion.

- Clarified docstrings.

- Improved code style (typing hints, naming) to align with repo conventions.

Author: akiels

maths/monte_carlo.py

@@ -8,35 +8,53 @@
 from statistics import mean
 
 
-def pi_estimator(iterations: int):
+def is_in_unit_circle(x: float, y: float) -> bool:
     """
-    An implementation of the Monte Carlo method used to find pi.
+    Checks whether the point (x, y) lies inside or on the boundary
+    of a unit circle centered at the origin.
+
+    >>> is_in_unit_circle(0.0, 0.0)
+    True
+    >>> is_in_unit_circle(1.0, 0.0)
+    True
+    >>> is_in_unit_circle(1.0, 1.0)
+    False
+    """
+    return sqrt(x**2 + y**2) <= 1
+
+
+def pi_estimator(iterations: int) -> float:
+    """
+    Estimate the value of pi using the Monte Carlo method.
+
+    The method simulates random points in the square and checks how
+    many fall within the inscribed unit circle.
+
+    Intuition:
     1. Draw a 2x2 square centred at (0,0).
-    2. Inscribe a circle within the square.
-    3. For each iteration, place a dot anywhere in the square.
-       a. Record the number of dots within the circle.
-    4. After all the dots are placed, divide the dots in the circle by the total.
-    5. Multiply this value by 4 to get your estimate of pi.
-    6. Print the estimated and numpy value of pi
+    2. Inscribe a unit circle within the square.
+    3. For each iteration, place a dot randomly within the square.
+    4. Divide the number of dots within the circle by the total number of dots.
+    5. Multiply this ratio by 4 to estimate pi.
+
+    Args:
+        iterations: Number of random points to generate.
+
+    Returns:
+        Estimated value of pi.
+
+    >>> round(pi_estimator(100000), 2) in [3.13, 3.14, 3.15]
+    True
     """
 
-    # A local function to see if a dot lands in the circle.
-    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.
-        return distance_from_centre <= 1
+    if iterations <= 0:
+        raise ValueError("Number of iterations must be positive")
 
-    # The proportion of guesses that landed in the circle
-    proportion = mean(
-        int(is_in_circle(uniform(-1.0, 1.0), uniform(-1.0, 1.0)))
+    inside_circle = sum(
+        is_in_unit_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(f"The estimated value of pi is {pi_estimate}")
-    print(f"The numpy value of pi is {pi}")
-    print(f"The total error is {abs(pi - pi_estimate)}")
+    return (inside_circle / iterations) * 4
 
 
 def area_under_curve_estimator(
@@ -46,83 +64,113 @@ def area_under_curve_estimator(
     max_value: float = 1.0,
 ) -> float:
     """
-    An implementation of the Monte Carlo method to find area under
-       a single variable non-negative real-valued continuous function,
-       say f(x), where x lies within a continuous bounded interval,
-       say [min_value, max_value], where min_value and max_value are
-       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. Evaluate f(x) at each of the drawn x values
-        c. Expected value = average of the function evaluations
-    4. Estimated value of integral = Expected value * (max_value - min_value)
-    5. Returns estimated value
+    Estimate the definite integral of a continuous, non-negative function f(x)
+    over a bounded interval [min_value, max_value] using the Monte Carlo method.
+
+    Intuition:
+    1. Randomly sample x values uniformly from [min_value, max_value].
+    2. Compute the average of f(x) at the sample points.
+       This estimates the expected value of f(x).
+    3. Multiply this average by the interval width to estimate the integral.
+
+    Args:
+        iterations: Number of random samples to draw.
+        function_to_integrate: A function f(x) to integrate over [min_value, max_value].
+        min_value: Lower bound of the interval.
+        max value: Upper bound of the interval.
+
+    Returns:
+        Estimated area under the curve (integral of f from min_value to max_value).
+
+    Raises:
+        ValueError: If iterations <= 0 or min_value >= max_value.
+
+    Example:
+    >>> def f(x): return x
+    >>> round(area_under_curve_estimator(100000, f), 2) in [0.49, 0.5, 0.51]
+    True
     """
+    if iterations <= 0:
+        raise ValueError("Number of iterations must be positive")
+    if min_value >= max_value:
+        raise ValueError("min_value must be less than max_value")
 
-    return mean(
+    samples = (
         function_to_integrate(uniform(min_value, max_value)) for _ in range(iterations)
-    ) * (max_value - min_value)
+    )
+    return mean(samples) * (max_value - min_value)
 
 
 def area_under_line_estimator_check(
     iterations: int, min_value: float = 0.0, max_value: float = 1.0
-) -> None:
+) -> tuple[float, float, float]:
     """
-    Checks estimation error for area_under_curve_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
+    Demonstrate the Monte Carlo integration method by estimating
+    the area under the line y = x over [min_value, max_value].
+
+    The exact integral of f(x) = x over [a, b] is:
+        (b**2 - a**2) / 2
+
+    This helps illustrate that the Monte Carlo method produces
+    a reasonable estimate for a simple, well-understood function.
+
+    Args:
+        iterations: Number of samples to use in estimation.
+        min_value: Lower bound of the integration interval.
+        max_value: Upper bound of the integration interval.
+
+    Returns:
+        A tuple of (estimated_value, expected_value, absolute_error).
+
+    >>> result = area_under_line_estimator_check(100000)
+    >>> round(result[0], 2) in [0.49, 0.5, 0.51]
+    True
     """
 
     def identity_function(x: float) -> float:
-        """
-        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
+    expected_value = (max_value**2 - min_value**2) / 2
+    error = abs(estimated_value - expected_value)
+    return estimated_value, expected_value, error
 
-    print("******************")
-    print(f"Estimating area under y=x where x varies from {min_value} to {max_value}")
-    print(f"Estimated value is {estimated_value}")
-    print(f"Expected value is {expected_value}")
-    print(f"Total error is {abs(estimated_value - expected_value)}")
-    print("******************")
 
-
-def pi_estimator_using_area_under_curve(iterations: int) -> None:
+def pi_estimator_using_area_under_curve(iterations: int) -> tuple[float, float, float]:
     """
-    Area under curve y = sqrt(4 - x^2) where x lies in 0 to 2 is equal to pi
+    Estimate the value of pi using Monte Carlo integration.
+
+    The area under the curve y = sqrt(4 - x**2) from x = 0 to 2 is equal to pi.
+    This function uses that fact to estimate pi.
+
+    Args:
+        iterations: Number of random samples to use for estimation.
+
+    Returns:
+        A tuple of (estimated_value, expected_value, absolute_error).
+
+    >>> result = pi_estimator_using_area_under_curve(100000)
+    >>> round(result[0], 2) in [3.13, 3.14, 3.15]
+    True
     """
 
-    def function_to_integrate(x: float) -> float:
+    def semicircle_function(x: float) -> float:
         """
         Represents semi-circle with radius 2
-        >>> [function_to_integrate(x) for x in [-2.0, 0.0, 2.0]]
+        >>> [semicircle_function(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
+        iterations, semicircle_function, 0.0, 2.0
     )
+    expected_value = pi
+    error = abs(estimated_value - expected_value)
 
-    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("******************")
+    return estimated_value, expected_value, error
 
 
 if __name__ == "__main__":