Ordinary Differential Equations (ODEs)

Differential equations are the language of physics and engineering, describing systems where the rate of change of a variable depends on the variable itself or independent parameters. An Ordinary Differential Equation (ODE) involves functions of a single variable and their derivatives.

1. Classification of ODEs

An $n$-th order ODE is an equation relating an independent variable $x$, a dependent variable $y(x)$, and its derivatives up to $y^{(n)}$.

Order and Linearity

• Order: The highest derivative present in the equation.
• Linearity: An ODE is linear if it can be written as: $a_n(x)y^{(n)}+a_{n-1}(x)y^{(n-1)}+\cdots +a_1(x)y^{\prime }+a_0(x)y=g(x)$ If $g(x)=0$, the equation is homogeneous; otherwise, it is non-homogeneous.

2. Existence and Uniqueness: Picard-Lindelöf

For a first-order initial value problem (IVP): $\frac{dy}{dx}=f(x,y),y(x_0)=y_0$

The Picard-Lindelöf Theorem states that if $f(x,y)$ is continuous and Lipschitz continuous in $y$ on some domain containing $(x_0,y_0)$, then there exists a unique solution $y(x)$ in some interval around $x_0$.

Failure of Lipschitz continuity often leads to non-uniqueness, such as in $y^{\prime }=y^{1/3}$ with $y(0)=0$, which has both $y(x)=0$ and $y(x)=(\frac{2}{3}x{)}^{3/2}$ as solutions.

3. First-Order Techniques

Separable Equations

If $y^{\prime }=g(x)h(y)$, we separate and integrate: $\int \frac{1}{h(y)}dy=\int g(x)dx$

Linear Equations (Integrating Factors)

For $y^{\prime }+p(x)y=q(x)$, we multiply by the integrating factor $\mu (x)=e^{\int p(x)dx}$: $\frac{d}{dx}[\mu (x)y]=\mu (x)q(x)⟹y(x)=\frac{1}{\mu (x)}\left(\int \mu (x)q(x)dx+C\right)$

Exact Equations

An equation $M(x,y)dx+N(x,y)dy=0$ is exact if $\frac{\partial M}{\partial y}=\frac{\partial N}{\partial x}$. This implies the existence of a potential function $\Psi (x,y)$ such that $\nabla \Psi =(M,N)$. The solution is $\Psi (x,y)=C$. This is a direct application of Poincaré’s Lemma in the context of differential forms.

4. Second-Order Linear Equations

Consider $a{y}^{\prime \prime }+b{y}^{\prime }+cy=g(x)$.

Homogeneous Solutions

For the homogeneous case ($g(x)=0$), we assume $y=e^{rx}$, leading to the auxiliary equation: $a{r}^2+br+c=0$

1. Distinct Real Roots ($r_1,r_2$): $y_h=c_1{e}^{r_{1}x}+c_2{e}^{r_{2}x}$
2. Repeated Roots ($r$): $y_h=c_1{e}^{rx}+c_{2}x{e}^{rx}$
3. Complex Roots ($\alpha \pm i\beta$): $y_h=e^{\alpha x}(c_1\cos \beta x+c_2\sin \beta x)$

Linear Independence and the Wronskian

Two solutions $y_1,y_2$ are linearly independent if their Wronskian is non-zero: $W(y_1,y_2)=\left|\begin{array}{cc}{y}_1& y_2\\ y_1^{\prime }& y_2^{\prime }\end{array}\right|=y_1{y}_2^{\prime }-y_2{y}_1^{\prime }\ne 0$

Non-homogeneous: Variation of Parameters

While Undetermined Coefficients works for simple $g(x)$, Variation of Parameters is universal. If $y_h=c_1{y}_1+c_2{y}_2$, then $y_p=u_1{y}_1+u_2{y}_2$ where: $u_1^{\prime }=\frac{-y_{2}g(x)}{W(y_1,y_2)},u_2^{\prime }=\frac{y_{1}g(x)}{W(y_1,y_2)}$

5. Power Series and the Method of Frobenius

When coefficients are functions (e.g., Bessel’s equation), we seek solutions as power series: $y(x)={\sum }_{n=0}^{\infty }{a}_n(x-x_0{)}^n$ For regular singular points, the Method of Frobenius assumes $y(x)=\sum a_n{x}^{n+r}$, where $r$ is determined by the indicial equation.

6. Modeling: The Harmonic Oscillator

The motion of a mass on a spring with damping and driving forces is modeled as: $m\frac{d^{2}x}{d{t}^2}+c\frac{dx}{dt}+kx=F(t)$

• Underdamped: $c^2<4mk$, oscillations decay exponentially.
• Critically Damped: $c^2=4mk$, Returns to equilibrium fastest without oscillation.
• Overdamped: $c^2>4mk$, Slow return to equilibrium.

7. Numerical Solution: The Nonlinear Pendulum

The exact equation for a pendulum is non-linear: ${\theta }^{\prime \prime }+\frac{g}{L}\sin \theta =0$. We use SciPy to solve this numerically.

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt

# Parameters: gravity g, length L
g, L = 9.81, 1.0

# System of first-order ODEs:
# Let y = [theta, omega] where omega = theta'
# y' = [omega, -g/L * sin(theta)]
def pendulum_dynamics(t, y):
    theta, omega = y
    return [omega, -(g/L) * np.sin(theta)]

# Initial conditions: 45 degrees, 0 velocity
y0 = [np.pi/4, 0.0]
t_span = (0, 10)
t_eval = np.linspace(0, 10, 500)

sol = solve_ivp(pendulum_dynamics, t_span, y0, t_eval=t_eval)

plt.plot(sol.t, sol.y[0], label='Angle (rad)')
plt.title("Nonlinear Pendulum Motion")
plt.xlabel("Time (s)")
plt.ylabel("Theta")
plt.grid(True)
plt.show()