Compare all models on this prompt →
Double pendulum · step-3.5-flash
Simulate a chaotic double pendulum on a canvas using the correct equations of motion, drawing the arms and a colorful fading trail of the tip.
⚠ no runnable output
This one-shot run produced no runnable artifact.
agent produced no index.html
The transcript below shows what the model did — read it to see why the one shot fell short.
transcript · 2 messages
user
Simulate a chaotic double pendulum on a canvas using the correct equations of motion, drawing the arms and a colorful fading trail of the tip.
assistant
thinking
The user wants me to create a single-page web app that simulates a chaotic double pendulum on a canvas. They want:
1. Correct equations of motion for a double pendulum
2. Drawing the arms
3. A colorful fading trail of the tip
I need to write this all into one self-contained file named index.html with CSS in a <style> tag and JS in a <script> tag. No external resources, CDNs, or network requests.
Let me think about what's needed:
1. A canvas element
2. CSS for styling
3. JavaScript that:
- Sets up the canvas
- Implements the physics of a double pendulum (Lagrangian mechanics)
- Animates the pendulum
- Draws the arms (lines)
- Draws a trail that fades out over time
- Makes it colorful
The double pendulum equations of motion are derived from Lagrangian mechanics. The system has two masses m1 and m2 attached by rigid rods of lengths l1 and l2. The angles θ1 and θ2 are measured from the vertical.
The equations are:
For θ1:
d/dt(∂L/∂θ1̇) - ∂L/∂θ1 = 0
For θ2:
d/dt(∂L/∂θ2̇) - ∂L/∂θ2 = 0
Where the Lagrangian L = T - V (kinetic minus potential energy).
The standard derived equations (from many sources) are:
Let's define:
- m1, m2: masses
- l1, l2: lengths
- g: gravity
- θ1, θ2: angles
- ω1 = θ1̇, ω2 = θ2̇: angular velocities
Then the accelerations α1 = θ1̈, α2 = θ2̈ are:
I need to recall or derive these. A common form is:
First compute some helper variables:
delta = θ2 - θ1
den1 = (m1 + m2)*l1 - m2*l1*cos(delta)*cos(delta)
den2 = (l2/l1) * den1
α1 = (m2*l1*ω1*ω1*sin(delta)*cos(delta) +
m2*g*sin(θ2)*cos(delta) +
m2*l2*ω2*ω2*sin(delta) -
(m1 + m2)*g*sin(θ1)) / den1
α2 = (-m2*l2*ω2*ω2*sin(delta)*cos(delta) +
(m1 + m2)*(g*sin(θ1)*cos(delta) +
l1*ω1*ω1*sin(delta) -
g*sin(θ2))) / den2
But I'm not fully sure about the signs. Let me recall a more standard form from physics:
The equations are often given as:
ω1̇ = ( -g*(2*m1 + m2)*sin(θ1) - m2*g*sin(θ1-2*θ2) - 2*sin(delta)*m2*(ω2*ω2*l2 + ω1*ω1*l1*cos(delta)) ) / ( l1*(2*m1 + m2 - m2*cos(2*delta)) )
ω2̇ = ( 2*sin(delta)*(ω1*ω1*l1*(m1 + m2) + g*(m1 + m2)*cos(θ1) + ω2*ω2*l2*m2*cos(delta)) ) / ( l2*(2*m1 + m2 - m2*cos(2*delta)) )
Yes, that's a common form. Let's check:
Let:
A = 2*m1 + m2
B = m2
C = 2*sin(delta)
D = ω1*ω1*l1
E = g*cos(θ1)
F = ω2*ω2*l2*cos(delta)
G = l2
Denom1 = l1*(A - B*cos(2*delta))
ω1̇ = ( -g*A*sin(θ1) - B*g*sin(θ1-2*θ2) - C*(D + F) ) / Denom1
Denom2 = l2*(A - B*cos(2*delta))
ω2̇ = ( 2*sin(delta)*(D*(m1+m2) + g*(m1+m2)*cos(θ1) + ω2*ω2*l2*m2*cos(delta)) ) / Denom2
Simplify: D*(m1+m2) = l1*(m1+m2)*ω1*ω1? Actually D = ω1*ω1*l1. So D*(m1+m2) = ω1*ω1*l1*(m1+m2). So the expression is:
2*sin(delta)*( l1*(m1+m2)*ω1*ω1 + g*(m1+m2)*cos(θ1) + ω2*ω2*l2*m2*cos(delta) ) / (l2*(2*m1+m2 - m2*cos(2*delta)))
Yes that seems right.
We need to implement this in JS. We'll use a time step and integrate using something like Runge-Kutta or just Euler with small steps. Since it's chaotic, small errors accumulate, but for visualization it's fine.
We'll maintain state:
- θ1, θ2 (angles in radians)
- ω1, ω2 (angular velocities)
- positions of the two masses and the tip
We'll also maintain a trail: an array of points (x,y) with timestamps or age, and we'll draw them with fading opacity and possibly color based on something (like HSV cycling). The trail should be colorful and fade out over time.
We'll use requestAnimationFrame for animation.
We need to handle canvas size and scaling. Typically we set a scale factor to map pendulum coordinates (in meters) to pixels.
We'll also need to set initial conditions. Maybe random or fixed that shows chaos.
Let's plan the HTML:
<!DOCTYPE html>
<html>
<head>
<title>Double Pendulum</title>
<style>
body { margin: 0; overflow: hidden; background: #000; }
canvas { display: block; }
</style>
</head>
<body>
<canvas id="canvas"></canvas>
<script>
// JavaScript code here
</script>
</body>
</html>
Now the JavaScript:
We'll set up canvas context, handle resizing.
Define constants:
- g = 9.81 (or maybe adjust for visual speed)
- m1, m2 = 1, 1 (or maybe different)
- l1, l2 = 100, 100 (pixel lengths, but careful: the equations use actual lengths in meters. We'll need to decide on a time scale. We can incorporate a scaling factor or just treat l1, l2 as "simulation lengths" and then scale for drawing.
Better: Use a dt (time step) that corresponds to real seconds. The equations assume SI units. If we set l1 = l2 = 1 (meter), then the natural frequencies are sqrt(g/l). With g=9.81, the period is ~2 seconds. That's fine. We'll draw at a scale of e.g., 200 pixels per meter.
So we define:
const g = 9.81;
const m1 = 1, m2 = 1;
const l1 = 1, l2 = 1; // in meters
Then position of mass1: (x1, y1) = (origin_x + l1*sin(θ1), origin_y + l1*cos(θ1))? Actually typical: y positive downward? We can set coordinate system with y increasing downward for convenience, so potential energy = m*g*y. Usually we set θ measured from vertical, with θ=0 when hanging down. Then:
x1 = l1 * sin(θ1)
y1 = l1 * cos(θ1)
Then for mass2:
x2 = x1 + l2 * sin(θ2)
y2 = y1 + l2 * cos(θ2)
But careful: If we use this, both angles are measured from vertical downward? Actually if θ=0 when pendulum hanging down (vertical), then cos(0)=1, sin(0)=0, so mass1 at (0, l1) relative to pivot. That's fine. The potential energy V = m1*g*y1 + m2*g*y2, with y increasing downward.
But sometimes we define θ from the vertical upward. Let's be consistent: use the standard equations that assume θ measured from vertical, with positive in one direction. The equations I recalled assume small angle approximations? No, they are exact. They should work for any θ.
Let's adopt:
θ1, θ2: angles from the vertical (downward), positive to the right (clockwise)? Actually the sign of sin and cos matters. In the equations, sin(θ) and cos(θ) appear. If we define θ such that when hanging straight down, θ=0. Then x = l*sin(θ), y = l*cos(θ). That gives positive x to the right, positive y downward.
In the equations, typically they have terms like -g*sin(θ) if θ measured from vertical and θ positive in the direction of displacement. Actually if θ=0 at bottom, the gravitational torque is restoring: torque = -m*g*l*sin(θ). That matches a term -g*sin(θ). So the equations with -g*A*sin(θ1) are consistent with that sign convention.
Thus we'll use:
θ1, θ2 initial maybe random or fixed.
We'll integrate using a small dt, say dt = 0.01 seconds. We'll compute derivatives and update.
We can use a simple Euler method: ω += α*dt; θ += ω*dt. But for better stability, we could use a semi-implicit Euler (velocity Verlet) or RK4. But for visual, Euler might suffice if dt small. We'll start with Euler and maybe switch to RK4 if too unstable.
But the system is chaotic and sensitive. Let's try to implement a more stable integrator like RK4. But that's more code. Alternatively, we can use a smaller dt and maybe add some damping? Not needed.
Given it's a visual simulation, we can try with a moderate dt like 0.016 (60 fps assuming 1 frame = 1 step? Actually we want to simulate physics at a higher rate than rendering? We can do fixed time step for physics, then interpolate. But simpler: each frame we do one or multiple physics steps. We could do e.g., 5 physics steps per frame with dt=0.002 to keep it stable. That's easier.
Let's plan: In each animation frame, we compute several physics iterations (substeps) to keep stable. For example, const SUB_STEPS = 5; and dt = 0.016/SUB_STEPS? Actually we want total elapsed per frame approx 16ms. So we can set a fixed physics dt = 0.01, and do as many steps as needed to catch up. Or simpler: we do a fixed number of steps per frame with a small dt.
Let's do this: const PHYSICS_DT = 0.01; // seconds per physics step
In each animation frame (approx 16ms), we perform N steps where N = frameDelta / PHYSICS_DT, rounding up, but to avoid spiral of death, we limit max steps. Or we do fixed steps: per frame, do maybe 5 steps with dt=0.003. That's about 0.015 seconds per frame.
I'll choose: per frame, do 4 physics steps with dt = 0.005. That gives effective 200Hz physics, stable.
Now the trail: we want a colorful fading trail of the tip (mass2). We'll store positions of mass2 over time. We can limit the trail length (max points). We'll draw each segment with color based on something like hue cycling over time, and opacity decreasing with age.
Simplify: store points with a lifetime or age. In each frame, we add the current tip position to an array. For each existing point, we increment its age. Remove those older than maxAge. Then draw as a path with varying stroke color along the path. But we want a fading effect: we can draw many small circles or line segments with varying opacity. Better: draw a continuous path but with gradient? Canvas supports creating a gradient along a path, but that's complex. Simpler: we draw many circles with small radius and decreasing size/opacity. That might be heavy if many points.
Alternative: we draw the trail as a polyline with a single color that fades along its length using a gradient stroke. But gradient along stroke is possible with LinearGradient but that's based on coordinate, not along path. There's no built-in path gradient in 2D canvas. However we can simulate by drawing many short line segments each with its own color and opacity. That's doable: for each consecutive pair of points, draw a line segment with color that fades based on its index in the trail (or age). That is O(n) per frame, with trail length maybe 500-1000 points. That's okay.
We'll store the trail as an array of {x, y, age} where age increments each step or frame. Or we could just store the points and use index to determine opacity. If we remove old points, we can set max length, say 500 points. Then we draw from oldest to newest? Actually we want newest to be bright, oldest to fade out. So we can draw segments in order of older to newer. For each segment between point i and i+1, we compute relative age = (i / trail.length) or (age / maxAge). Let's use age: each point has an age (in frames or seconds). When adding, set age=0. Then each frame we increment all ages. Remove if age > maxAge (e.g., 2 seconds). Then for drawing, opacity = 1 - age/maxAge. Color: we can use HSL with hue cycling over time, maybe based on angle or velocity? Or simply rainbow: hue = (some function of index or time) mod 360. To make it colorful, we can change hue over time of each point: assign hue based on the time when it was added (or a global counter). So each point gets a hue based on its creation's timestamp. Then as it ages, its color remains same but fades in opacity. That yields a rainbow trail.
Implementation:
- global variable trail = []
- global lastTime for hue cycling? Or we can compute hue based on frame count: hue = (frameCount * 2) % 360. When adding a point, store its hue.
- Each frame: add new point at tip position with current hue.
- Increment age of each point by dt (or 1 per frame? We'll use frame-based aging if we want to limit by seconds, we need dt). We can store age in frames or seconds. Simpler: store age in seconds, add dt each physics update? But we add point once per frame, so age should increment per frame. We can do: at each frame, after adding, loop through trail and age += frameDeltaTime (or just 1 per frame if we assume 60fps? But actual frame time varies). We'll compute deltaTime from the animation loop.
Actually easier: we store points as {x, y, hue, age: 0} and in each frame after adding, we increase age by frameDeltaTime (in seconds). Then we remove if age > maxAge (e.g., 2 seconds). Then for drawing, opacity = 1 - age / maxAge.
But careful: frameDeltaTime is time between frames. However, our physics steps are multiple per frame, but trail update only once per frame (we want continuous line, not too dense). So we add point once per frame at the tip position after physics updates. That's fine.
Now drawing:
- Clear canvas (with some opacity to create a fading background effect? The user asked for "colorful fading trail of the tip", not necessarily a fading background. We can clear fully each frame to avoid smear. But the trail itself should fade over time. So we'll clear completely and redraw all trail points each frame.
- Draw the pendulum arms: two lines from pivot to mass1, mass1 to mass2. Use a color (white or bright).
- Draw masses: maybe circles at the ends.
- Draw the trail: iterate over trail points. For each point, we can draw a small circle with fillStyle = `hsla(${hue}, 100%, 50%, ${opacity})`. But that would be many circles; better to draw connected lines to make smooth trail? But circles might be okay if radius small. However, thousands of circles could be heavy. Let's approximate: draw line segments between consecutive points, with stroke style based on the average hue and opacity of that segment. Since we have stored hue for each point, we can interpolate? Actually we want each segment to have the color of the older point? Or the newer? To have proper gradient, we could vary color along segment, but that's expensive. Simpler: for each segment, draw a line with strokeStyle = `hsla(${hue}, 100%, 50%, ${opacity})`. But which hue? The segment is between point i and i+1. They might have different hues if hue changes over time. We could use the hue of the older point or the newer; a gradual change might be smoother if we draw many short segments and the hue changes a little each frame. So we can set the segment's color to the hue of the point i (or i+1) and it will create a gradient effect if we draw in order. That's acceptable.
But to reduce number of draw calls, we could batch segments with same hue? Not needed.
Alternatively, we could use a Path2D and set stroke style with a gradient, but gradient would be linear, not along path. So stick with segment by segment.
However, we might compute the trail polyline as one path and then use a gradient that is linear but maybe we can map the gradient to the path's bounding box? That would not give a proper fade from old to new; it would be spatial. Since the trail may loop around spatially, that wouldn't work. So we need segment-level coloring.
We'll do:
if (trail.length > 1) {
for (let i = 0; i < trail.length - 1; i++) {
const p1 = trail[i];
const p2 = trail[i+1];
const opacity = 1 - p1.age / maxAge;
if (opacity <= 0) continue;
ctx.strokeStyle = `hsla(${p1.hue}, 100%, 50%, ${opacity})`;
ctx.lineWidth = 2; // maybe vary width?
ctx.beginPath();
ctx.moveTo(p1.x, p1.y);
ctx.lineTo(p2.x, p2.y);
ctx.stroke();
}
}
That's fine.
Now physics implementation:
We need to compute α1 and α2 given the current angles and angular velocities.
Let's define function computeDerivs(state) that returns [dθ1/dt=ω1, dθ2/dt=ω2, dω1/dt=α1, dω2/dt=α2].
State: [θ1, θ2, ω1, ω2].
Compute:
sin1 = Math.sin(θ1), cos1 = Math.cos(θ1)
sin2 = Math.sin(θ2), cos2 = Math.cos(θ2)
sinDelta = Math.sin(θ2 - θ1)
cosDelta = Math.cos(θ2 - θ1)
Compute denominator common: den = 2*m1 + m2 - m2*Math.cos(2*delta) ??? Actually from the equations: Denom1 = l1*(2*m1 + m2 - m2*cos(2*delta)). And Denom2 = l2*(2*m1 + m2 - m2*cos(2*delta)). So denominator factor D = (2*m1 + m2 - m2*cos(2*delta)).
But careful: 2*delta is 2*(θ2-θ1). So cos(2*delta) = Math.cos(2*delta). So D = 2*m1 + m2 - m2*Math.cos(2*delta).
Then:
ω1̇ numerator:
- g*(2*m1 + m2)*sin(θ1)
- m2*g*Math.sin(θ1 - 2*θ2) ??? Actually the term is m2*g*sin(θ1 - 2*θ2)? Wait: the term is -m2*g*sin(θ1 - 2*θ2) but I wrote sin(θ1-2*θ2) but earlier I wrote sin(θ1-2*θ2)? Actually it's sin(θ1-2*θ2). Let's check: The term is -m2*g*sin(θ1-2*θ2). But some sources have sin(θ1 - 2*θ2), others have sin(θ2? I need to verify.
Standard equations (from Wikipedia, "Double pendulum"):
Let δ = θ2 - θ1
Then:
θ1̈ = ( -g*(2*m1+m2)*sin θ1 - m2*g*sin(θ1-2θ2) - 2*sin δ * m2*(θ2̇²*l2 + θ1̇²*l1*cos δ) ) / ( l1*(2*m1 + m2 - m2*cos(2δ)) )
θ2̈ = ( 2*sin δ * (θ1̇²*l1*(m1+m2) + g*(m1+m2)*cos θ1 + θ2̇²*l2*m2*cos δ) ) / ( l2*(2*m1 + m2 - m2*cos(2δ)) )
Yes that matches.
Check numerator for θ2̈: includes +g*(m1+m2)*cos θ1? Actually it's g*(m1+m2)*cos θ1? Wait: the term is 2*sin δ * ( θ1̇²*l1*(m1+m2) + g*(m1+m2)*cos θ1 + θ2̇²*l2*m2*cos δ ). That is positive in numerator. That matches a common form.
So let's code exactly:
let delta = theta2 - theta1;
let sinDelta = Math.sin(delta);
let cosDelta = Math.cos(delta);
let sin2delta = Math.sin(2*delta);
let cos2delta = Math.cos(2*delta);
let denom1 = l1 * (2*m1 + m2 - m2*cos2delta);
let denom2 = l2 * (2*m1 + m2 - m2*cos2delta);
// But careful: denom1 and denom2 are separate but same factor times l1 or l2.
let num1 = -g*(2*m1 + m2)*Math.sin(theta1)
- m2*g*Math.sin(theta1 - 2*theta2)
- 2*sinDelta * (omega2*omega2*l2 + omega1*omega1*l1*cosDelta);
let alpha1 = num1 / denom1;
let num2 = 2*sinDelta * ( omega1*omega1*l1*(m1+m2)
+ g*(m1+m2)*Math.cos(theta1)
+ omega2*omega2*l2*m2*cosDelta );
let alpha2 = num2 / denom2;
Note: Some sources have the term (θ2̇²*l2*m2*cos δ) inside the parentheses, but then it's multiplied by 2*sin δ. That's what we have.
But we should double-check signs: In the theta2 double dot numerator, some texts have -g*(m1+m2)*cos θ1? Let's derive quickly from Lagrangian to be sure. But given the known equations from many sources, the above is standard.
We'll trust it.
Now we need to be careful about potential division by zero. The denominator could become zero if cos(2δ) = (2*m1+m2)/m2? That's >1 for typical m1,m2 positive? 2*m1+m2 / m2 = 2*m1/m2 + 1, which can be >1, so it's possible denominator zero? Actually cos(2δ) is between -1 and 1. The expression 2*m1 + m2 - m2*cos(2δ) is minimized when cos(2δ)=1 => 2*m1 + m2 - m2*1 = 2*m1, >0. Maximized when cos(2δ)=-1 => 2*m1 + m2 + m2 = 2*m1+2*m2 >0. So denominator always positive if m1,m2>0. Good.
Now we have α1, α2.
Integration:
Given dt small, use either:
- Euler: ω += α*dt; θ += ω*dt.
But that can cause energy drift. Better: semi-implicit Euler: first compute α based on current state, then update ω: ω += α*dt; then update θ: θ += ω*dt. That's okay.
Or use RK4 for better accuracy. Let's try RK4? That would be more code but might be more stable for chaotic system. However, with small dt and substeps, Euler might be okay. We'll start with Euler, if we see instability we can improve.
But I recall that double pendulum simulation using Euler can quickly explode due to numerical errors. So maybe we should implement RK4 or at least a midpoint method. However, given the constraints, we can try Euler with enough substeps. But the user expects "correct equations of motion". Euler is not perfect but it's an approximation. It's acceptable for a visual simulation as long as it's stable.
Alternatively, we can use velocity Verlet: store previous accelerations? But I think we can do RK4.
Let's weigh: Code complexity: Euler is 4 lines. RK4 is more but doable. Since we need to produce a single file with no external deps, and we want it to look good, I'd like to use a more stable integrator. But maybe it's overkill. Many online double pendulum simulations use Euler with small dt and it's fine. I'll try Euler with dt = 0.002 and do maybe 5 sub-steps per frame (total dt=0.01 per frame). That might be stable.
We'll see.
Now the trail drawing: we want "colorful fading trail". We'll use HSL with varying hue over time. We'll increment a global hue counter each frame and assign to new point.
Now overall structure:
HTML with canvas.
JS:
const canvas = document.getElementById('canvas');
const ctx = canvas.getContext('2d');
let width, height;
function resize() {
width = canvas.width = window.innerWidth;
height = canvas.height = window.innerHeight;
}
window.addEventListener('resize', resize);
resize();
Simulation parameters:
const g = 9.81;
const m1 = 1, m2 = 1;
const l1 = 1, l2 = 1; // in meters (for simulation)
const pixelsPerMeter = 200; // scale
const pivotX = width/2;
const pivotY = height/2; // actually we need to offset? We'll set pivot at center.
But note: y increases downward. So position formulas:
x1 = pivotX + l1 * Math.sin(θ1) * ppm
y1 = pivotY + l1 * Math.cos(θ1) * ppm
x2 = x1 + l2 * Math.sin(θ2) * ppm
y2 = y1 + l2 * Math.cos(θ2) * ppm
Now state variables:
let theta1 = Math.PI/2; // start somewhere
let theta2 = Math.PI/2;
let omega1 = 0;
let omega2 = 0;
We can randomize: theta1 = Math.random()*2*Math.PI; etc.
Trail:
let trail = [];
const maxTrailAge = 2.0; // seconds, aging based on frame delta? Actually we'll compute age in seconds. We'll add age increment = frameDeltaTime.
let lastTime = 0;
But careful: We add a point each frame (after physics). So we need to track time to compute delta for aging.
In animation loop:
function animate(time) {
const dt = (time - lastTime) / 1000; // seconds
lastTime = time;
// physics substeps
const subSteps = 5;
const subDt = 0.005; // or maybe dt? Actually we want consistent physics regardless of frame rate. So we should use a fixed physics step and perform as many as needed to catch up to real time. But simpler: we'll use fixed subSteps with fixed subDt, so the simulation speed is independent of frame rate? Actually that would make simulation run faster on faster machines if we just do fixed number steps per frame. But we also have dt for aging the trail. So we need to ensure physics time increments match real time. If we do 5 steps of 0.005 each, that's 0.025 seconds of simulation per frame. If the frame rate is 60fps (~0.0167s per frame), we are simulating more time than real, causing speed-up. Not good.
Better: Use a fixed time step for physics and accumulate time. Standard approach:
let accumulator = 0;
const fixedDt = 0.01; // 10ms physics
In each frame:
accumulator += dt;
while (accumulator >= fixedDt) {
stepPhysics(fixedDt);
accumulator -= fixedDt;
}
Then add trail point after all steps? Or add after each step? We want trail at each rendered frame, so we add point once per frame after physics updates. That point represents the tip after all physics steps for that frame.
That is more robust.
Let's implement:
let lastTime = performance.now();
let accumulator = 0;
const fixedDt = 0.01; // seconds
function stepPhysics(dt) {
// compute derivatives and update state using Euler or RK4
// We'll implement Euler:
let { alpha1, alpha2 } = computeAlpha(theta1, theta2, omega1, omega2);
omega1 += alpha1 * dt;
omega2 += alpha2 * dt;
theta1 += omega1 * dt;
theta2 += omega2 * dt;
// Normalize angles to keep them in [-π, π] to avoid large numbers? Not necessary but can help with precision. But we might want to keep them as is.
// Could wrap: theta1 = (theta1 + Math.PI) % (2*Math.PI) - Math.PI;
}
But computeAlpha function returns alpha1, alpha2.
We'll implement computeAlpha.
Now in animate:
function animate(time) {
const deltaTime = (time - lastTime) / 1000;
lastTime = time;
accumulator += deltaTime;
while (accumulator >= fixedDt) {
stepPhysics(fixedDt);
accumulator -= fixedDt;
}
// after physics steps, update trail
// add current tip position
const hue = (frameCount * 2) % 360; // or increment by some amount
trail.push({x: tipX, y: tipY, hue, age: 0});
// age all points by deltaTime
for (let p of trail) {
p.age += deltaTime;
}
// remove old
trail = trail.filter(p => p.age < maxTrailAge);
// render
draw();
frameCount++;
requestAnimationFrame(animate);
}
But careful: The physics steps we do inside the while loop: we are using fixedDt, so we simulate fixed increments. However, the tip position after the while loop is after the last fixed step. That's fine.
Now drawing:
Clear canvas: ctx.fillStyle = 'rgba(0,0,0,1)'; ctx.fillRect(0,0,width,height); Or use clearRect.
We'll draw:
- Arms: lines from pivot to mass1, mass1 to mass2. Color: maybe white or light.
- Masses: circles.
- Trail: as described.
Potential improvements: We could also draw a small circle at pivot.
Now the code for computeAlpha:
function computeAlpha(theta1, theta2, omega1, omega2) {
const m1_const = m1, m2_const = m2; // avoid globals? but we'll use global constants.
const g = 9.81;
const l1_const = l1, l2_const = l2;
const delta = theta2 - theta1;
const sinDelta = Math.sin(delta);
const cosDelta = Math.cos(delta);
const sin2delta = Math.sin(2*delta);
const cos2delta = Math.cos(2*delta);
const denomFactor = 2*m1_const + m2_const - m2_const*cos2delta;
const denom1 = l1_const * denomFactor;
const denom2 = l2_const * denomFactor;
const term1 = -g*(2*m1_const + m2_const)*Math.sin(theta1);
const term2 = -m2_const*g*Math.sin(theta1 - 2*theta2);
const term3 = -2*sinDelta * (omega2*omega2*l2_const + omega1*omega1*l1_const*cosDelta);
const num1 = term1 + term2 + term3;
const alpha1 = num1 / denom1;
const term4 = 2*sinDelta * ( omega1*omega1*l1_const*(m1_const+m2_const)
+ g*(m1_const+m2_const)*Math.cos(theta1)
+ omega2*omega2*l2_const*m2_const*cosDelta );
const alpha2 = term4 / denom2;
return {alpha1, alpha2};
}
Check term2: sin(theta1 - 2*theta2) = Math.sin(theta1 - 2*theta2). That's fine.
But note: Some derivations have sin(θ1 - 2θ2) and also the denominator uses cos(2δ) = cos(2(θ2-θ1)) = cos(2θ2-2θ1). That's fine.
But there might be a sign error? Let's verify using known special cases: when θ2=0 (second pendulum hanging below first), then δ = -θ1. Then the equations should reduce to a single pendulum? Not exactly. It's complicated. We'll trust known formula.
However, I've seen another formulation:
Let's check another common source:
https://www.myphysicslab.com/pendulum/double-pendulum-en.html
They use:
dθ1/dt = ω1
dθ2/dt = ω2
dω1/dt = (-g*(2*m1+m2)*sin(θ1) - m2*g*sin(θ1-2*θ2) - 2*sin(θ2-θ1)*m2*(ω2²*l2 + ω1²*l1*cos(θ2-θ1))) / (l1*(2*m1+m2 - m2*cos(2*(θ2-θ1))))
dω2/dt = (2*sin(θ2-θ1)*(ω1²*l1*(m1+m2) + g*(m1+m2)*cos(θ1) + ω2²*l2*m2*cos(θ2-θ1))) / (l2*(2*m1+m2 - m2*cos(2*(θ2-θ1))))
That matches ours except the term3 includes ω2²*l2 + ω1²*l1*cos(delta) and multiplied by 2*sin(delta). Yes.
So it's correct.
Now about the term for θ1: they have -2*sin(θ2-θ1)*m2*(ω2²*l2 + ω1²*l1*cos(θ2-θ1)). Since δ = θ2-θ1, that's -2*sinδ * m2*(...). That's exactly our term3: -2*sinDelta * (omega2*omega2*l2 + omega1*omega1*l1*cosDelta) multiplied by m2? Wait we missed factor m2? In term3 we have -2*sinDelta * (omega2*omega2*l2 + omega1*omega1*l1*cosDelta). But the equation I just wrote includes a factor m2. Let's check the numerator for dω1/dt:
= (-g*(2m1+m2)*sinθ1 - m2*g*sin(θ1-2θ2) - 2*sin(θ2-θ1)*m2*(ω2²*l2 + ω1²*l1*cos(θ2-θ1))) / (...)
So indeed there is an m2 factor multiplying the third term. Our term3 lacks that m2. So we need to multiply term3 by m2_const.
Let's correct:
Term1: -g*(2m1+m2)*sinθ1
Term2: -m2*g*sin(θ1-2θ2)
Term3: -2*sinδ * m2 * (ω2²*l2 + ω1²*l1*cosδ)
Our earlier derived formula had: - 2*sin(δ)*m2*(ω2²*l2 + ω1²*l1*cosδ). Yes we omitted m2. So we need to add m2.
Thus:
const term3 = -2 * sinDelta * m2_const * (omega2*omega2*l2_const + omega1*omega1*l1_const*cosDelta);
Now term4 numerator:
2*sinδ * ( ω1²*l1*(m1+m2) + g*(m1+m2)*cosθ1 + ω2²*l2*m2*cosδ )
That's correct. It has m2 inside the parentheses with ω2²*l2*cosδ.
So our term4 is: 2*sinDelta * ( omega1*omega1*l1_const*(m1_const+m2_const) + g*(m1_const+m2_const)*Math.cos(theta1) + omega2*omega2*l2_const*m2_const*cosDelta );
Yes.
So code:
const { alpha1, alpha2 } = computeAlpha(...);
Now let's ensure we handle large angles and that the integration works.
Now trail point addition:
We need to compute tipX, tipY each frame after physics steps. We'll compute them using the current state.
function getPositions() {
const x1 = pivotX + l1 * Math.sin(theta1) * pixelsPerMeter;
const y1 = pivotY + l1 * Math.cos(theta1) * pixelsPerMeter;
const x2 = x1 + l2 * Math.sin(theta2) * pixelsPerMeter;
const y2 = y1 + l2 * Math.cos(theta2) * pixelsPerMeter;
return {x1, y1, x2, y2};
}
Then tip is (x2,y2).
Now we need to set pivotX, pivotY based on canvas center. We'll update on resize.
Now we need to handle the trail aging: we add deltaTime (the real time elapsed) to each point's age. This deltaTime is the time between frames, not the physics time. That's fine because we want trail to fade in real time.
Potential issue: if physics substepping makes simulation time advance more than real time? Actually we use accumulator to simulate physics with fixedDt steps that sum to the real elapsed time. So after the while loop, the simulation time advanced by deltaTime (because we consumed exactly that much from accumulator). So physics and real time are synchronized. So the tip position corresponds to real time now.
Thus we use deltaTime (in seconds) to age the trail.
Now we need to ensure we add a point only once per frame. That's fine.
Now we need to set initial state. Let's pick a non-trivial initial condition to see chaos. For example: θ1 = Math.PI/2, θ2 = Math.PI/2, ω1=0, ω2=0. That's not very chaotic because both are at same angle, they'll swing in sync? Actually they will be in phase? Might be regular. To get chaos, we need different initial angles or velocities. For instance θ1 = Math.PI/2, θ2 = Math.PI/2 + 0.1, or give ω1 a small push. Or random: theta1 = Math.random()*Math.PI; etc. We'll do random.
Let's do:
theta1 = Math.PI/2 + (Math.random() - 0.5) * 0.5;
theta2 = Math.PI/2 + (Math.random() - 0.5) * 0.5;
Or we could have a "start" button? Not needed. Just start.
Now we need to consider performance: trail length could be up to maxTrailAge / frameAge? If we add per frame at 60fps, with maxAge=2 seconds, we have 120 points. That's fine. Even if we have 500, okay.
But if we want longer trail, we could increase maxTrailAge. 2 seconds is fine.
Now drawing: We'll draw the trail with line segments. But line segments between points might be too fine if points are sparse? But we add each frame, so points are roughly at screen resolution if 60fps. The tip moves quickly, so distance between points could be large. But line segments will still be okay because they connect them directly. However, if the tip moves fast, the trail might appear as a continuous line because points are close enough? With 60fps, typical speed: tip speed maybe few m/s, scaled to 200 pixels/m => say 2m/s => 400 pixels/s => 6.67 pixels per frame at 60fps. That's okay, line segments will be about 6.7 pixels each, which is fine.
But if dt is 0.01, tip speed maybe 4 m/s? Actually the tip can go quite fast. But still okay.
Now we need to draw arms. We'll draw lines with lineWidth 2 or 3, color white.
We might also draw circles at masses: radius maybe 5-10 pixels.
Let's add:
ctx.strokeStyle = '#fff';
ctx.lineWidth = 3;
ctx.beginPath();
ctx.moveTo(pivotX, pivotY);
ctx.lineTo(x1, y1);
ctx.lineTo(x2, y2);
ctx.stroke();
Then draw circles:
ctx.fillStyle = '#fff';
ctx.beginPath();
ctx.arc(x1, y1, 5, 0, Math.PI*2);
ctx.fill();
ctx.beginPath();
ctx.arc(x2, y2, 5, 0, Math.PI*2);
ctx.fill();
We could also add a small circle at pivot.
Now the trail drawing: as described.
Make sure to set globalCompositeOperation maybe to add glow? Not needed.
Now we should also consider the possibility of numerical instability: The angles might become very large or NaN if denominator zero? We'll add a check if denom1 or denom2 is 0, but unlikely. Also, if values become infinite, we could reset simulation. But not likely with small dt.
Now let's write the complete code.
But before final, let's test mentally: If we start with random angles, we should see chaotic motion. However, the trail will be colorful and fade.
Also note: The equations use meters and g=9.81. If we set l1=l2=1, then natural period ~2s, suitable. With dt=0.01, Euler might be okay but could drift energy. We'll implement RK4 to be safe? Since it's a code generation, we want it to work well. But RK4 is more code but doable. However, the user asked for "correct equations of motion", which we are using. The integration method is separate; we could use a more accurate method. I'll consider using RK4 to ensure stability and accuracy, making the simulation more impressive.
Given that we have to write it all in one file, and code length is not a constraint, I can implement a generic RK4 integrator. But simpler: we can use the "acceleration" function that returns derivatives.
Define function acceleration(state) returns [dθ1/dt=ω1, dθ2/dt=ω2, dω1/dt=α1, dω2/dt=α2].
Then RK4:
k1 = f(state)
k2 = f(state + 0.5*dt*k1)
k3 = f(state + 0.5*dt*k2)
k4 = f(state + dt*k3)
newState = state + (dt/6)*(k1 + 2*k2 + 2*k3 + k4)
We'll need to handle arrays.
We'll maintain state as [θ1, θ2, ω1, ω2].
Implement:
function derivs(state) {
const [theta1, theta2, omega1, omega2] = state;
const {alpha1, alpha2} = computeAlpha(theta1, theta2, omega1, omega2);
return [omega1, omega2, alpha1, alpha2];
}
function stepRK4(dt) {
const k1 = derivs(state);
const k2 = derivs(state.map((s,i) => s + dt*0.5*k1[i]));
const k3 = derivs(state.map((s,i) => s + dt*0.5*k2[i]));
const k4 = derivs(state.map((s,i) => s + dt*k3[i]));
for (let i=0; i<4; i++) {
state[i] += (dt/6)*(k1[i] + 2*k2[i] + 2*k3[i] + k4[i]);
}
}
That's a bit heavy per step but with small number of points it's fine. We'll do RK4 with fixedDt=0.02 maybe? RK4 can handle larger dt. But we can still use fixedDt=0.02 and it should be stable and accurate. That would reduce substepping.
But we need to be careful: RK4 with dt=0.02 might be okay. However, double pendulum is sensitive; I'd probably use dt=0.01 or 0.02 with RK4. Let's try dt=0.02.
But we have an accumulator that sums real time and then does multiple RK4 steps until caught up. So per frame we might do 1-2 steps if using dt=0.02 and frame time ~0.016, then sometimes 2 steps (if frame took longer). That's fine.
Let's adopt RK4 for better accuracy. However, the code will be longer. But it's acceptable.
Alternatively, we could use a semi-implicit Euler with small dt and it's probably fine. Many online simulations use Euler. I'll go with RK4 because it's more correct.
But note: The user said "simulate a chaotic double pendulum on a canvas using the correct equations of motion". That implies using the correct physical equations, not necessarily a perfect integrator. So we could use Euler and it's okay.
Given simplicity, I'm leaning toward Euler with fixedDt=0.005 and 4 substeps per frame to keep simulation time close to real. But to be safe against instability, maybe I'll go with RK4.
Let's decide: I'll implement RK4 to ensure stability and accuracy, making the simulation look good for longer. It will also showcase skill.
Thus, state as array, and stepPhysics(dt) will call stepRK4(dt) if we choose RK4.
But then we need to ensure we handle the accumulation properly. We'll do:
let state = [theta1, theta2, omega1, omega2];
function stepPhysics(dt) {
// RK4 integration
const k1 = derivs(state);
const k2 = derivs(state.map((s,i) => s + 0.5*dt*k1[i]));
const k3 = derivs(state.map((s,i) => s + 0.5*dt*k2[i]));
const k4 = derivs(state.map((s,i) => s + dt*k3[i]));
for (let i=0; i<4; i++) {
state[i] += (dt/6)*(k1[i] + 2*k2[i] + 2*k3[i] + k4[i]);
}
}
But we have separate variables? We can keep state array and then unpack when needed for drawing. That is fine.
So initially:
let state = [theta1, theta2, omega1, omega2];
After integration, we have new values.
Now drawing: we need to get theta1 = state[0], theta2 = state[1].
That's fine.
But careful: computeAlpha uses the global constants m1, m2, l1, l2, g. So we can keep those.
Now we need to ensure that when we compute positions, we use state[0], state[1].
Now the color trail: hue = (frameCount * 1) % 360; We can also use something like hue = (frameCount * 2) % 360 to cycle quickly.
Now we need to handle resize: update width, height, and recompute pivotX, pivotY. That's it.
Now we need to handle the starting time: we set lastTime = performance.now(); accumulator = 0; frameCount=0;
Now we need to define maxTrailAge = 2.0; // seconds
Now we also need to add a point each frame after physics steps.
But note: When we first start, we should add a point at initial position. So before the loop or inside the first iteration we add a point. In our animate function, after the physics steps (which on first frame will be the initial state? Actually accumulator initially 0, deltaTime from time 0 to first frame maybe nonzero. We'll add that time to accumulator and then step. Then add point. That's fine.
But we must ensure that trail includes initial position. That's okay.
Now we need to also update pivotX, pivotY on resize.
Now we also need to think about performance: The trail drawing with individual line segments might be fine for 100-200 segments. But if we set maxTrailAge larger, say 5 seconds, that's 300 segments. Still fine.
Now we should also consider that we might want the trail to be smooth. If we draw line segments with lineCap='round' and lineJoin='round', it might appear smoother. We can set ctx.lineCap = 'round'; ctx.lineJoin = 'round'; before drawing trail.
We'll do:
ctx.lineCap = 'round';
ctx.lineJoin = 'round';
Now let's write the entire HTML with script.
We'll also add some styling to make canvas fullscreen.
Potential issues: The simulation might be too fast or too slow. We can adjust g or l. But we set g=9.81, l=1, dt=0.02 for RK4. The natural period is about 2 seconds for a single pendulum, so double pendulum will be faster but still okay.
But we might want to make it slower to see the motion. We could reduce g to e.g., 5? But that would not be correct. Maybe we can increase lengths? Actually we can treat l1 and l2 as 1.5 meters? That would lower frequency. But the visual scale is pixelsPerMeter; we can adjust to make arms longer. But we want to see the whole pendulum. If l1=l2=1, and pivot at center, the tip can go about 2 meters away from pivot (max distance = l1+l2 = 2m) => 400 pixels radius. That's fine.
But if canvas is smaller, it might go off-screen. We'll set ppm to something like Math.min(width, height)/4? Actually pivot at center, and max reach = (l1+l2)*ppm. If we set ppm=200, max reach = 400 pixels. That fits on most screens. But if user resizes to small window, it might clip. We could adjust ppm dynamically: ppm = min(width, height) / 4; That ensures it fits. Let's do that: let ppm = Math.min(width, height) * 0.4; Actually 0.4 of half? Let's compute: half of min dimension is half height/width. The maximum distance from pivot to tip is l1+l2 = 2. So we need ppm * 2 <= min(width,height)/2? Actually pivot at center, so distance from center to edge = min(width,height)/2. So we need 2*ppm ≤ min(width,height)/2? That would be too small. We want the pendulum to be visible within the canvas. We can set ppm = (min(width,height)/2 - margin) / (l1+l2). If we want a margin of say 20 pixels. So let's compute:
function updateScale() {
const minDim = Math.min(width, height);
const maxReach = l1 + l2; // 2 meters
ppm = (minDim/2 - 20) / maxReach;
}
But we also want to ensure not too small. But okay.
Better: We'll set ppm = width/4 maybe. But let's keep it simple: pixelsPerMeter = 200 constant. If the canvas is smaller than 800, it may clip. But the user will likely view on decent size. Since it's a demo, we can adapt. I'll implement dynamic scaling: compute ppm based on min(width,height) each frame? Actually it's constant unless resized. So in resize, we set ppm accordingly.
Let's set:
const pivotX = width/2;
const pivotY = height/2;
pixelsPerMeter = Math.min(width, height) / 4; // because max arm length sum = 2, so 2*ppm = half of min dim? That gives tip up to half min dim. Actually if minDim=800, ppm=200, max reach=400, pivot at center, so distance to edge is 400, exactly fits. Good.
So in resize(): width=canvas.width=window.innerWidth; height=canvas.height=window.innerHeight; pivotX = width/2; pivotY = height/2; pixelsPerMeter = Math.min(width, height) / 4;
That ensures it fits.
Now the trail: we add point each frame. But if we add at the very beginning (frame 0) before any physics? Actually we will add after physics steps. But initially state is set. So first frame after accumulating deltaTime (which will be maybe a few ms), we will step physics (maybe one step if accumulator>=fixedDt). Then add point. That's fine.
But we also might want to add the initial point before starting the loop so that there's a point at time zero. We could add it at start. In animate, we add after physics. So before requestAnimationFrame is called, we might not have any point. That's fine; after first draw, there will be a point.
Now let's write the full HTML.
We'll structure:
<!DOCTYPE html>
<html>
<head>
<meta charset="UTF-8">
<title>Double Pendulum Simulation</title>
<style>
* { margin: 0; padding: 0; box-sizing: border-box; }
body { overflow: hidden; background: #000; }
canvas { display: block; }
</style>
</head>
<body>
<canvas id="canvas"></canvas>
<script>
// all JS code
</script>
</body>
</html>
Now JS:
- Get canvas and context.
- Set width/height and pivot, ppm.
- Constants: g, m1, m2, l1, l2.
- State: array [θ1, θ2, ω1, ω2]. Initialize with random or fixed values.
- Trail array, maxTrailAge, frameCount, lastTime, accumulator, fixedDt.
- Functions:
computeAlpha(state) returns alphas. But actually computeAlpha uses the angles and omegas separately. We can have a function that takes theta1, theta2, omega1, omega2.
We'll have derivative function that takes state array and returns derivatives array.
stepPhysics(dt): integrate using RK4 (or Euler? We'll do RK4). I'll implement RK4.
- getPositions(state): returns {x1,y1,x2,y2} in pixels.
- draw(): clear canvas, draw arms, masses, trail.
- animate(time): compute delta, accumulate, while loop stepPhysics(fixedDt), after loop, add trail point, age points, filter, draw, requestAnimationFrame.
Now I'll write the code accordingly.
But before writing, let's check if RK4 is okay with the fact that fixedDt may be larger than we'd like for RK4? Typically RK4 works well with moderate dt. For chaotic systems, still need small dt. But we can set fixedDt=0.01 or 0.02. I'll start with 0.01. That means per frame (at 60fps, delta ~0.016), we will do either 1 or 2 steps. That's fine.
Alternatively, we could do multiple smaller steps per accumulator consumption to reduce error. But we'll just do while (accumulator >= fixedDt) and subtract. If the frame rate is high, accumulator may be less than fixedDt and we skip physics? That would cause physics to run slower than real time. We need to ensure we step physics with fixedDt until accumulator is negative? Actually the typical approach:
accumulator += deltaTime;
while (accumulator >= fixedDt) {
stepPhysics(fixedDt);
accumulator -= fixedDt;
}
But if accumulator is less than fixedDt, we don't step. That means physics time advances in discrete jumps of fixedDt each frame, but if frames are fast, we might step only occasionally. That results in simulation time lagging real time by less than fixedDt. That's acceptable; it's the standard variable time step with fixed physics.
But if we want to add a trail point each frame exactly at real time, the physics state might not be exactly at the real time but up to fixedDt behind. That's fine.
But note: if we do not step because accumulator < fixedDt, the state remains same, but we add a trail point based on that state, which is okay. The simulation will eventually catch up.
But there is a risk of spiral of death if the accumulator builds up when frame rate drops. We can limit max steps per frame, but unlikely needed.
Alternatively, we could have a more nuanced approach: we could have the trail point added based on the actual simulated time, not per frame. But simpler to add per frame after physics catch-up.
Now I'll code.
One more thing: We should use high DPI? Not necessary.
Now the trail: We'll store points with {x, y, hue, age}. We'll compute hue based on a global hueCycle that increments each frame. But we want smooth color variation, so we could increment by a small amount each frame, not full 360. We'll do: let globalHue = 0; each frame, globalHue = (globalHue + 0.5) % 360; Then assign that hue to new point. That will give a slowly cycling rainbow.
Now drawing the trail: we'll loop and draw line segments. But to improve performance, we could use a Path2D and draw with stroke once per segment? That's still many draw calls. But we can batch segments with same hue? Not needed. Instead, we can draw the entire trail as one path but with gradient? Not possible. We'll do per segment.
But we can reduce number of segments by not drawing segments for very old points that have low opacity? We already filter.
Alternatively, we could draw circles for each point: but bars? I'll stick with line segments.
Now implement.
Potential pitfalls: In the computeAlpha, using global m1, m2, l1, l2, g. We'll define them as const.
But in the derivative function, we need to use these constants. We'll just compute inside.
We must be careful about the order of parameters. The derivative function expects state = [θ1, θ2, ω1, ω2].
Now we need to convert meters to pixels: In getPositions, we use pivotX, pivotY, pixelsPerMeter.
Now I'll write the code.
I'll include comments for clarity, but not too many.
Now final file: index.html.
Let's craft.
But before finalizing, consider whether RK4 might be overkill and more bug-prone. I'll test the Euler version mentally for a few steps? Might be okay but could explode. I recall double pendulum is sensitive to integration errors; even Euler with small dt can cause energy gain and chaos. But for short time it's okay. RK4 is better.
Given we have time, I'll implement RK4.
Now we need to ensure that the state array is mutated correctly.
Implement:
function derivs(state) {
const [theta1, theta2, omega1, omega2] = state;
const {alpha1, alpha2} = computeAlpha(theta1, theta2, omega1, omega2);
return [omega1, omega2, alpha1, alpha2];
}
But computeAlpha might be called with those values.
Now stepPhysics using RK4:
function stepPhysics(dt) {
const k1 = derivs(state);
const k2 = derivs(state.map((s,i) => s + dt*0.5*k1[i]));
const k3 = derivs(state.map((s,i) => s + dt*0.5*k2[i]));
const k4 = derivs(state.map((s,i) => s + dt*k3[i]));
for (let i = 0; i < state.length; i++) {
state[i] += (dt/6)*(k1[i] + 2*k2[i] + 2*k3[i] + k4[i]);
}
}
We need to ensure we don't modify state inside derivs? No.
Now computeAlpha:
function computeAlpha(theta1, theta2, omega1, omega2) {
const delta = theta2 - theta1;
const sinDelta = Math.sin(delta);
const cosDelta = Math.cos(delta);
const sin2delta = Math.sin(2*delta);
const cos2delta = Math.cos(2*delta);
const denomFactor = 2*m1 + m2 - m2*cos2delta;
const denom1 = l1 * denomFactor;
const denom2 = l2 * denomFactor;
const term1 = -g*(2*m1 + m2)*Math.sin(theta1);
const term2 = -m2*g*Math.sin(theta1 - 2*theta2);
const term3 = -2 * sinDelta * m2 * (omega2*omega2*l2 + omega1*omega1*l1*cosDelta);
const num1 = term1 + term2 + term3;
const alpha1 = num1 / denom1;
const term4 = 2 * sinDelta * (omega1*omega1*l1*(m1+m2) + g*(m1+m2)*Math.cos(theta1) + omega2*omega2*l2*m2*cosDelta);
const alpha2 = term4 / denom2;
return {alpha1, alpha2};
}
We must define g, m1, m2, l1, l2 as constants.
Now initial state: let's do:
let theta1 = Math.PI/2 + (Math.random() - 0.5) * 0.5;
let theta2 = Math.PI/2 + (Math.random() - 0.5) * 0.5;
let omega1 = 0;
let omega2 = 0;
But could also give small random velocities? Not needed.
State array: let state = [theta1, theta2, omega1, omega2];
Now getPositions:
function getPositions(state) {
const [theta1, theta2] = state;
const x1 = pivotX + l1 * Math.sin(theta1) * ppm;
const y1 = pivotY + l1 * Math.cos(theta1) * ppm;
const x2 = x1 + l2 * Math.sin(theta2) * ppm;
const y2 = y1 + l2 * Math.cos(theta2) * ppm;
return {x1, y1, x2, y2};
}
Note: l1 and l2 are in meters, ppm is pixels per meter.
Now draw function:
function draw() {
ctx.fillStyle = 'rgba(0,0,0,1)';
ctx.fillRect(0,0,width,height);
const {x1, y1, x2, y2} = getPositions(state);
// Draw arms
ctx.strokeStyle = '#fff';
ctx.lineWidth = 3;
ctx.beginPath();
ctx.moveTo(pivotX, pivotY);
ctx.lineTo(x1, y1);
ctx.lineTo(x2, y2);
ctx.stroke();
// Draw masses
ctx.fillStyle = '#fff';
ctx.beginPath(); ctx.arc(x1, y1, 6, 0, Math.PI*2); ctx.fill();
ctx.beginPath(); ctx.arc(x2, y2, 6, 0, Math.PI*2); ctx.fill();
// Draw trail
if (trail.length > 1) {
ctx.lineCap = 'round';
ctx.lineJoin = 'round';
for (let i = 0; i < trail.length - 1; i++) {
const p1 = trail[i];
const p2 = trail[i+1];
const opacity = 1 - p1.age / maxTrailAge;
if (opacity <= 0.01) continue;
ctx.strokeStyle = `hsla(${p1.hue}, 100%, 50%, ${opacity})`;
ctx.lineWidth = 2;
ctx.beginPath();
ctx.moveTo(p1.x, p1.y);
ctx.lineTo(p2.x, p2.y);
ctx.stroke();
}
}
}
Now animation loop:
let lastTime = performance.now();
let accumulator = 0;
const fixedDt = 0.01; // seconds
let trail = [];
const maxTrailAge = 2.0; // seconds
let globalHue = 0;
let frameCount = 0;
function animate(time) {
const deltaTime = (time - lastTime) / 1000;
lastTime = time;
accumulator += deltaTime;
while (accumulator >= fixedDt) {
stepPhysics(fixedDt);
accumulator -= fixedDt;
}
// Add trail point (once per frame)
const {x2} = getPositions(state);
const tipX = x2;
const tipY = getPositions(state).y2; // but we already computed x2 and y2? Actually getPositions returns both, we should store.
// Better:
const pos = getPositions(state);
trail.push({x: pos.x2, y: pos.y2, hue: globalHue, age: 0});
// Age all points
for (let p of trail) {
p.age += deltaTime;
}
// Remove old
trail = trail.filter(p => p.age < maxTrailAge);
// Update global hue
globalHue = (globalHue + 0.5) % 360;
draw();
frameCount++;
requestAnimationFrame(animate);
}
Start: requestAnimationFrame(animate);
Now handle resize:
function resize() {
width = canvas.width = window.innerWidth;
height = canvas.height = window.innerHeight;
pivotX = width/2;
pivotY = height/2;
ppm = Math.min(width, height) / 4;
}
window.addEventListener('resize', resize);
resize();
Now initial state and start.
But note: stepPhysics uses state array which is initially defined. We need to define state after constants? Yes.
Now we also need to consider that if we set fixedDt=0.01 and use RK4, we may get good results. But we should test: With RK4, the integration is order 4, so error is O(dt^4). Good.
Potential performance: RK4 does 4 derivative evaluations per step, each derivative does a bunch of trig. That's fine for 2D.
Now one more improvement: The trail might look like a string of connected lines. But we could also draw the trail with variable width along the path to give a fading tail effect. But not necessary.
Now we might also add a glow effect by using shadowBlur? Not needed.
Now let's consider that the double pendulum might be too energetic and the tip might exceed screen bounds. Our scaling ensures that the maximum reach (l1+l2) times ppm is at most half the min dimension, so tip stays within the canvas if it doesn't exceed 2m from pivot. Actually if the pendulum swings, the tip's distance from pivot is sqrt(l1^2 + l2^2 + 2*l1*l2*cos(δ)? Actually the maximum distance is l1+l2 (when aligned straight outward) and minimum is |l1-l2|. So our scaling ensures max distance fits. So safe.
Now potential improvement: We could also draw the pendulum arms with a gradient or color based on angle? Not requested.
Now we should also handle the case where the trail array becomes large; we filter each frame by creating new array. That's okay for <1000.
Now code structure: all inside script tag.
Now we need to ensure that we use 'use strict'; maybe not needed.
Now let's write the final code.
But before finalizing, we should test mentally: On first frame, deltaTime maybe ~16ms. accumulator becomes 0.016, which is >= fixedDt=0.01, so we do one step, accumulator becomes 0.006. Then we add trail point based on new state. Good.
Now if frame rate is lower, say 30fps (33ms), accumulator becomes 0.033, while loop: step (accumulator>=0.01), do step, accumulator=0.023, still >=, do second step, accumulator=0.013, do third step? Actually subtract 0.01 each time: after 3 steps, accumulator = 0.003. So we do 3 steps. That's okay; physics simulation will be accurate.
Now ensure that we add trail point after all steps, representing the tip at the end of all steps for that frame. That's fine.
Now we need to consider the effect of multiple steps on the trail point: we only add one point per frame, so the trail spacing corresponds to frame time, not physics steps. That's fine.
Now we also need to consider that the state might become NaN if something goes wrong. We could add a sanity check: if (isNaN(state[0])) reset or something. Not needed for demo.
Now we also might want to add a small damping to prevent it from gaining energy uncontrollably due to integration errors? RK4 should be fine.
Now about the trail colors: we use globalHue that increments by 0.5 each frame, so it cycles in ~720 frames (12 seconds at 60fps). That's okay.
Alternatively, we could base hue on something like velocity or angle to make it more interesting. But not required.
Now let's write the HTML.
We'll write the code in a clean manner.
Potential improvements: Use requestAnimationFrame with time parameter. That's fine.
Now we must ensure all variables are declared with let/const.
I'll write:
const canvas = document.getElementById('canvas');
const ctx = canvas.getContext('2d');
let width, height, pivotX, pivotY, ppm;
// Physical constants
const g = 9.81;
const m1 = 1, m2 = 1;
const l1 = 1, l2 = 1; // meters
// Simulation state
let state; // [θ1, θ2, ω1, ω2]
let trail = [];
const maxTrailAge = 2.0; // seconds
let lastTime = performance.now();
let accumulator = 0;
const fixedDt = 0.01;
let globalHue = 0;
let frameCount = 0;
function resize() { ... }
function computeAlpha(theta1, theta2, omega1, omega2) { ... }
function derivs(state) { ... }
function stepPhysics(dt) { ... }
function getPositions(state) { ... }
function draw() { ... }
function animate(time) { ... }
resize();
// Initialize state with random angles
const rand = () => Math.random() * 0.4 - 0.2; // small offset
state = [Math.PI/2 + rand(), Math.PI/2 + rand(), 0, 0];
requestAnimationFrame(animate);
But careful: The random offset should be small enough but can be anything.
Alternatively, we can pick random in [-π, π]? That might be too extreme and cause very fast spinning, but it's fine.
I'll do: state = [Math.PI/2, Math.PI/2, 0, 0]; Actually that's symmetric and might be regular. So add a small offset: e.g., state = [Math.PI/2, Math.PI/2 + 0.1, 0, 0]; That should break symmetry. But we can randomize: state = [Math.PI/2 + (Math.random()-0.5)*1, Math.PI/2 + (Math.random()-0.5)*1, 0, 0]; That gives up to 0.5 rad offset each. That's fine.
Now we need to be careful about the while loop: use while (accumulator >= fixedDt). That's fine.
Now also we need to ensure that accumulator doesn't become negative due to floating errors. We can do accumulator %= fixedDt? But not needed.
Now we also need to consider that if accumulator is very large (like if tab inactive for a long time), the while loop could take many steps and freeze. To prevent spiral of death, we can cap the number of steps per frame, e.g., maxSteps = 10. Then while (accumulator >= fixedDt && steps < maxSteps) { ... }. But we'll ignore for simplicity.
Now we should also note that we use l1 and l2 as constants. In computeAlpha we used l1, l2, m1, m2, g. Ensure they are defined in scope.
Now we need to check the term2: Math.sin(theta1 - 2*theta2). That's correct.
But wait: Some sources have sin(θ1 - 2θ2) with a positive sign in the numerator for the second term? Actually the term is -m2*g*sin(θ1 - 2θ2). So we have negative.
But I've also seen: -m2*g*sin(θ1 - 2θ2) appears in numerator of ω1̇. Yes.
But is it sin(θ1 - 2θ2) or sin(2θ2 - θ1)? sin(θ1 - 2θ2) = -sin(2θ2 - θ1). So sign matters. Our term2 is negative: -m2*g*sin(theta1 - 2*theta2). This equals +m2*g*sin(2*theta2 - theta1) if we flip. But we'll trust as is.
Let's cross-check with a known implementation. I'll quickly recall from memory: Many sources provide:
dω1 = (-g*(2*m1+m2)*sin(θ1) - m2*g*sin(θ1-2θ2) - 2*sin(θ2-θ1)*m2*(ω2^2*l2 + ω1^2*l1*cos(θ2-θ1))) / (l1*(2*m1+m2 - m2*cos(2*(θ2-θ1)))
Yes.
So we are good.
Now one more thing: In term3 we have -2 * sinDelta * m2 * (omega2^2 * l2 + omega1^2 * l1 * cosDelta). That matches.
Now term4: 2*sinDelta*(omega1^2*l1*(m1+m2) + g*(m1+m2)*cosθ1 + omega2^2*l2*m2*cosDelta). Yes.
Now we should also consider that the equations use lengths l1 and l2, and we're using l1=1, l2=1. That's fine.
Now we need to think about potential divide-by-zero if denomFactor is zero. That would require cos2delta = (2*m1+m2)/m2. Since RHS > 1 for positive masses, impossible because cos2delta ≤1. So safe.
Now about the trail drawing: We need to ensure we use the correct opacity: p1.age / maxTrailAge gives fraction; we subtract from 1. But we must ensure opacity is between 0 and 1. So if age > maxTrailAge, we already filtered out.
Now we might want to also vary line width along the trail? Not required.
Now we also want the trail to be "colorful". Using HSLA with varying hue works.
Now also note: The trail could be drawn with a small blur for glow? Not needed.
Now we need to consider performance: The trail drawing loop does:
ctx.beginPath(); ctx.moveTo...; ctx.lineTo...; ctx.stroke(); for each segment. That's many draw calls. Could optimize by grouping segments with same hue into a single path? But hue changes per segment. Could use a gradient path? Not easily.
Maybe we can draw the trail as a set of circles with radius based on age, to give a fading dotted effect. That would also be many draw calls.
But for 200 segments, it's fine.
Now we should also consider that the trail might be too long? We set maxTrailAge=2 seconds, so at 60fps that's up to 120 points. That's fine.
Now we need to test that the code is syntactically correct.
Now let's write the file.
We'll include everything.
One more thing: In the derivative function, we use l1 and l2 as distances in meters, and g. That's fine.
But we also have to be careful about mutating state in stepPhysics: we use state.map to create arrays for k2 and k3, but state.map returns a new array. That's fine.
But we need to ensure that in stepPhysics, we don't modify state until the final update. That's correct.
Now we also need to ensure that we are using the correct variable names: l1, l2, m1, m2, g.
Now also note: In computeAlpha, we use the same m1, m2, l1, l2, g. They are constants defined outside. So they are in closure.
Now for the colors: we use `hsla(${p1.hue}, 100%, 50%, ${opacity})`. That's fine.
Now we also should set ctx.lineWidth for trail segments maybe 2 or 3.
Now maybe we want the trail to be slightly transparent overall so it fades nicely.
Now let's think about the initial point: we add trail point each frame after physics. So the first point added will be after the first physics step. That's fine. But the trail may appear to start at the first moved position, not the exact starting position. To include the starting point, we could add a point before the loop. We could push the initial tip position to trail before requestAnimationFrame. Let's do that: after initializing state, compute initial pos and push to trail with age 0. Then in animate, we add each frame.
Thus:
state = [θ1,...];
const initPos = getPositions(state);
trail.push({x: initPos.x2, y: initPos.y2, hue: globalHue, age: 0});
But globalHue is 0 initially. That's okay.
Then we start animation.
But then inside animate we add a new point each frame. That's fine.
Now we need to set globalHue and increment per frame.
Now we also need to ensure that the initial point ages correctly when we start aging. At first frame, we age all points by deltaTime (which is time from initialization to first frame? Actually we start requestAnimationFrame after we push initial point. The first animate call will happen at some time after that, so deltaTime = time - lastTime. We set lastTime = performance.now() after we set state? Actually we should set lastTime before calling animate. Let's structure:
- Set up everything.
- Initialize state.
- Compute initial trail point.
- Set lastTime = performance.now().
- Then requestAnimationFrame(animate);
We'll do that.
Now code:
window.addEventListener('load', () => {
// all code here
} but not necessary; we can place script at end of body, it will run after DOM. So we can just put script with code directly, and it'll run when parsed. That's fine.
Now final code:
<!DOCTYPE html>
<html>
<head>
<meta charset="UTF-8">
<title>Double Pendulum</title>
<style>
* { margin: 0; padding: 0; box-sizing: border-box; }
body { overflow: hidden; background: #000; }
canvas { display: block; }
</style>
</head>
<body>
<canvas id="canvas"></canvas>
<script>
// All JS code
</script>
</body>
</html>
Now inside script:
// Get context and dimensions
const canvas = document.getElementById('canvas');
const ctx = canvas.getContext('2d');
let width, height, pivotX, pivotY, ppm;
// Physical constants
const g = 9.81;
const m1 = 1, m2 = 1;
const l1 = 1, l2 = 1; // meters
// Simulation variables
let state; // [θ1, θ2, ω1, ω2]
let trail = [];
const maxTrailAge = 2.0; // seconds
let lastTime = 0;
let accumulator = 0;
const fixedDt = 0.01; // seconds per physics step
let globalHue = 0;
// Resize handling
function resize() {
width = canvas.width = window.innerWidth;
height = canvas.height = window.innerHeight;
pivotX = width / 2;
pivotY = height / 2;
ppm = Math.min(width, height) / 4; // pixels per meter
}
window.addEventListener('resize', resize);
resize();
// Initialize state with slightly perturbed equilibrium
function randomOffset() {
return (Math.random() - 0.5) * 0.5; // ±0.25 rad
}
theta1_initial = Math.PI/2 + randomOffset();
theta2_initial = Math.PI/2 + randomOffset();
state = [theta1_initial, theta2_initial, 0, 0];
// Compute accelerations from current angles and angular velocities
function computeAlpha(theta1, theta2, omega1, omega2) {
const delta = theta2 - theta1;
const sinDelta = Math.sin(delta);
const cosDelta = Math.cos(delta);
const cos2delta = Math.cos(2*delta);
const denomFactor = 2*m1 + m2 - m2*cos2delta;
const denom1 = l1 * denomFactor;
const denom2 = l2 * denomFactor;
const term1 = -g * (2*m1 + m2) * Math.sin(theta1);
const term2 = -m2 * g * Math.sin(theta1 - 2*theta2);
const term3 = -2 * sinDelta * m2 * (omega2*omega2*l2 + omega1*omega1*l1*cosDelta);
const num1 = term1 + term2 + term3;
const alpha1 = num1 / denom1;
const term4 = 2 * sinDelta * (
omega1*omega1*l1*(m1 + m2)
+ g*(m1 + m2)*Math.cos(theta1)
+ omega2*omega2*l2*m2*cosDelta
);
const alpha2 = term4 / denom2;
return { alpha1, alpha2 };
}
// Compute derivatives for RK4
function derivs(s) {
const [theta1, theta2, omega1, omega2] = s;
const {alpha1, alpha2} = computeAlpha(theta1, theta2, omega1, omega2);
return [omega1, omega2, alpha1, alpha2];
}
// Perform one RK4 step
function stepPhysics(dt) {
const k1 = derivs(state);
const k2 = derivs(state.map((v,i) => v + 0.5*dt*k1[i]));
const k3 = derivs(state.map((v,i) => v + 0.5*dt*k2[i]));
const k4 = derivs(state.map((v,i) => v + dt*k3[i]));
for (let i = 0; i < 4; i++) {
state[i] += (dt/6) * (k1[i] + 2*k2[i] + 2*k3[i] + k4[i]);
}
}
// Convert state to pixel coordinates
function getPositions(s) {
const [theta1, theta2] = s;
const x1 = pivotX + l1 * Math.sin(theta1) * ppm;
const y1 = pivotY + l1 * Math.cos(theta1) * ppm;
const x2 = x1 + l2 * Math.sin(theta2) * ppm;
const y2 = y1 + l2 * Math.cos(theta2) * ppm;
return {x1, y1, x2, y2};
}
// Draw everything
function draw() {
// Clear
ctx.fillStyle = 'rgba(0,0,0,1)';
ctx.fillRect(0,0,width,height);
const {x1, y1, x2, y2} = getPositions(state);
// Arms
ctx.strokeStyle = '#fff';
ctx.lineWidth = 3;
ctx.beginPath();
ctx.moveTo(pivotX, pivotY);
ctx.lineTo(x1, y1);
ctx.lineTo(x2, y2);
ctx.stroke();
// Masses
ctx.fillStyle = '#fff';
ctx.beginPath();
ctx.arc(x1, y1, 6, 0, Math.PI*2);
ctx.fill();
ctx.beginPath();
ctx.arc(x2, y2, 6, 0, Math.PI*2);
ctx.fill();
// Trail
if (trail.length > 1) {
ctx.lineCap = 'round';
ctx.lineJoin = 'round';
for (let i = 0; i < trail.length - 1; i++) {
const p1 = trail[i];
const p2 = trail[i+1];
const opacity = 1 - p1.age / maxTrailAge;
if (opacity <= 0.01) continue;
ctx.strokeStyle = `hsla(${p1.hue}, 100%, 50%, ${opacity})`;
ctx.lineWidth = 2;
ctx.beginPath();
ctx.moveTo(p1.x, p1.y);
ctx.lineTo(p2.x, p2.y);
ctx.stroke();
}
}
}
// Animation loop
function animate(time) {
if (!lastTime) lastTime = time;
const deltaTime = (time - lastTime) / 1000; // seconds
lastTime = time;
accumulator += deltaTime;
// Fixed timestep simulation
while (accumulator >= fixedDt) {
stepPhysics(fixedDt);
accumulator -= fixedDt;
}
// Add a trail point for the current tip
const pos = getPositions(state);
trail.push({x: pos.x2, y: pos.y2, hue: globalHue, age: 0});
// Age all points
for (let p of trail) {
p.age += deltaTime;
}
// Remove old points
trail = trail.filter(p => p.age < maxTrailAge);
// Cycle hue
globalHue = (globalHue + 0.5) % 360;
draw();
requestAnimationFrame(animate);
}
// Add initial trail point
const initialPos = getPositions(state);
trail.push({x: initialPos.x2, y: initialPos.y2, hue: globalHue, age: 0});
// Start animation
requestAnimationFrame(animate);
</script>
</body>
</html>
That should be the complete file.
But we need to check that all variables are defined. We used randomOffset function but defined as theta1_initial etc. That's fine.
But there is a subtlety: In stepPhysics, we used state.map((v,i) => ...). That's okay because state is an array.
Also, we used l1 and l2 in computeAlpha, which are meters. That's fine.
Now, is the RK4 implementation correct? Yes, standard.
But note: In RK4, we compute k1 = f(state). Then k2 = f(state + 0.5*dt*k1). That's correct. Then k3, k4. Then new state = state + dt/6*(k1 + 2*k2 + 2*k3 + k4). That's correct.
But we must ensure that we use the same dt for k2, k3, k4. Yes.
Now we should consider that state might become very large angles; but that's okay.
Now about the trail: we add a point each frame, but the physics steps may have advanced the state many steps. The point added corresponds to the state after all those steps. That's fine.
But the age of points is incremented by deltaTime, which is real time between frames. That's fine.
Now potential issue: If we set maxTrailAge = 2.0, but deltaTime could be large on first frame (e.g., if the page loads and the first call to animate happens after a long delay due to browser throttling? Usually the first frame is after load, it might be ~0.1s or so. That's okay; we add the initial point with age 0. Then we age by deltaTime; if deltaTime > maxTrailAge on the first frame, the initial point would be removed before being drawn. But deltaTime on first frame is the time between page load and first animation frame, which is typically around 16ms or maybe a bit more due to setup. Not a problem.
Now test mental: The code seems fine.
Now I'm thinking about the trail: We draw line segments between p1 and p2. But p1 and p2 are consecutive points. The opacity of the segment is based on p1's age. So the segment from p1 to p2 will have constant opacity along its length. That's okay. However, the segment from p2 to p3 will have opacity based on p2's age. So there is a step change in opacity at each point. But since points are close, it's okay.
Alternative: We could draw the trail as a polyline with a single stroke but with a gradient along the path using a canvas gradient. But that would be based on coordinates, not age. Not needed.
Now also we could draw the trail with a slight blur for aesthetics? Not required.
Now we might want to add a "bloom" effect by using a larger shadowBlur and shadowColor. That could look nice. But we'll keep simple.
Now we should ensure that the code is within a single HTML file and no external dependencies.
Now I'll produce the final answer: the content of index.html.
But before final, let's check for any syntax errors:
- In stepPhysics, we have `state.map((v,i) => v + 0.5*dt*k1[i])`. That's valid.
- In draw, we used `ctx.fillStyle = 'rgba(0,0,0,1)';` That's okay.
- In computeAlpha, we used `Math.cos(2*delta)`. That's fine.
- In animation, we used `if (!lastTime) lastTime = time;` but we already set lastTime before calling requestAnimationFrame? Actually we didn't set lastTime initially. We can either set let lastTime = 0; and then in animate, if lastTime===0, set it to time. But we also add initial point before starting animation. We could set lastTime = performance.now() just before requesting first frame. That might be easier: after we add initial trail point, set lastTime = performance.now(); then requestAnimationFrame(animate). Then in animate, we compute deltaTime as (time - lastTime)/1000, and then set lastTime = time. That's straightforward.
So I'll modify:
// Initial setup
let lastTime = performance.now();
requestAnimationFrame(animate);
And inside animate, we start with:
const deltaTime = (time - lastTime) / 1000;
lastTime = time;
We should also set lastTime before the first call. That's done.
So we'll define let lastTime = performance.now(); right before requestAnimationFrame.
Now also we need to define accumulator = 0 initially.
Now one more thing: In the while loop, we might have accumulator be slightly negative due to floating error after subtraction. So we could do accumulator = Math.max(0, accumulator); but not critical.
Now finalize.
But before final answer, I'll test the code in my head for a few steps:
Initial state: θ1=π/2=1.5708, θ2=π/2+0.1=1.6708, ω1=ω2=0.
Compute alpha:
delta = 0.1, sinDelta=0.0998, cosDelta=0.995, cos2delta = cos(0.2)=0.98007.
denomFactor = 2*1+1 -1*0.98007 = 3 - 0.98007 = 2.01993.
denom1 = 1*2.01993 = 2.01993.
denom2 = 1*2.01993 = 2.01993.
term1 = -9.81*3*sin(1.5708)= -9.81*3*1 = -29.43.
term2 = -1*9.81*sin(1.5708 - 2*1.6708)= -9.81*sin(1.5708 - 3.3416)= -9.81*sin(-1.7708)= -9.81*(-0.978)= +9.81*0.978? Actually sin(-1.7708) = -sin(1.7708) = -0.978. So term2 = -m2*g*(-0.978) = +9.81*0.978 = 9.59. Wait carefully: term2 = -m2*g*Math.sin(theta1 - 2*theta2). theta1 - 2θ2 = 1.5708 - 3.3416 = -1.7708. sin(-1.7708) ≈ -0.978. So term2 = -1*9.81*(-0.978) = +9.81*0.978 = 9.59.
term3 = -2*sinDelta*m2*(ω2^2*l2 + ω1^2*l1*cosDelta) = 0 because ω1=ω2=0.
num1 = -29.43 + 9.59 = -19.84. alpha1 = -19.84/2.01993 ≈ -9.82 rad/s^2.
So θ1 will start accelerating downward (negative since positive is say rightward? Actually θ1 positive is to right? If pivot at center, θ=0 is vertical down. Positive θ means pendulum swings to right? sin(θ) positive => x positive, so positive θ swings right. Our initial θ1=90° (π/2) means pendulum to right? Actually at π/2, sin=1, so x positive, y = cos(π/2)=0, so horizontal to right? That would be horizontal? Wait, if θ measured from vertical down, then at θ=0, pendulum hanging down: x=0, y=l. At θ=π/2, sin=1, cos=0, so x=l (right), y=0 (pivot height). So that's horizontal to the right. That's an unstable equilibrium? Actually a pendulum at horizontal is not equilibrium; it will swing down. The acceleration should be negative (toward zero). Alpha1 negative means ω1 decreases (since ω1=0, becomes negative) meaning θ1 will decrease from π/2 toward 0. That makes sense: it will swing downward to the right? Actually negative means clockwise? Let's think: If θ1 is positive (to right), then negative angular acceleration means it will decrease θ, i.e., rotate back towards zero (downward). So that's correct.
Now alpha2: term4 = 2*sinDelta * (ω1^2*l1*(m1+m2) + g*(m1+m2)*cosθ1 + ω2^2*l2*m2*cosDelta). The first and third terms zero because ω=0. So = 2*0.0998 * (g*2*cosθ1). cosθ1 = cos(π/2)=0. So term4 = 0. So alpha2 = 0. So second pendulum initially has no angular acceleration? That seems plausible because both masses are at rest and symmetric? Actually if θ2=θ1+0.1, and both masses at rest, there is a torque? But cosθ1=0 so gravity term zero? Wait: cosθ1 is cos(π/2)=0, so the term g*(m1+m2)*cosθ1 is zero. So indeed alpha2=0 initially. So second pendulum initially has zero angular acceleration, but after θ1 changes, it will get acceleration. So that seems okay.
Thus initial step seems fine.
Now the code should work.
Now about the trail: we add initial point before loop. Then each frame we add new point after physics steps. That's fine.
Now we should also consider that we may want to clear the trail on reset if we reset state? Not needed.
Now final code: we'll produce the HTML.
But wait: The user said "draw the arms and a colorful fading trail of the tip." That's what we did.
Now we need to ensure that the simulation is "chaotic". The double pendulum is chaotic for many initial conditions. Our random perturbation should give chaotic behavior.
Now I'll write the final file content.
Let's ensure proper indentation and no missing semicolons.
One more thing: In computeAlpha, we used const delta = theta2 - theta1; But what if delta is large? That's fine.
Now also, the equations assume l1 and l2 are in meters. We treat them as 1. That's fine.
Now about the RK4: We compute state.map for k2 and k3. That creates new arrays each time. That's okay.
But we need to be careful: state.map((v,i) => v + 0.5*dt*k1[i]) uses k1[i]. That is the derivative for component i. Since state has 4 elements, k1 has 4.
Yes.
Now we also need to define the variable `frameCount`? We didn't actually use it. We can remove it. I used globalHue instead. So I'll remove frameCount.
Now also we don't use randomOffset function after we set initial state. That's fine.
Now I'll produce the final code.
Potential improvement: Use `performance.now()` inside requestAnimationFrame directly, but that's fine.
Now we need to include the <script> tag with type="javascript"? Not required.
Now final file:
I'll output exactly the HTML content.
Let's write it.```html <!DOCTYPE html> <html> <head> <meta charset="UTF-8"> <title>Double Pendulum</title> <style> * { margin: 0; padding: 0; box-sizing: border-box; } body { overflow: hidden; background: #000; } canvas { display: block; } </style> </head> <body> <canvas id="canvas"></canvas> <script> // Get canvas and context const canvas = document.getElementById('canvas'); const ctx = canvas.getContext('2d'); let width, height, pivotX, pivotY, ppm; // Physical constants (SI units) const g = 9.81; const m1 = 1, m2 = 1; const l1 = 1, l2 = 1; // lengths in meters // Simulation state: [θ1, θ2, ω1, ω2] let state; let trail = []; const maxTrailAge = 2.0; // seconds let lastTime = 0; let accumulator = 0; const fixedDt = 0.01; // physics timestep in seconds let globalHue = 0; // Resize canvas and scale function resize() { width = canvas.width = window.innerWidth; height = canvas.height = window.innerHeight; pivotX = width / 2; pivotY = height / 2; ppm = Math.min(width, height) / 4; // pixels per meter } window.addEventListener('resize', resize); resize(); // Random offset for initial angles (to break symmetry) function randomOffset() { return (Math.random() - 0.5) * 0.5; } // Initial state: slightly perturbed from vertical (π/2) state = [ Math.PI/2 + randomOffset(), Math.PI/2 + randomOffset(), 0, 0 ]; // Compute angular accelerations (α1, α2) from current state function computeAlpha(theta1, theta2, omega1, omega2) { const delta = theta2 - theta1; const sinDelta = Math.sin(delta); const cosDelta = Math.cos(delta); const cos2delta = Math.cos(2*delta); const denomFactor = 2*m1 + m2 - m2*cos2delta; const denom1 = l1 * denomFactor; const denom2 = l2 * denomFactor; // Numerator for α1 const term1 = -g * (2*m1 + m2) * Math.sin(theta1); const term2 = -m2 * g * Math.sin(theta1 - 2*theta2); const term3 = -2 * sinDelta * m2 * (omega2*omega2*l2 + omega1*omega1*l1*cosDelta); const num1 = term1 + term2 + term3; const alpha1 = num1 / denom1; // Numerator for α2 const term4 = 2 * sinDelta * ( omega1*omega1*l1*(m1 + m2) + g*(m1 + m2)*Math.cos(theta1) + omega2*omega2*l2*m2*cosDelta ); const alpha2 = term4 / denom2; return { alpha1, alpha2 }; } // Derivatives for RK4: [dθ1/dt, dθ2/dt, dω1/dt, dω2/dt] function derivs(s) { const [theta1, theta2, omega1, omega2] = s; const {alpha1, alpha2} = computeAlpha(theta1, theta2, omega1, omega2); return [omega1, omega2, alpha1, alpha2]; } // Perform one Runge-Kutta 4 integration step function stepPhysics(dt) { const k1 = derivs(state); const k2 = derivs(state.map((v,i) => v + 0.5*dt*k1[i])); const k3 = derivs(state.map((v,i) => v + 0.5*dt*k2[i])); const k4 = derivs(state.map((v,i) => v + dt*k3[i])); for (let i = 0; i < 4; i++) { state[i] += (dt/6) * (k1[i] + 2*k2[i] + 2*k3[i] + k4[i]); } } // Convert simulation state to pixel coordinates function getPositions(s) { const [theta1, theta2] = s; const x1 = pivotX + l1 * Math.sin(theta1) * ppm; const y1 = pivotY + l1 * Math.cos(theta1) * ppm; const x2 = x1 + l2 * Math.sin(theta2) * ppm; const y2 = y1 + l2 * Math.cos(theta2) * ppm; return { x1, y1, x2, y2 }; } // Draw everything function draw() { // Clear canvas ctx.fillStyle = 'rgba(0,0,0,1)'; ctx.fillRect(0, 0, width, height); const { x1, y1, x2, y2 } = getPositions(state); // Draw arms ctx.strokeStyle = '#fff'; ctx.lineWidth = 3; ctx.beginPath(); ctx.moveTo(pivotX, pivotY); ctx.lineTo(x1, y1); ctx.lineTo(x2, y2); ctx.stroke(); // Draw masses ctx.fillStyle = '#fff'; ctx.beginPath(); ctx.arc(x1, y1, 6, 0, Math.PI*2); ctx.fill(); ctx.beginPath(); ctx.arc(x2, y2, 6, 0, Math.PI*2); ctx.fill(); // Draw colorful fading trail if (trail.length > 1) { ctx.lineCap = 'round'; ctx.lineJoin = 'round'; for (let i = 0; i < trail.length - 1; i++) { const p1 = trail[i]; const p2 = trail[i+1]; const opacity = 1 - p1.age / maxTrailAge; if (opacity <= 0.01) continue; ctx.strokeStyle = `hsla(${p1.hue}, 100%, 50%, ${opacity})`; ctx.lineWidth = 2; ctx.beginPath(); ctx.moveTo(p1.x, p1.y); ctx.lineTo(p2.x, p2.y); ctx.stroke(); } } } // Animation loop function animate(time) { const deltaTime = (time - lastTime) / 1000; // seconds lastTime = time; accumulator += deltaTime; // Fixed timestep simulation while (accumulator >= fixedDt) { stepPhysics(fixedDt); accumulator -= fixedDt; } // Add current tip to trail const pos = getPositions(state); trail.push({ x: pos.x2, y: pos.y2, hue: globalHue, age: 0 }); // Age trail points and remove old ones for (let p of trail) p.age += deltaTime; trail = trail.filter(p => p.age < maxTrailAge); // Cycle hue for next segment globalHue = (globalHue + 0.5) % 360; draw(); requestAnimationFrame(animate); } // Add initial trail point (starting tip position) const initPos = getPositions(state); trail.push({ x: initPos.x2, y: initPos.y2, hue: globalHue, age: 0 }); // Start animation lastTime = performance.now(); requestAnimationFrame(animate); </script> </body> </html> ```
TAKE IT TO CODE
This run produced no artifact — reuse the prompt with a different model or your own pi-agent run.
RUN
modelstep-3.5-flash
vendorStepFun
harnesspi-agent
generated2026-08-02
METRICS
runtime353.4s
tokens in1.6k
tokens out24k
cached0
cost / run$0.0073