What is a Differential Equation

Understanding differential equations, order, linearity, and mathematical modeling of continuous dynamic laws.

What it is

A differential equation (DE) is a functional relationship linking one or more independent variables, one or more unknown functions, and one or more derivatives of those functions with respect to the independent variables: $F\left(x, y(x), \frac{dy}{dx}, \frac{d^2 y}{dx^2}, \dots, \frac{d^n y}{dx^n}\right) = 0$ The order of the equation is determined by the highest derivative present in the expression. The degree corresponds to the algebraic power of that highest derivative once cleared of radicals and fractional exponents.

They are primarily classified into:

1. Ordinary Differential Equations (ODE): The unknown function depends on a single independent variable (commonly time $t$ or spatial coordinate $x$).
2. Partial Differential Equations (PDE): The unknown function depends on multiple independent variables (e.g., space and time $\psi(x,y,z,t)$), involving partial derivatives $\frac{\partial \psi}{\partial t}$, $\nabla^2 \psi$.

Furthermore, an equation is linear if the unknown function and all its derivatives appear to the first power without products between the function and its derivatives: $a_n(x) \frac{d^n y}{dx^n} + a_{n-1}(x) \frac{d^{n-1} y}{dx^{n-1}} + \dots + a_1(x) \frac{dy}{dx} + a_0(x) y = g(x)$ If the coefficients $a_i(x)$ are constants, the differential operator is translationally invariant. If $g(x) = 0$, the equation is homogeneous and the principle of linear superposition holds.

Why it appears in physics

The fundamental laws of nature rarely dictate the instantaneous absolute value of an observable; they dictate how that state changes in response to interactions and fields:

• Newton's second law does not dictate position $x(t)$, but relates acceleration (second time derivative) to net force: $m \frac{d^2 x}{dt^2} = F(x, \dot{x}, t)$.
• Newton's law of cooling relates the rate of temperature change to thermal gradient with the surroundings: $\frac{dT}{dt} = -k(T - T_{\text{env}})$.
• Maxwell-Faraday induction links the curl of the electric field to the time rate of change of magnetic flux: $\nabla \times \mathbf{E} = -\frac{\partial \mathbf{B}}{\partial t}$.
• The Schrödinger equation governs quantum state evolution via the Hamiltonian operator: $i\hbar \frac{\partial \Psi}{\partial t} = \hat{H}\Psi$.

Solving the physical problem means integrating these local instantaneous laws to reconstruct global trajectories, field profiles, or wavefunctions throughout the spacetime domain.

Visual or geometric intuition

Geometrically, a first-order ODE $\frac{dy}{dx} = f(x, y)$ defines a slope field (direction field) in the $(x, y)$ plane. At every coordinate point $(x, y)$, the function $f(x, y)$ specifies an infinitesimal line element with tangent slope $m = f(x, y)$. Solving the differential equation corresponds to tracing integral curves: smooth curves $y(x)$ whose tangent at every point matches the local slope field element. For a one-parameter family of curves, the general solution $y(x; C)$ foliates the plane. Imposing an initial condition $(x_0, y_0)$ selects the unique integral curve passing through that point.

Minimal symbolic language

