From c71c280726fb4e9487833993042e54598fe94fd9 Mon Sep 17 00:00:00 2001 From: Ravi Kumar <119737193+ravi-ivar-7@users.noreply.github.com> Date: Thu, 26 Oct 2023 12:50:28 +0530 Subject: [PATCH] added runge kutta gills method to maths/ numerical_analysis (#10967) * added runge kutta gills method * Apply suggestions from code review --------- Co-authored-by: Tianyi Zheng --- maths/numerical_analysis/runge_kutta_gills.py | 89 +++++++++++++++++++ 1 file changed, 89 insertions(+) create mode 100644 maths/numerical_analysis/runge_kutta_gills.py diff --git a/maths/numerical_analysis/runge_kutta_gills.py b/maths/numerical_analysis/runge_kutta_gills.py new file mode 100644 index 000000000..2bd9cd612 --- /dev/null +++ b/maths/numerical_analysis/runge_kutta_gills.py @@ -0,0 +1,89 @@ +""" +Use the Runge-Kutta-Gill's method of order 4 to solve Ordinary Differential Equations. + +https://www.geeksforgeeks.org/gills-4th-order-method-to-solve-differential-equations/ +Author : Ravi Kumar +""" +from collections.abc import Callable +from math import sqrt + +import numpy as np + + +def runge_kutta_gills( + func: Callable[[float, float], float], + x_initial: float, + y_initial: float, + step_size: float, + x_final: float, +) -> np.ndarray: + """ + Solve an Ordinary Differential Equations using Runge-Kutta-Gills Method of order 4. + + args: + func: An ordinary differential equation (ODE) as function of x and y. + x_initial: The initial value of x. + y_initial: The initial value of y. + step_size: The increment value of x. + x_final: The final value of x. + + Returns: + Solution of y at each nodal point + + >>> def f(x, y): + ... return (x-y)/2 + >>> y = runge_kutta_gills(f, 0, 3, 0.2, 5) + >>> y[-1] + 3.4104259225717537 + + >>> def f(x,y): + ... return x + >>> y = runge_kutta_gills(f, -1, 0, 0.2, 0) + >>> y + array([ 0. , -0.18, -0.32, -0.42, -0.48, -0.5 ]) + + >>> def f(x, y): + ... return x + y + >>> y = runge_kutta_gills(f, 0, 0, 0.2, -1) + Traceback (most recent call last): + ... + ValueError: The final value of x must be greater than initial value of x. + + >>> def f(x, y): + ... return x + >>> y = runge_kutta_gills(f, -1, 0, -0.2, 0) + Traceback (most recent call last): + ... + ValueError: Step size must be positive. + """ + if x_initial >= x_final: + raise ValueError( + "The final value of x must be greater than initial value of x." + ) + + if step_size <= 0: + raise ValueError("Step size must be positive.") + + n = int((x_final - x_initial) / step_size) + y = np.zeros(n + 1) + y[0] = y_initial + for i in range(n): + k1 = step_size * func(x_initial, y[i]) + k2 = step_size * func(x_initial + step_size / 2, y[i] + k1 / 2) + k3 = step_size * func( + x_initial + step_size / 2, + y[i] + (-0.5 + 1 / sqrt(2)) * k1 + (1 - 1 / sqrt(2)) * k2, + ) + k4 = step_size * func( + x_initial + step_size, y[i] - (1 / sqrt(2)) * k2 + (1 + 1 / sqrt(2)) * k3 + ) + + y[i + 1] = y[i] + (k1 + (2 - sqrt(2)) * k2 + (2 + sqrt(2)) * k3 + k4) / 6 + x_initial += step_size + return y + + +if __name__ == "__main__": + import doctest + + doctest.testmod()