Second-Order ODE Solver

Solve second-order linear ODEs with constant coefficients, homogeneous or not. Derives the characteristic equation, classifies the damping as over, critical or under-damped, gives the closed form and plots the phase-plane trajectory.

Free to useNo sign-up requiredInstant Results

Second-Order ODE SolverTry it now — free ▼

Second-order linear ODE a·y″ + b·y′ + c·y = g(x)

Quick examples:

Undamped SHM Underdamped Critically damped Overdamped Driven (const) Decaying force Unstable

Coefficients

a·y″

b·y′

c·y

The equation is a·y″ + b·y′ + c·y = g(x) with real constants a, b, c. Coefficient a must be non-zero.

Forcing term (right-hand side)

g(x) =

Closed-form y_p is derived for constants, polynomials (up to degree 2) and A·exp(k·x). Any other expression still gets a high-accuracy numerical solution.

Initial conditions & plot range

x₀

y₀

y′₀

x min

x max

Two initial conditions are required: y(x₀) = y₀ and y′(x₀) = y′₀. The solver integrates outward from x₀ in both x-directions.

Embed Second-Order ODE Solver Widget

About Second-Order ODE Solver

The Second-Order ODE Solver takes a linear ordinary differential equation of the form a·y″ + b·y′ + c·y = g(x) with constant real coefficients, automatically derives its characteristic equation, classifies the damping regime (overdamped, critically damped, underdamped, undamped or unstable), and produces both a symbolic closed-form solution and a high-accuracy numerical solution. The interactive output pairs a dual-curve time plot of y(x) and y′(x) with a phase-plane trajectory of (y, y′) — a view that exposes the regime at a glance: spiral-in for underdamped, node-in for overdamped, closed loop for undamped, spiral-out for unstable.

What Is a Second-Order Linear ODE With Constant Coefficients?

A second-order linear ordinary differential equation with real constant coefficients is an equation of the form

a · y″(x) + b · y′(x) + c · y(x) = g(x)

where a ≠ 0, b, c are real constants and g(x) is the forcing term. Two initial conditions y(x₀) = y₀ and y′(x₀) = y′₀ turn this into an initial-value problem with a unique solution on a neighbourhood of x₀ — this follows from the Picard-Lindelöf theorem applied to the equivalent first-order system.

If g(x) = 0 the equation is homogeneous. Otherwise it is non-homogeneous, and the full solution decomposes as

y(x) = y_h(x) + y_p(x)

where y_h is the general solution of the associated homogeneous equation (contains two free constants) and y_p is any particular solution of the full equation. Applying the two initial conditions fixes the two free constants.

The Characteristic Equation

Guessing y = e^(r·x) in the homogeneous equation gives the characteristic (or auxiliary) equation

a · r² + b · r + c = 0

a quadratic whose discriminant Δ = b² − 4ac controls the entire qualitative behaviour:

Three Root Cases & the Damping Regime

Discriminant Δ = b² − 4ac	Roots of a·r² + b·r + c = 0	Homogeneous solution y_h(x)	Physical regime
Δ > 0	Two distinct real roots r₁, r₂	C₁·e^(r₁·x) + C₂·e^(r₂·x)	Overdamped — slow monotonic return to equilibrium (both roots negative).
Δ = 0	One repeated real root r = −b/(2a)	(C₁ + C₂·x)·e^(r·x)	Critically damped — fastest return without oscillation.
Δ < 0	Complex conjugate roots α ± β·i with α = −b/(2a), β = √(−Δ)/(2a)	e^(α·x)·(C₁·cos(β·x) + C₂·sin(β·x))	Underdamped (α < 0) — decaying oscillation; undamped (α = 0) — pure sinusoidal motion; unstable (α > 0) — growing oscillation.

Method of Undetermined Coefficients (Non-homogeneous Case)

When g(x) takes one of the following simple forms, the method of undetermined coefficients supplies a particular solution by assuming a trial of the same form with unknown coefficients and solving for them:

Constant g(x) = k. Trial: y_p = K. If c = 0 multiply by x; if b = 0 also, multiply by x again.
Polynomial of degree n. Trial: general polynomial of degree n. Multiply by x or x² if the constant or linear term resonates.
Exponential g(x) = A·e^(k·x). Trial: y_p = K·e^(k·x). If k coincides with a characteristic root, multiply by x (simple root) or x² (double root) — this is resonance.
Sinusoidal g(x) = A·cos(ω·x) + B·sin(ω·x). Trial: y_p = K₁·cos(ω·x) + K₂·sin(ω·x). Multiply by x if iω is a root (pure-frequency resonance).
Products and sums follow by linearity and the product rule.

Reading the Phase Plane

The equivalent first-order system is u = y, v = y′ with u′ = v and v′ = (g(x) − b·v − c·u)/a. Plotting v against u parametrically in x gives the phase-plane trajectory. For homogeneous autonomous systems (no x in g), orbits are uniquely determined by their starting point (y₀, y′₀) and reveal the regime at a glance:

Underdamped: trajectory spirals inward toward the origin.
Overdamped: trajectory approaches origin along an invariant line (slow eigenvector).
Critically damped: degenerate node, trajectory tangent to the single eigenvector.
Undamped: closed ellipse encircling the origin — perpetual oscillation.
Unstable: trajectory spirals or runs outward to infinity.

Worked Example: Driven Damped Harmonic Oscillator

