Ordinary Differential Equation-2: 2nd Order ODE (RK4 Method) - Computational Physics Formulation & Simulation

Program 1

2nd Order ODE (RK4 Method)

Python 3.11 Ready

import numpy as np
import matplotlib.pyplot as plt

def f(w, t, omega_n, zeta, F_0, omega):
    x, v = w
    dx_dt = v
    dv_dt = -2 * zeta * omega_n * v - omega_n**2 * x + F_0 * np.sin(omega * t)
    return np.array([dx_dt, dv_dt])

# Parameters
omega_n = 1.0       # Natural frequency
zeta = 0.1          # Damping ratio
F_0 = 1.0           # Amplitude of forcing term
omega = 0.5         # Frequency of forcing term

t0 = 0
x0 = np.array([1, 0])   #Initial condition
t_end = 20
h = 0.1
t = np.arange(t0, t_end, h)
N = len(t)

# RK4 function
def rk4(f, x0, omega_n, zeta, F_0, omega, t, h):
    x = np.zeros((N, 2)) 
    x[0] = x0     
    for i in range(N - 1):
        k1 = f(x[i], t[i], omega_n, zeta, F_0, omega)
        k2 = f(x[i] + h / 2 * k1, t[i] + h / 2, omega_n, zeta, F_0, omega)
        k3 = f(x[i] + h / 2 * k2, t[i] + h / 2, omega_n, zeta, F_0, omega)
        k4 = f(x[i] + h * k3, t[i] + h, omega_n, zeta, F_0, omega)
        x[i + 1] = x[i] + (h / 6) * (k1 + 2 * k2 + 2 * k3 + k4)
    
    return x

sol = rk4(f, x0, omega_n, zeta, F_0, omega, t, h)

plt.plot(t, sol[:, 0], label='Position vs Time (RK4)', color='blue')
plt.plot(t, sol[:, 1], label='Velocity vs Time (RK4)', color='red')
plt.xlabel('Time (s)')
plt.ylabel('Position (m) / Velocity (m/s)')
plt.title('Forced Damped Harmonic Oscillator')
plt.legend()
plt.savefig('plot.png')
plt.show()

Click "Run Code" to execute script and render figures.

Mathematical Problem Formulation

Fourth-Order Runge-Kutta (RK4) Method

The RK4 method is used to solve the differential equation of a forced, damped harmonic oscillator represented by the following second-order ordinary differential equation:

$$ m \frac{d^2 x}{dt^2} + b \frac{dx}{dt} + kx = F_0 \sin(\omega t) $$

The equation can be rewritten in terms of the damping ratio $ \zeta $ and the natural frequency $ \omega_n $ as follows:

$$ \frac{d^2 x}{dt^2} = -2 \zeta \omega_n \frac{dx}{dt} - \omega_n^2 x + \frac{F_0}{m} \sin(\omega t) $$

where:

$$ \omega_n = \sqrt{\frac{k}{m}} $$ is the natural frequency,
$$ \zeta = \frac{b}{2\sqrt{km}} $$ is the damping ratio.

Given Parameters:

$$ \omega_n = 1.0 \quad \zeta = 0.1 \quad F_0 = 1.0 \quad \omega = 0.5 $$

Initial Conditions:

$$ x(0) = 1 \quad \frac{dx}{dt} \bigg|_{t=0} = 0$$

Theoretical Principles & Computational Physics Guide

2nd Order ODE (RK4 Method): Differential Equations & Dynamical Simulation

Dr. Alorika Chatterjee & Dr. Anirban Shaw Science & Physics Learners Python 3 / NumPy / SciPy

1. Core Principles & Theoretical Framework

Dynamical Evolution: Solves differential equations governing kinematics, planetary orbits, electrical circuits, and wave propagation.

Step-by-Step Stepping: Advances continuous state vectors forward in discrete time intervals $\Delta t$ using slope estimates.

Conservation Laws: Evaluates preservation of physical invariants (total energy, angular momentum) across long integration horizons.

2. Mathematical Formulation & Algorithmic Pipeline

Vector initial value problem and numerical state advancement:

$$ \frac{d\mathbf{y}}{dt} = \mathbf{f}(t, \mathbf{y}), \quad \mathbf{y}(t_0) = \mathbf{y}_0 \implies \mathbf{y}_{n+1} = \mathbf{y}_n + \Delta t \cdot \mathbf{\Phi}(t_n, \mathbf{y}_n, \Delta t) $$

3. Numerical Mechanics & Implementation Logic

Coupled Reduction: Convert $n$-th order ODEs into systems of $n$ coupled first-order equations.
Intermediate Slopes: Compute weighted intermediate derivatives ($k_1, k_2, k_3, k_4$ in RK4).
State Update: Advance state variables simultaneously and verify energy/phase invariants.

4. Relevant Computational Practice Exercises

Exercise 1

Solve the differential equation for 2nd Order ODE (RK4 Method) over 1000 time steps and plot state coordinates vs time.

Exercise 2

Compare numerical solution stability for small vs large step sizes $h = \Delta t$.

Exercise 3

Plot phase space trajectory demonstrating physical stability or damping attractors.

Standard Literature & References:

John R. Taylor, "Classical Mechanics", University Science Books (2005).
Mark Newman, "Computational Physics", University of Michigan (2013).

Ch. 13.4 2nd Order ODE (RK4 Method)