Van der Pol Oscillator

DEFINITION

ẍ − μ(1 − x²)ẋ + x = 0
drawn as the phase portrait (x, ẋ)
μ → 0 a circle · μ large a relaxation oscillation

NOTES

Damping that changes sign: inside |x| < 1 the term feeds the oscillation, outside it drains one. Whatever the starting state, the system ends on the same closed orbit — the limit cycle, and van der Pol coined "relaxation oscillation" for what it becomes as μ grows, when the motion splits into long slow crawls and sudden jumps. He got it out of a triode circuit; the same equation later turned up in heartbeats. The sweep walks μ, so the plate morphs from the near-circle to the stiff cornered loop and back.

PROVENANCE

Origin

B. van der Pol, "On ‘relaxation-oscillations’", The London, Edinburgh and Dublin Philosophical Magazine and Journal of Science, ser. 7, 2(11), 1926, 978–992

Standing

Public domain — a differential equation

Constants

μ is swept from the runner. The transient is integrated off before anything is recorded: the limit cycle is the exhibit, and the spiral onto it is not part of the figure

Source

doi:10.1080/14786442608564127

HOUDINI · VEX

The same published mathematics as a Detail Wrangle body. Paste it into a Wrangle with Run Over set to Detail; every constant is the published value plus a tweak channel, so Create Spare Parameters gives a slider that starts where the paper does.

// FORMA — PL. 76 · VAN DER POL OSCILLATOR — Balthasar van der Pol, 1926
//   ẍ − μ(1 − x²)ẋ + x = 0
//   drawn as the phase portrait (x, ẋ)
//   μ → 0 a circle · μ large a relaxation oscillation
// Paste into a Detail Wrangle (Run Over: Detail), no inputs needed.
// Written from the published mathematics, not adapted from any code.
// Constants arrive at their published values. Press the node's Create
// Spare Parameters button and every tweak becomes a slider — starting
// at 0, the published figure, and moving in the constant's own units.
// https://forma-gen.com/#plate=vanderpol

float p_mu = 1.6 + chf('mu_tweak');     // μ — nonlinearity · live 0.2 .. 6

// The plate's own colour: FORMA's ATTRACTORS accent as a cosine ramp,
// brightest near t = 0 and t = 1, near-black around t = 0.5.
vector forma_ramp(float t){
  return set(
    0.0956 + 0.1039 * cos(6.28318530718 * (t + 0)),
    0.4041 + 0.4392 * cos(6.28318530718 * (t + 0.05)),
    0.46 + 0.5 * cos(6.28318530718 * (t + 0.1)));
}

// The limit cycle in the phase plane (x, ẋ): settle 120 time units onto the
// cycle, then record from one upward crossing of ẋ = 0 to the next, so the
// loop closes on a real point of the orbit rather than wherever recording
// happened to stop. Decimated, as the plate is, so the polyline carries the
// cycle rather than the integrator's step count.
int   forma_settle = 60000;
int   forma_max    = 40000;
int   forma_n      = 900;
float forma_dt     = 0.002;

float x = 2.0, y = 0.0;
for (int i = 0; i < forma_settle; i++){
    float dy = p_mu * (1.0 - x * x) * y - x;
    x += y * forma_dt;  y += dy * forma_dt;
}
// walk to an upward crossing of ẋ = 0
for (int i = 0; i < forma_max; i++){
    float py = y, dy = p_mu * (1.0 - x * x) * y - x;
    x += y * forma_dt;  y += dy * forma_dt;
    if (py < 0.0 && y >= 0.0) break;
}
// record one full cycle, crossing to crossing
float rx[], ry[];
push(rx, x);  push(ry, y);
for (int i = 0; i < forma_max; i++){
    float py = y, dy = p_mu * (1.0 - x * x) * y - x;
    x += y * forma_dt;  y += dy * forma_dt;
    push(rx, x);  push(ry, y);
    if (i > 200 && py < 0.0 && y >= 0.0) break;
}
// the phase plane is drawn ẋ-up already; no flip needed
int last = len(rx) - 1;
int prim = addprim(0, "polyline");
for (int i = 0; i <= forma_n; i++){
    int j = min(last, int(rint(float(i) / float(forma_n) * float(last))));
    int pt = addpoint(0, set(rx[j], ry[j], 0.0));
    float u = float(i) / float(forma_n);
    setpointattrib(0, "Cd", pt, forma_ramp(0.8 + 0.3 * u));
    addvertex(0, prim, pt);
}

