Exploring the Runge-Kutta method in JavaScript

Playing with differential equations via the RK4 method

Differential equations often show up in models tied to physical systems that store energy, such as mathematical models with springs, capacitors, pneumatics, etc. But finding the solution of a differential equation is not trivial, and many of them still have no known algebraic solution, which is why these mathematical models are an active field of study. Still, not having an algebraic solution doesn’t stop us from obtaining their values — iterative methods were developed for exactly that.

Before looking at some of these methods in JavaScript, below you can graphically explore differential equations solved by the Runge-Kutta method implemented here.

Use JavaScript format. The Math namespace may be omitted.

More examples to explore

The buttons above fill in the equation and parameters for four classic first-order cases:

• RL circuit: (5-4700*y)/0.0001 — current rising until it settles, the same circuit as in the text. Note the tiny step (nanoseconds): this circuit is stiff; with a large step the method genuinely blows up — try raising h to see.
• Capacitor charging: (5-y)/0.001 — capacitor voltage approaching 5V.
• Logistic growth: 0.5*y*(1-y/100) — population growing fast and saturating at 100.
• Cooling: -(y-25)/5 — a body at 90 cooling down to room temperature 25.

Play with h: too small is slow, too large loses stability and the chart explodes — you can see in practice why the step matters.

Euler’s method

The first iterative method for solving differential equations was created by the famous mathematician Euler; known precisely as Euler’s method. It makes only a pointwise estimate of the next value and continues from the last estimate. Since its estimate is weak, errors accumulate each iteration, making the values useless as they explode or oscillate. But being simple and easy to implement, Euler’s method is usually the starting example when studying these algorithms.

Like every differential equation modeling a physical system, these methods start from a user-defined initial condition, and each iteration estimates the next value, repeating the process from the previous one.

Since this site is mostly about electronics, let’s start from the RL circuit example below, the topic of a previous article, to explore some methods:

For a 5V DC source the circuit is modeled by the following differential equation

$5 = 4700 i{\left(t \right)} + 0.0001 \frac{d}{d t} i{\left(t \right)}$

Isolating $di(t) \over {dt}$ we get

${di(t) \over dt} = \frac{5 - 4700 * i(t)}{0.0001}$

Euler’s method is simple and didactic. First we need to define the function with the differential equation that will be called each iteration to get the next value. The didt(t,i) function defined below is the equation of the RL circuit above. Initial conditions come through the y0 and x0 parameters. h is the step taken each iteration — think of it as an amount of incremented time, and that matters, because in a fast system, say the RL circuit with small inductance, values happen quickly, while in a very slow model, such as weather models, time is on the order of hours.

function didt(t, i) {
return (5 - 4700 * i) / 0.0001;
}

function Euler(y0, x0, h, final, f) {
const r = [];
let x = x0;
let y = y0;

r[0] = y;
while (x < final) {
    y = y + h * f(x, y);
    x = x + h;
    r.push(y);
}
return r;
}

r = Euler(0, 0, 0.0001, 0.1, didt);

console.log(r);

The circuit’s initial conditions at the first instant (t = 0): current equals 0, so y0 = 0 and x0 = 0. Since the algorithm is iterative, each step represents an amount of time, so h is the time step size. As such a circuit settles fast, h should be on the order of microseconds or nanoseconds. The final parameter tells the algorithm when to stop. Euler’s method is purely didactic. It makes terrible next-value estimates that accumulate, returning unrealistic, oscillatory values. Indeed, in the program above the values become unrealistic and explosive within a few iterations.

Runge-Kutta method

Correction (2026-09-21): this page’s solver carried a subtle bug for almost two years — the step h was multiplied twice (into the ks and again in y +=), so everything evolved in slow motion. For smooth systems the chart even looked right, which is the worst kind of bug. Fixed in the code and examples above; if you copied the old version, replace it with this one.

Widely used, the Runge-Kutta method converges properly, doesn’t accumulate much error, and rarely oscillates. Unlike Euler’s method, estimates come from a weighted arithmetic mean of 4 values (K_n).

function didt(t, i) {
return (5 - 4700 * i) / 0.0001;
}

function RK(y0, x0, h, final, f) {
const r = [];
let x = x0;
let y = y0;
let k1, k2, k3, k4;
const h2 = h / 2;
const h6 = h / 6; // precomputed: 1 multiplication per step instead of 5

r.push(y);
while (x < final) {
    k1 = f(x, y);
    k2 = f(x + h2, y + h2 * k1);
    k3 = f(x + h2, y + h2 * k2);
    k4 = f(x + h, y + h * k3);
    y = y + h6 * (k1 + 2 * k2 + 2 * k3 + k4);
    x = x + h;

    r.push(y);
}
return r;
}

r = RK(0, 0, 0.00000001, 0.000001, didt);

console.log(r);

The array can certainly be discarded, with values used dynamically to plot a chart, or only the final value may be of interest.