• Explicit $n$-th order ODE: $\frac{d^n y}{dx^n} = f\left(x, y, y', \dots, y^{(n-1)}\right)$
• Linear differential operator $L$: $L[y] = \sum_{k=0}^n a_k(x) \frac{d^k y}{dx^k} = g(x)$
• General solution: Function family containing $n$ independent integration constants $C_1, C_2, \dots, C_n$: $y(x) = y_h(x; C_1, \dots, C_n) + y_p(x)$
• Particular solution: Specific function obtained by fixing constants using initial or boundary data.
• Singular solution: A curve satisfying the ODE everywhere but unobtainable by setting finite constant values in the general solution (often an envelope curve).

What the result means

The solution of a differential equation is a continuous state field or dynamic evolution trajectory.

1. Causal determinism: Given complete initial data (position and velocity in classical mechanics), the unique solution predicts past and future system states unambiguously.
2. Integration constants: Reflect experimental preparation and physical boundary constraints, independent of the governing dynamic law.
3. Asymptotic behaviour: Reveals whether the system reaches a stationary equilibrium ($t \to \infty$), sustains a limit-cycle oscillation, or diverges unstably.

When not to use it

• Discrete jumps without smooth derivatives: Impact mechanics with instantaneous momentum transfer require impulse formulations or difference equations.
• Non-local memory systems: Viscoelastic materials remembering entire stress histories require integro-differential equations or fractional calculus.
• Physical boundary domain violations: Mathematical solutions predicting negative absolute temperatures or faster-than-light velocities must be bounded within physical validity limits.

Operational summary

1. Classify the equation: ODE or PDE? Order? Linear or non-linear? Homogeneous or non-homogeneous?
2. Expect exactly $n$ independent boundary/initial conditions for an $n$-th order system.
3. Check for symmetries: Autonomous (missing independent variable)? Missing dependent variable?
4. Choose the appropriate analytical method (separation of variables, integrating factor, undetermined coefficients) or qualitative/numerical methods (Runge-Kutta).

Procedure: Analysis and Solution of a Differential Equation

Before calculating

Before attempting any algebraic manipulation or integration:

1. Identify the variables: Clearly distinguish between the independent variable (typically time $t$ or coordinate $x$) and the unknown dependent function (position $x(t)$, charge $q(t)$, temperature $T(t)$).
2. Determine order and linearity: Locate the highest derivative. Check whether the unknown function or any derivative appears with exponents other than 1 or multiplied together.
3. Check dimensional homogeneity: Every term in the differential equation must have identical physical dimensions. If $x$ is in meters and $t$ in seconds, $\frac{d^2 x}{dt^2}$ has units $\text{m}/\text{s}^2$.

Step 1: Identify the mathematical object

Write the equation in standard normalized canonical form.

• For first order: $\frac{dy}{dx} + P(x) y = Q(x) \quad \text{or} \quad M(x, y)\,dx + N(x, y)\,dy = 0$
• For second order linear: $y'' + p(x) y' + q(x) y = g(x)$ Record whether the coefficients are constants or variable functions.

Step 2: Choose the operation

Select the solving technique according to the recognized structure:

1. Separable variables: If $\frac{dy}{dx} = g(x) h(y)$, move all $y$-terms with $dy$ and all $x$-terms with $dx$: $\int \frac{dy}{h(y)} = \int g(x)\,dx$.
2. First-order linear: Apply the integrating factor method $\mu(x) = \exp\left(\int P(x)\,dx\right)$.
3. Exact equation: If $\frac{\partial M}{\partial y} = \frac{\partial N}{\partial x}$, find the potential function $\Phi(x, y) = C$.
4. Higher-order linear with constant coefficients: Set up the characteristic polynomial equation $a r^n + \dots + b r + c = 0$ to span the homogeneous nullspace $y_h$.
5. Inhomogeneous forcing term: Use undetermined coefficients (for polynomials, exponentials, harmonic sines/cosines) or variation of parameters (general Wronskian method).

Step 3: Conserve units and signs

• Any integrating factor $\mu(x)$ must have a dimensionless exponent; verify that $\int P(x)\,dx$ has no net units.
• Pay close attention to negative signs in restoring or dissipative forces: in friction $F = -\gamma v$, a sign error turns stable decay into an unphysical exponential explosion.

Step 4: Calculate

1. Carry out quadratures (analytical integrations).
2. Explicitly include integration constants ($+C$ for first order, $C_1, C_2$ for second order). Never append them casually at the end.
3. Solve for the dependent variable $y(x)$ algebraically if possible; otherwise express the solution as an implicit level curve $F(x, y) = C$.

Step 5: Interpret the result

1. Match initial/boundary data: Evaluate $y(x_0) = y_0$ and solve the algebraic system for the constants $C_i$.
2. Asymptotic limits: Check behaviour as $t \to \infty$ or $x \to 0$. Does it agree with physical intuition?
3. Limiting cases: If damping or forcing is set to zero, do you recover the simple harmonic oscillator or uniform linear motion?

Self-check checklist

• Does the count of integration constants equal the differential order of the equation?
• Have I substituted the computed solution back into the original ODE to verify identity?
• Are the dimensions and units of all integration constants consistent with physical variables?
• Does the solution contain any poles/singularities within the physical operational interval?
• Does the sign of the response preserve physical stability and energy conservation?

Worked Examples: Differential Equations

Step-by-step guided example

Problem

A projectile of mass $m = 2\,\text{kg}$ falls vertically through a viscous medium subject to gravity ($g = 9.8\,\text{m}/\text{s}^2$) and linear hydrodynamic drag $F_d = -\gamma v$ with drag coefficient $\gamma = 0.4\,\text{kg}/\text{s}$. Given that it starts from rest $v(0) = 0\,\text{m}/\text{s}$:

1. Classify and formulate the governing differential equation of motion.
2. Find the analytic solution for velocity $v(t)$.
3. Determine the asymptotic terminal velocity $v_L$.

Detailed solution

Step 1: Formulation and classification Applying Newton's second law along the downward vertical axis: $m \frac{dv}{dt} = m g - \gamma v \implies \frac{dv}{dt} + \frac{\gamma}{m} v = g$ Substituting numerical values: $\frac{dv}{dt} + \frac{0.4}{2} v = 9.8 \implies \frac{dv}{dt} + 0.2 v = 9.8$

• Type: Ordinary Differential Equation (ODE).
• Order: First order.
• Linearity: Linear with constant coefficients.
• Homogeneity: Inhomogeneous (constant source term $g = 9.8$).

Step 2: Solution via integrating factor Identify $P(t) = 0.2$ and $Q(t) = 9.8$. The integrating factor is: $\mu(t) = \exp\left(\int 0.2\,dt\right) = e^{0.2 t}$ Multiplying the ODE by $\mu(t)$: $e^{0.2 t} \frac{dv}{dt} + 0.2 e^{0.2 t} v = 9.8 e^{0.2 t} \implies \frac{d}{dt}\left(v(t) e^{0.2 t}\right) = 9.8 e^{0.2 t}$ Integrating with respect to $t$: $v(t) e^{0.2 t} = \int 9.8 e^{0.2 t}\,dt = \frac{9.8}{0.2} e^{0.2 t} + C = 49 e^{0.2 t} + C$ Solving for $v(t)$: $v(t) = 49 + C e^{-0.2 t}$

Step 3: Initial condition application The body starts from rest: $v(0) = 0$: $0 = 49 + C e^0 \implies C = -49$ The particular solution to the initial value problem is: $v(t) = 49\left(1 - e^{-0.2 t}\right)\,\text{m}/\text{s}$

Step 4: Asymptotic behavior and terminal velocity As $t \to \infty$, the transient exponential term $e^{-0.2 t} \to 0$: $v_L = \lim_{t \to \infty} v(t) = 49\,\text{m}/\text{s}$ Direct physical check: terminal velocity occurs when net force vanishes: $m g - \gamma v_L = 0 \implies v_L = \frac{m g}{\gamma} = \frac{2 \times 9.8}{0.4} = 49\,\text{m}/\text{s}$ Matches perfectly.

Inverse example

Problem

From laboratory measurements of an electrical transient response, the capacitor charge is found to follow: $q(t) = Q_0 e^{-t/\tau} \cos(\omega t)$ where $Q_0, \tau, \omega$ are positive physical constants. Reconstruct the second-order homogeneous linear ODE with constant coefficients that governs this system and identify the physical circuit components.

Detailed solution

Step 1: Identifying characteristic roots The observed transient waveform has the form: $q(t) = e^{\alpha t} \left(A \cos(\beta t) + B \sin(\beta t)\right)$ with $\alpha = -\frac{1}{\tau}$ and $\beta = \omega$. Therefore, the roots of the associated characteristic polynomial are the complex conjugate pair: $r_{1,2} = \alpha \pm i\beta = -\frac{1}{\tau} \pm i\omega$

Step 2: Reconstructing the characteristic equation The monic polynomial with roots $r_{1,2}$ is: $(r - r_1)(r - r_2) = \left(r + \frac{1}{\tau}\right)^2 + \omega^2 = 0 \implies r^2 + \frac{2}{\tau} r + \left(\frac{1}{\tau^2} + \omega^2\right) = 0$

Step 3: Formulating the differential equation Associating $r^k$ with the derivative $\frac{d^k q}{dt^k}$: $\frac{d^2 q}{dt^2} + \frac{2}{\tau} \frac{dq}{dt} + \left(\omega^2 + \frac{1}{\tau^2}\right) q = 0$

Step 4: Physical circuit identification The standard equation for an unforced series $RLC$ circuit is: $L \frac{d^2 q}{dt^2} + R \frac{dq}{dt} + \frac{1}{C} q = 0 \implies \frac{d^2 q}{dt^2} + \frac{R}{L} \frac{dq}{dt} + \frac{1}{L C} q = 0$ Matching coefficients: $\frac{R}{L} = \frac{2}{\tau} \implies \tau = \frac{2 L}{R}$ $\frac{1}{L C} = \omega^2 + \frac{1}{\tau^2} = \omega^2 + \left(\frac{R}{2 L}\right)^2 = \omega_0^2$ The damped frequency is $\omega = \sqrt{\frac{1}{L C} - \left(\frac{R}{2 L}\right)^2}$, proving this is an underdamped series $RLC$ resonator.

Common Errors: Differential Equations

Error 1: Confusing differential order with algebraic degree

A frequent misconception is classifying an equation such as $\left(\frac{dy}{dx}\right)^2 + y = 0$ as "second order".

• Cause: Conflating the algebraic power of a derivative with the number of times the unknown function has been differentiated.
• Correction: Order is the highest derivative present (first order here, since only $y'$ appears). Degree is the algebraic power of that highest derivative (degree 2 here).

Error 2: Appending the integration constant as an afterthought at the end

In non-linear ODEs or when integrating factor methods are applied, casually adding "$+C$" to an already isolated expression rather than at the exact quadrature step.

• Erroneous path: Writing $\frac{dy}{dx} = 2 x y \implies \ln y = x^2 \implies y = e^{x^2} + C$.
• Correct path: $\int \frac{dy}{y} = \int 2x\,dx \implies \ln |y| = x^2 + C_1 \implies y(x) = e^{x^2 + C_1} = C e^{x^2}$.
• Consequence: The appended "$+C$" function fails to satisfy the original differential equation.

Error 3: Assuming linear superposition holds for non-linear systems

Adding two known solutions $y_1$ and $y_2$ to form a general solution $y = c_1 y_1 + c_2 y_2$ in equations containing terms like $y^2$, $y y'$, or $\sin(y)$.

• Correction: The superposition principle applies exclusively to linear operators satisfying $L[c_1 y_1 + c_2 y_2] = c_1 L[y_1] + c_2 L[y_2]$. For the full non-linear pendulum $\ddot{\theta} + \frac{g}{L} \sin\theta = 0$, the sum of two solutions is not a solution.

Error 4: Losing singular/equilibrium solutions when dividing by $g(y)$

In separation of variables $\frac{dy}{dx} = f(x) g(y)$, dividing by $g(y)$ without testing values $y = y^*$ where $g(y^*) = 0$.

• Consequence: These roots represent equilibrium solutions (or singular solutions) that are omitted from the integration family. For instance, in $y' = y^{2/3}$, the constant solution $y(x) = 0$ is lost if division by zero is ignored.

Error 5: Overlooking dimensional units of coefficients and constants

Treating integration constants as pure dimensionless numbers.

• Correction: If $x(t)$ represents meters and $t$ seconds, in $\ddot{x} + \omega_0^2 x = 0$, parameter $\omega_0$ has units $\text{s}^{-1}$, and constants $C_1, C_2$ in $x(t) = C_1 \cos(\omega_0 t) + C_2 \sin(\omega_0 t)$ must have physical dimensions of meters.

Quick control rule

To verify if your solution $y(x)$ is mathematically sound: differentiate it up to the order of the ODE and substitute it into the original equation. If it does not reduce identically to $0 = 0$, an algebraic mistake was made.

Physical Applications: Differential Equations

1. Classical dynamics and celestial mechanics

Newtonian mechanics is formulated directly via second-order ordinary differential equations. For a mass in an inverse-square gravitational field: $m \frac{d^2 \mathbf{r}}{dt^2} = -\frac{G M m}{r^3} \mathbf{r}$ Expressing this vector system in polar coordinates yields Binet's equation: $\frac{d^2 u}{d\theta^2} + u = \frac{G M}{h^2}, \quad u = \frac{1}{r}$ whose general solution yields Keplerian conic sections (ellipses, parabolas, and hyperbolas), governing planets, comets, and spacecraft orbital mechanics.

2. Electrodynamics and circuit theory (Kirchhoff's laws)

Passive electrical networks comprising resistors ($R$), inductors ($L$), and capacitors ($C$) are modeled by linear ODEs derived from Kirchhoff's voltage loop rule: $L \frac{d^2 q}{dt^2} + R \frac{dq}{dt} + \frac{1}{C} q = V(t)$ where charge $q(t)$ and current $i(t) = \frac{dq}{dt}$ determine the frequency response of tuning filters, RF antennas, and electrical grid transient stability.

3. Chemical kinetics and reaction networks

In reactive systems, concentration rates of change follow the law of mass action: $\frac{d[A]}{dt} = -k_1 [A][B] + k_{-1}[C]$ In biochemistry, Michaelis-Menten enzyme kinetics and Turing reaction-diffusion systems for biological pattern formation rely on coupled nonlinear systems of differential equations.

4. Heat transport and diffusion

Transient one-dimensional thermal conduction across a solid medium is governed by Fourier's parabolic heat PDE: $\frac{\partial T}{\partial t} = \alpha \frac{\partial^2 T}{\partial x^2}$ where $\alpha = \frac{k}{\rho c_p}$ denotes thermal diffusivity. Solving this equation guides metallurgical quenching, architectural insulation, and semiconductor microchip cooling design.

5. Quantum mechanics and the Schrödinger equation

Matter wave dynamics at microscopic scales are governed by the fundamental partial differential equation of Schrödinger: $i\hbar \frac{\partial \Psi(\mathbf{r}, t)}{\partial t} = \left(-\frac{\hbar^2}{2m} \nabla^2 + V(\mathbf{r})\right) \Psi(\mathbf{r}, t)$ Separating time for stationary potentials yields the spatial eigenvalue problem $H \psi_n = E_n \psi_n$, whose discrete spectrum explains atomic orbitals and electronic energy quantization.

Interdisciplinary connection

Differential equations form the mathematical infrastructure of dynamic systems: from SIR epidemiological propagation models and Lotka-Volterra predator-prey population cycles in ecology, to the Black-Scholes stochastic PDE used for financial derivative valuation in quantitative finance.

Guided Practice: Differential Equations

Exercise 1: Radioactive decay and Carbon-14 dating

A fossil wood sample retains $25\%$ of the initial concentration of Carbon-14 ($^{14}\text{C}$) found in living organisms. The nuclear decay rate is directly proportional to the number of active nuclei: $\frac{dN}{dt} = -\lambda N$, with half-life $t_{1/2} = 5730\,\text{years}$.

1. Solve the differential equation by separation of variables subject to $N(0) = N_0$.
2. Compute the decay constant $\lambda$.
3. Determine the age of the fossil.

Solution

1. Separating variables: $\frac{dN}{N} = -\lambda dt \implies \ln N = -\lambda t + C_1 \implies N(t) = N_0 e^{-\lambda t}$
2. Using the half-life condition $N(t_{1/2}) = \frac{N_0}{2}$: $\frac{1}{2} = e^{-\lambda t_{1/2}} \implies \lambda = \frac{\ln 2}{t_{1/2}} = \frac{0.69315}{5730} \approx 1.2097 \times 10^{-4}\,\text{yr}^{-1}$
3. For $N(t) = 0.25 N_0 = (1/2)^2 N_0$: $0.25 = e^{-\lambda t} \implies t = \frac{-\ln 0.25}{\lambda} = \frac{2 \ln 2}{\lambda} = 2 t_{1/2} = 11460\,\text{years}$

Exercise 2: Reactor cooling via Newton's law

A temperature sensor in an industrial reactor reads $T(0) = 100^\circ\text{C}$ in an ambient chamber held at constant temperature $T_{\text{env}} = 20^\circ\text{C}$. After $10\,\text{minutes}$, the sensor reads $60^\circ\text{C}$.

1. Formulate and solve Newton's cooling ODE $\frac{dT}{dt} = -k(T - T_{\text{env}})$.
2. Calculate the cooling constant $k$ in $\text{min}^{-1}$.
3. How much additional time is required for the reactor to cool to $25^\circ\text{C}$?

Solution

1. Using shifted variable $\theta = T - T_{\text{env}}$: $\frac{d\theta}{dt} = -k\theta \implies T(t) = T_{\text{env}} + (T_0 - T_{\text{env}}) e^{-k t}$ With $T_0 = 100^\circ\text{C}$ and $T_{\text{env}} = 20^\circ\text{C}$: $T(t) = 20 + 80 e^{-k t}$
2. At $t = 10\,\text{min}$, $T(10) = 60^\circ\text{C}$: $60 = 20 + 80 e^{-10 k} \implies 40 = 80 e^{-10 k} \implies e^{-10 k} = 0.5$ $k = \frac{\ln 2}{10} \approx 0.06931\,\text{min}^{-1}$
3. Setting $T(t) = 25^\circ\text{C}$: $25 = 20 + 80 e^{-k t} \implies 5 = 80 e^{-k t} \implies e^{-k t} = \frac{1}{16} = 2^{-4}$ $k t = 4 \ln 2 \implies t = \frac{4 \ln 2}{k} = 40\,\text{minutes}$ The additional elapsed time beyond the 10-minute mark is $40 - 10 = 30\,\text{minutes}$.

Exercise 3: Torricelli's law for tank drainage

A cylindrical water tank with constant cross-sectional area $A_c = 1.5\,\text{m}^2$ is filled to initial height $h_0 = 4\,\text{m}$. An orifice of area $A_o = 0.01\,\text{m}^2$ is opened at the bottom. Torricelli's law gives the exit speed as $v = \sqrt{2 g h}$ ($g = 9.8\,\text{m}/\text{s}^2$).

1. Derive the differential equation governing the liquid height $h(t)$.
2. Solve the equation for $h(t)$.
3. Find the exact time required to drain the tank completely.

Solution

1. By mass conservation (volume rate of loss = exit flow rate): $A_c \frac{dh}{dt} = -A_o \sqrt{2 g h} \implies \frac{dh}{dt} = -\frac{A_o}{A_c} \sqrt{2 g} \sqrt{h}$
2. Separating variables: $\frac{dh}{\sqrt{h}} = -\alpha dt, \quad \alpha = \frac{A_o}{A_c} \sqrt{2 g} = \frac{0.01}{1.5} \sqrt{19.6} \approx 0.02951\,\text{s}^{-1}$ Integrating: $2\sqrt{h} - 2\sqrt{h_0} = -\alpha t \implies h(t) = \left(\sqrt{h_0} - \frac{\alpha}{2} t\right)^2$
3. Complete drainage occurs when $h(t_{\text{drain}}) = 0$: $t_{\text{drain}} = \frac{2 \sqrt{h_0}}{\alpha} = \frac{2 A_c \sqrt{h_0}}{A_o \sqrt{2 g}} = \frac{2 \times 1.5 \times 2}{0.01 \times 4.427} \approx 135.5\,\text{seconds} \approx 2.26\,\text{minutes}$