diff --git a/DIRECTORY.md b/DIRECTORY.md
index f0a34a553946..8854369894a3 100644
--- a/DIRECTORY.md
+++ b/DIRECTORY.md
@@ -719,6 +719,7 @@
   * [Pollard Rho](maths/pollard_rho.py)
   * [Polynomial Evaluation](maths/polynomial_evaluation.py)
   * Polynomials
+    * [Legendre](maths/polynomials/legendre.py)
     * [Single Indeterminate Operations](maths/polynomials/single_indeterminate_operations.py)
   * [Power Using Recursion](maths/power_using_recursion.py)
   * [Prime Check](maths/prime_check.py)
diff --git a/maths/monte_carlo.py b/maths/monte_carlo.py
index d174a0b188a2..29c25a67930c 100644
--- a/maths/monte_carlo.py
+++ b/maths/monte_carlo.py
@@ -18,6 +18,10 @@ def pi_estimator(iterations: int):
     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
+    >>> pi_estimator(1000)
+    The estimated value of pi is 3.145
+    The numpy value of pi is 3.141592653589793
+    The total error is 0.003
     """
 
     # A local function to see if a dot lands in the circle.
@@ -61,8 +65,11 @@ def area_under_curve_estimator(
         c. Expected value = average of the function evaluations
     4. Estimated value of integral = Expected value * (max_value - min_value)
     5. Returns estimated value
+    >>> def test_function(x):
+    >>>     return x * x
+    >>> area_under_curve_estimator(1000, test_function)
+    0.334 (estimated value should be close to 1/3)
     """
-
     return mean(
         function_to_integrate(uniform(min_value, max_value)) for _ in range(iterations)
     ) * (max_value - min_value)
@@ -77,12 +84,19 @@ def area_under_line_estimator_check(
     1. Calls "area_under_curve_estimator" function
     2. Compares with the expected value
     3. Prints estimated, expected and error value
+    >>> area_under_line_estimator_check(1000)
+    ******************
+    Estimating area under y=x where x varies from 0.0 to 1.0
+    Estimated value is 0.332
+    Expected value is 0.5
+    Total error is 0.168
+    ******************
     """
 
     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]]
+        >>> [identity_function(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
@@ -103,6 +117,13 @@ def identity_function(x: float) -> float:
 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
+    >>> pi_estimator_using_area_under_curve(1000)
+    ******************
+    Estimating pi using area_under_curve_estimator
+    Estimated value is 3.141
+    Expected value is 3.141592653589793
+    Total error is 0.0004
+    ******************
     """
 
     def function_to_integrate(x: float) -> float:
diff --git a/maths/polynomials/legendre.py b/maths/polynomials/legendre.py
new file mode 100644
index 000000000000..d5cbb2564c59
--- /dev/null
+++ b/maths/polynomials/legendre.py
@@ -0,0 +1,55 @@
+# Imports de bibliothèques standard
+from math import factorial
+
+# Imports de bibliothèques tierces
+import pytest
+from numpy.polynomial import Polynomial
+
+
+def legendre(n: int) -> list[float]:
+    """
+    Compute the coefficients of the nth Legendre polynomial.
+
+    The Legendre polynomials are solutions to Legendre's differential equation
+    and are widely used in physics and engineering.
+
+    Parameters:
+        n (int): The order of the Legendre polynomial.
+
+    Returns:
+        list[float]: Coefficients of the polynomial in ascending order of powers.
+    """
+    legendre_polynomial = (1 / (factorial(n) * (2**n))) * (Polynomial([-1, 0, 1]) ** n)
+    return legendre_polynomial.deriv(n).coef.tolist()
+
+
+def test_legendre_0() -> None:
+    """Test the 0th Legendre polynomial."""
+    assert legendre(0) == [1.0], "The 0th Legendre polynomial should be [1.0]"
+
+
+def test_legendre_1() -> None:
+    """Test the 1st Legendre polynomial."""
+    assert legendre(1) == [0.0, 1.0], "The 1st Legendre polynomial should be [0.0, 1.0]"
+
+
+def test_legendre_2() -> None:
+    """Test the 2nd Legendre polynomial."""
+    assert legendre(2) == [-0.5, 0.0, 1.5]
+    "The 2nd Legendre polynomial should be [-0.5, 0.0, 1.5]"
+
+
+def test_legendre_3() -> None:
+    """Test the 3rd Legendre polynomial."""
+    assert legendre(3) == [0.0, -1.5, 0.0, 2.5]
+    "The 3rd Legendre polynomial should be [0.0, -1.5, 0.0, 2.5]"
+
+
+def test_legendre_4() -> None:
+    """Test the 4th Legendre polynomial."""
+    assert legendre(4) == pytest.approx([0.375, 0.0, -3.75, 0.0, 4.375])
+    "The 4th Legendre polynomial should be [0.375, 0.0, -3.75, 0.0, 4.375]"
+
+
+if __name__ == "__main__":
+    pytest.main()