I'm trying to make a simplified simulation of a 1-dimensional wave with a chain of harmonic oscillators. The differential equation for the position x[i]
of the i-th oscillator (assuming equilibrium at x[i]=0
) turns out to be - according to textbooks - this one:
m*x[i]''=-k(2*x[i]-x[i-1]-x[i+1])
(the derivative is wrt the time) so i tried to numerically compute the dynamics with the following algorithm. Here osc[i]
is the i-th oscillator as an object with attributes loc
(absolute location), vel
(velocity), acc
(acceleration) and bloc
(equilibrium location), dt
is the time increment:
for every i:
osc[i].vel+=dt*osc[i].acc;
osc[i].loc+=dt*osc[i].vel;
osc[i].vel*=0.99;
osc[i].acc=-k*2*(osc[i].loc-osc[i].bloc);
if(i!=0){
osc[i].acc+=+k*(osc[i-1].loc-osc[i-1].bloc);
}
if(i!=N-1){
osc[i].acc+=+k*(osc[i+1].loc-osc[i+1].bloc);
}
You can see the algorithm in action here, just click to give an impulse to the 6-th oscillator. You can see that not only it doesn't generate a wave at all but it also generate a motion with growing total energy (even if I added the dampening!). What am I doing wrong?