dim y(100) as double 'string height
dim yp(100) as double 'velocity
dim ypp(100) as double 'acceleration
dim x as integer 'index for y's
dim t as double 'time step
dim ymax as double
dim d as double 'delta and spring constant
dim k as double
d = 1
k = 1
ymax = 10
y(0) = 0 'one of the boundary conditions
'this sets up the initial condition
'the string is held like an asymmetric triangle
'and let go
'It is done asymmetrically to allow even harmonics
'to exist
'"held and let go" implies initial velocity
'is zero
for x = 1 to 99
'y(x) = ymax*((x/100)^3-(x/100))
'y(x) = ymax*(sin(3.1415926535897932384626*x/100)+sin(2*3.1415926535897932384626*x/100))
y(x) = ymax*(.5*sin(3.1415926535897932384626*x/100)+2*sin(20*3.1415926535897932384626*x/100))
next
for x = 1 to 99
yp(x) = 0
next
'we defined everything as "1", so the fundamental
'frequency is 1/2 (period = 2).
'On the other hand, the shortest possible wavelength is
'2/100, so the highest possible frequency is 50 (period = 1/50).
'Consequently, we will set time step = 1/100 and let this run for
'20.48 periods, which, as a power of 2, Excel can do a FFT on
open "stringreal-10d-4.csv" for output as #1
for t = 1 to 4096
for x = 0 to 100
print #1, y(x);",";
next
print #1,
for x = 1 to 99 'only going to 99 rather than 100 because y(100) = 0 is fixed
'top line: real
'ypp(x) = -100*100*k*((y(x)-y(x-1))*(d/100+sqr(1/100^2+(y(x)-y(x-1))^2)-1/100)/sqr(1/100^2+(y(x)-y(x-1))^2) +(y(x)-y(x+1))*(d/100+sqr(1/100^2+(y(x)-y(x+1))^2)-1/100)/sqr(1/100^2+(y(x)-y(x+1))^2))
'bottom line: small angle approximation
ypp(x) = -100*100*k*((y(x)-y(x-1))*(d/100)/sqr(1/100^2) +(y(x)-y(x+1))*(d/100)/sqr(1/100^2))
next
for x = 1 to 99
yp(x) = yp(x) + ypp(x)/100
next
for x = 1 to 99
y(x) = y(x) + yp(x)/100
next
next
close #1