Merge branch 'add-ode-numerical-methods' of https://github.com/sauriopqno/Python into add-ode-numerical-methods

Edit theordinary_differential_ecuations.py file

Author: sauriopqno

maths/numerical_analysis/ordinary_differential_equations.py

@@ -1,27 +1,31 @@
 """
-Numerical methods for solving Ordinary Differential Equations (ODEs) of the form y' = f(x, y).
+Numerical methods for solving Ordinary Differential Equations (ODEs)
+of the form y' = f(x, y).
 This module implements the Euler method and Heun's method.
 
 References:
 - Euler Method: https://en.wikipedia.org/wiki/Euler_method
 - Heun's Method: https://en.wikipedia.org/wiki/Heun%27s_method
-
 """
 
-from typing import Callable, List, Tuple
+from collections.abc import Callable
 
 
 def euler_method(
-    f: Callable[[float, float], float], x0: float, y0: float, h: float, steps: int
-) -> Tuple[List[float], List[float]]:
+    differential_function: Callable[[float, float], float],
+    x_initial: float,
+    y_initial: float,
+    step_size: float,
+    total_steps: int,
+) -> tuple[list[float], list[float]]:
     """
     Approximates the solution of a first-order ODE y' = f(x, y) using Euler's method.
 
-    :param f: The differential equation function f(x, y).
-    :param x0: Initial value of x.
-    :param y0: Initial value of y (y(x0) = y0).
-    :param h: Step size.
-    :param steps: Number of iterations to perform.
+    :param differential_function: The differential equation function f(x, y).
+    :param x_initial: Initial value of x.
+    :param y_initial: Initial value of y (y(x_initial) = y_initial).
+    :param step_size: Step size.
+    :param total_steps: Number of iterations to perform.
     :return: A tuple containing lists of x and y coordinates.
 
     >>> def f_ode(x: float, y: float) -> float: return y - x**2 + 1
@@ -31,56 +35,69 @@ def euler_method(
     >>> [round(i, 4) for i in y_vals]
     [0.5, 0.875, 1.3281]
     """
-    x = [0.0] * (steps + 1)
-    y = [0.0] * (steps + 1)
-    x[0], y[0] = x0, y0
+    x_estimates = [0.0] * (total_steps + 1)
+    y_estimates = [0.0] * (total_steps + 1)
+    x_estimates[0], y_estimates[0] = x_initial, y_initial
 
-    for i in range(steps):
-        y[i + 1] = y[i] + f(x[i], y[i]) * h
-        x[i + 1] = x[i] + h
+    for i in range(total_steps):
+        y_estimates[i + 1] = (
+            y_estimates[i]
+            + differential_function(x_estimates[i], y_estimates[i]) * step_size
+        )
+        x_estimates[i + 1] = x_estimates[i] + step_size
 
-    return x, y
+    return x_estimates, y_estimates
 
 
 def heuns_method(
-    f: Callable[[float, float], float], x0: float, y0: float, h: float, steps: int
-) -> Tuple[List[float], List[float]]:
+    differential_function: Callable[[float, float], float],
+    x_initial: float,
+    y_initial: float,
+    step_size: float,
+    total_steps: int,
+) -> tuple[list[float], list[float]]:
     """
-    Approximates the solution of a firs -order ODE y' = f(x, y) using Heun's method
+    Approximates the solution of a first-order ODE y' = f(x, y) using Heun's method
     (also known as the improved Euler method).
 
-    :param f: The differential equation function f(x, y).
-    :param x0: Initial value of x.
-    :param y0: Initial value of y.
-    :param h: Step size.
-    :param steps: Number of iterations to perform.
+    :param differential_function: The differential equation function f(x, y).
+    :param x_initial: Initial value of x.
+    :param y_initial: Initial value of y.
+    :param step_size: Step size.
+    :param total_steps: Number of iterations to perform.
     :return: A tuple containing lists of x and y coordinates.
 
     >>> def f_ode(x: float, y: float) -> float: return y - x**2 + 1
     >>> x_vals, y_vals = heuns_method(f_ode, 0.0, 0.5, 0.25, 2)
     >>> [round(i, 4) for i in y_vals]
     [0.5, 0.9141, 1.4114]
     """
-    x = [0.0] * (steps + 1)
-    y = [0.0] * (steps + 1)
-    x[0], y[0] = x0, y0
-
-    for i in range(steps):
-        y_predictor = y[i] + f(x[i], y[i]) * h
-        x[i + 1] = x[i] + h
-        y[i + 1] = y[i] + ((f(x[i], y[i]) + f(x[i + 1], y_predictor)) / 2.0) * h
-
-    return x, y
+    x_estimates = [0.0] * (total_steps + 1)
+    y_estimates = [0.0] * (total_steps + 1)
+    x_estimates[0], y_estimates[0] = x_initial, y_initial
+
+    for i in range(total_steps):
+        y_predictor = (
+            y_estimates[i]
+            + differential_function(x_estimates[i], y_estimates[i]) * step_size
+        )
+        x_estimates[i + 1] = x_estimates[i] + step_size
+        y_estimates[i + 1] = (
+            y_estimates[i]
+            + (
+                (
+                    differential_function(x_estimates[i], y_estimates[i])
+                    + differential_function(x_estimates[i + 1], y_predictor)
+                )
+                / 2.0
+            )
+            * step_size
+        )
+
+    return x_estimates, y_estimates
 
 
 if __name__ == "__main__":
     import doctest
 
-    doctest.testmod()
-
-    # example
-    def example_ode(x: float, y: float) -> float:
-        return y - x**2 + 1
-
-    x_res, y_res = heuns_method(example_ode, 0.0, 0.5, 0.25, 4)
-    print(f"Final heuns methot y value at x={x_res[-1]}: {round(y_res[-1], 4)}")
+    doctest.testmod()