Consider the equation y″ + 2·y′ + 5·y = 10 with y(0) = 0, y′(0) = 0 — a driven, underdamped system.

Characteristic equation: r² + 2r + 5 = 0 → Δ = 4 − 20 = −16 → r = −1 ± 2i.
Homogeneous solution: y_h = e^(−x)·(C₁·cos 2x + C₂·sin 2x).
Particular solution for constant forcing g = 10: try y_p = K, so 5K = 10, giving y_p = 2.
Apply ICs: y(0) = 0 → C₁ + 2 = 0 → C₁ = −2. y′(0) = 0 → −C₁ + 2C₂ = 0 → C₂ = −1.
Final answer: y(x) = 2 − e^(−x)·(2·cos 2x + sin 2x) — oscillates with decaying envelope and limit y → 2.

How to Use This Calculator

Enter the coefficients a, b, c in the top row. a must be non-zero (otherwise the equation is first-order).
Type the forcing term g(x), or leave it as 0 for a homogeneous problem. Closed-form particular solutions are derived for constants, polynomials up to degree 2, and single exponentials A·e^(k·x) including the resonance case.
Provide initial conditions (x₀, y₀, y′₀). Both y and y′ at x₀ must be specified because the equation is second-order.
Pick the x range for the plots. The solver integrates outward from x₀ in both x-directions using RK4.
Click Solve & Visualize. You get the characteristic equation with its roots on the complex plane, the damping-regime classification, the homogeneous and particular closed-form solutions, a dual-curve time plot of y and y′, and the phase-plane trajectory.

Common Applications

Mechanical spring-mass-damper systems: m·x″ + c·x′ + k·x = F(t). Overdamped, critically damped, and underdamped correspond to different damping ratios ζ = c/(2·√(m·k)).
RLC electrical circuits: series RLC circuits obey L·Q″ + R·Q′ + Q/C = V(t) — identical structure, different symbols.
Pendulum (small-angle): θ″ + (g/L)·θ = 0 gives simple harmonic motion; adding air resistance gives damped oscillation.
Building response to earthquakes: single-degree-of-freedom structure with base acceleration as the forcing term.
PID-controlled servo systems: closed-loop error dynamics reduce to a second-order ODE whose damping ratio governs overshoot.
Population models with inertia: economic growth with capital-accumulation lag, or ecological models with delayed response.

Numerical Method — Classical Runge-Kutta (RK4) on the 2D System

The tool reduces a·y″ + b·y′ + c·y = g(x) to the first-order system

u′ = v v′ = ( g(x) − b·v − c·u ) / a

with u(x₀) = y₀, v(x₀) = y′₀. Four-stage Runge-Kutta is then applied to the vector state (u, v). RK4 has local truncation error O(h⁵) and global error O(h⁴); the default 400 sub-steps in each direction give roughly six-digit accuracy for non-stiff problems.

Frequently Asked Questions

What is a second-order linear ODE with constant coefficients?

A second-order linear ODE with constant coefficients has the form a·y″ + b·y′ + c·y = g(x), where a, b, c are real constants and g(x) is the forcing (non-homogeneous) term. With two initial conditions y(x₀) = y₀ and y′(x₀) = y′₀ the solution is unique. The homogeneous case g(x) = 0 always admits a closed-form solution via the characteristic equation a·r² + b·r + c = 0; the non-homogeneous case is solved as y(x) = y_h(x) + y_p(x).

What is the characteristic equation?

For a·y″ + b·y′ + c·y = 0, substituting the ansatz y = e^(r·x) yields a·r² + b·r + c = 0 — the characteristic or auxiliary equation. Its roots determine the form of the homogeneous solution: two distinct real roots give y_h = C₁·e^(r₁·x) + C₂·e^(r₂·x); a repeated root r gives y_h = (C₁ + C₂·x)·e^(r·x); complex conjugate roots α ± β·i give y_h = e^(α·x)·(C₁·cos(β·x) + C₂·sin(β·x)).

What does under-, critically, and overdamped mean?

The terminology comes from the spring-mass-damper model m·x″ + c·x′ + k·x = 0. Overdamped (discriminant > 0, two real roots) means the system returns to equilibrium slowly without oscillation. Critically damped (discriminant = 0, repeated root) is the fastest return without overshoot. Underdamped (discriminant < 0, complex roots) gives decaying oscillation. Undamped (b = 0, c/a > 0) gives pure sinusoidal oscillation forever.

What is the method of undetermined coefficients?

For simple forcing g(x) — constants, polynomials, exponentials, sines, cosines, and their products — the particular solution y_p is assumed to have the same form as g with unknown coefficients, which are determined by substituting into the ODE and matching terms. The trial must be multiplied by x (or x² for double roots) when g(x) resonates with a characteristic root.

What is a phase plane?

For a second-order equation reduced to the 2D system (y, y′) the phase plane plots y′ versus y as x advances. Solution curves in the phase plane reveal regime at a glance: decaying spirals for underdamped, inward nodes for overdamped, closed ellipses for undamped harmonic motion, and outward spirals for unstable oscillation. It is the geometric counterpart of the characteristic-equation root diagram.

What numerical method does this tool use?

The classical fourth-order Runge-Kutta (RK4) method is applied to the equivalent first-order system u = y, v = y′, with u′ = v and v′ = (g(x) − b·v − c·u)/a. RK4 has local truncation error O(h⁵) and the default 400 sub-steps per direction give roughly six-digit accuracy for non-stiff equations over the chosen window.