AFTER EFFECTS · EXPRESSION

The same published mathematics as a Shape Layer path expression. Paste it onto a Path property; every constant is the published value plus a Slider Control named _tweak, so a bare paste already draws the figure and each slider moves one constant in its own units. Trim Paths is the comet.

// FORMA — PL. 76 · VAN DER POL OSCILLATOR — Balthasar van der Pol, 1926
//   ẍ − μ(1 − x²)ẋ + x = 0
//   drawn as the phase portrait (x, ẋ)
//   μ → 0 a circle · μ large a relaxation oscillation
// After Effects port — paste onto a Shape Layer's Path property
// (Contents › Shape › Path). Written from the published mathematics, not
// adapted from any code. Constants arrive at their published values; add a
// Slider Control (Effect › Expression Controls) named <k>_tweak and that
// constant moves in its own units, starting at 0 — the published figure.
// The plate's comet and its reveal are Trim Paths; the stroke colour is
// FORMA's ATTRACTORS accent, #35E0FF. Animation runs on time.
// https://forma-gen.com/#plate=vanderpol

// A missing slider reads 0, so a bare paste already draws the figure.
function forma_tweak(n){ try { return effect(n)("Slider"); } catch (e){ return 0; } }
var p_mu = 1.6 + forma_tweak("mu_tweak");   // μ — nonlinearity · live 0.2 .. 6

// The frame: the plate's W × H canvas is this comp, with the origin at the
// layer's anchor; canvas y already runs down, as After Effects' does.
var forma_W = thisComp.width, forma_H = thisComp.height, forma_t = time;
var forma_phase = 0.2815996913705021;   // this plate's own fixed phase, as the page has it
function forma_pt(x, y){ return [x - forma_W / 2, y - forma_H / 2]; }

// The Van der Pol limit cycle: ẍ − μ(1 − x²)ẋ + x = 0 integrated by Euler at
// dt = 0.002 from (2, 0) — 60,000 steps to settle onto the cycle, then one
// full lap recorded between upward crossings of ẋ = 0 and resampled to 900
// points, scaled so the cycle fills 0.40 of the frame. The page's head laps
// every 9 s: Trim Paths, End 0 → 100%.
var dt = 0.002, x = 2, y = 0;
function forma_dydt(x, y){ return p_mu * (1 - x * x) * y - x; }
for (var i = 0; i < 60000; i++){ var d0 = forma_dydt(x, y); x += y * dt; y += d0 * dt; }
for (var j = 0; j < 40000; j++){
  var py = y, d1 = forma_dydt(x, y);
  x += y * dt; y += d1 * dt;
  if (py < 0 && y >= 0) break;
}
var raw = [[x, y]];
for (var k = 0; k < 40000; k++){
  var py2 = y, d2 = forma_dydt(x, y);
  x += y * dt; y += d2 * dt;
  raw.push([x, y]);
  if (k > 200 && py2 < 0 && y >= 0) break;
}
var N = 900, samp = [], ax = 0, ay = 0;
for (var m = 0; m <= N; m++){
  var q = raw[Math.min(raw.length - 1, Math.round(m / N * (raw.length - 1)))];
  samp.push(q);
  if (Math.abs(q[0]) > ax) ax = Math.abs(q[0]);
  if (Math.abs(q[1]) > ay) ay = Math.abs(q[1]);
}
var s = Math.min(forma_W * 0.40 / (ax || 1), forma_H * 0.40 / (ay || 1));
var pts = [];
for (var n = 0; n < samp.length; n++) pts.push(forma_pt(forma_W / 2 + samp[n][0] * s, forma_H / 2 - samp[n][1] * s));
createPath(pts, [], [], true);