Euler methods for ordinary differential equations

Euler methods are the simplest class of algorithms to numerically solve ordinary differential equations (ODE) of form y' = f(t,y), y(t_0) = y_0. Developed by Leonhard Euler in 1768 (Institutiones Calculi Integralis), they are the starting point of ODE numerical analysis.

Explicit Euler (Forward Euler)

Formula: y_{n+1} = y_n + h f(t_n, y_n) where h = (T - t_0)/N is integration step. Global convergence order 1: error = O(h). Fast and cheap but inaccurate. Unstable for stiff problems or oscillating ODEs.

Implicit Euler (Backward Euler)

Formula: y_{n+1} = y_n + h f(t_{n+1}, y_{n+1}). Requires solving implicit equation each step (fixed-point or Newton iteration). Double cost per iter, but A-stable: no limit on h for stability. Ideal for stiff ODEs.

Heun's method (Improved Euler)

Predictor-corrector: explicit Euler predictor y* = y_n + h f(t_n, y_n), trapezoidal corrector y_{n+1} = y_n + h/2 (f(t_n, y_n) + f(t_{n+1}, y*)). Order 2: error = O(h^2). Cost: 2 evaluations per step. Optimal accuracy/cost compromise.

Midpoint (Modified Euler)

Formula: k1 = f(t_n, y_n), y_{n+1} = y_n + h f(t_n + h/2, y_n + h/2 k1). Order 2 like Heun. Stable for oscillating ODEs. Variant of two-stage Runge-Kutta.

Global convergence and order

Local truncation error for explicit Euler is -h^2/2 y''(xi) for some xi. Summing N = (T-t0)/h steps, global error = O(h). For Heun and midpoint, global error = O(h^2). Classical RK4 reaches O(h^4). Tradeoff: higher order = more f evaluations per step, but larger h possible.

Stability

For y' = lambda y (test problem), explicit Euler stable only if |1 + h lambda| < 1, i.e. h < 2/|lambda| for negative real lambda. Implicit Euler A-stable: |1/(1 - h lambda)| < 1 for Re(lambda) < 0 always. For stiff ODEs (large eigenvalue ratio), implicit is essential.

Stiff problems

Stiff ODEs have very different time scales (e.g., chemical dynamics with fast and slow reactions). Explicit Euler requires tiny h for stability, even if solution varies little. Implicit doesn't have this limit. Example: Robertson equation y'_1 = -0.04 y_1 + 1e4 y_2 y_3 (stiff for very different constants).

Notable ODEs

• Exponential decay: y' = -k y, solution y = y_0 exp(-k t).
• Logistic growth: y' = r y (1 - y/K), saturation at K.
• Newton cooling: T' = -k (T - T_amb), tendency to T_amb.
• Riccati: y' = y^2 - t, no elementary solution.
• Lotka-Volterra (system): predator-prey.
• Harmonic oscillator: y'' + y = 0 (2x2 system).

Numerical comparison

For y' = -y, y(0) = 1, T = 5 (exact y(5) = e^-5 ~ 0.006738):

• Explicit Euler h=0.1 (50 steps): y_50 ~ 0.005154 (24% error)
• Explicit Euler h=0.01 (500 steps): y_500 ~ 0.006568 (2.5% error)
• Heun h=0.1: y_50 ~ 0.006740 (0.03% error)
• RK4 h=0.1: y_50 ~ 0.006738 (1e-7% error)

When to use which method

• Explicit Euler: didactic, non-stiff ODEs, small h acceptable.
• Implicit Euler: stiff ODEs (chemistry, fluids, equilibria).
• Heun: educational code, medium accuracy.
• RK4: practical standard for non-stiff ODEs.
• Adaptive (RK45, RK23): ODEs with variable speeds.

History

Leonhard Euler 1768 introduces method in "Institutiones Calculi Integralis". Cauchy 1840 rigorously proves convergence. Heun 1900 and Kutta 1901 extend to higher orders. Dahlquist 1956 introduces A-stability concept. Modern software (LSODE, scipy.integrate.solve_ivp) includes dozens of specialized methods.

Local and global error: where the order comes from

The distinction between local and global error explains why explicit Euler is called a first-order method. The local truncation error is the error committed in a single step, assuming the starting value exact: expanding the solution in a Taylor series shows that Euler neglects the term (h^2/2) y'', so its local error is proportional to h^2. The global error, accumulated over the whole interval divided into about 1/h steps, is instead proportional to h: halving the step halves the final error.

Second-order methods, such as Heun and the midpoint, also cancel the quadratic term and have local error h^3 and global error h^2: halving h reduces the error to a quarter. A practical way to estimate the error without knowing the exact solution is Richardson extrapolation, which compares the results obtained with two different steps (h and h/2) to deduce the effective order and a quantitative error estimate.

Toward Runge-Kutta: the ladder of one-step methods

The Euler methods are the first rungs of a wider hierarchy, that of one-step Runge-Kutta methods. Explicit Euler is the first-order Runge-Kutta; the midpoint and Heun's method are second-order variants that evaluate the slope at intermediate points to double the accuracy. The famous RK4, of order four, combines four evaluations per step and offers an excellent compromise between accuracy and cost, and is the workhorse of numerical ODE solving.

The next step is adaptive step control: methods such as Runge-Kutta-Fehlberg (RKF45) simultaneously compute a fourth-order and a fifth-order estimate, and from their difference deduce the local error, automatically adapting h to keep it below a threshold. Thus the solver uses large steps where the solution is smooth and small steps where it varies rapidly, optimizing the number of evaluations: an idea that descends directly from the simple Euler methods.

Updated July 2026