Die Menschen drängen sich zum Lichte, nicht um besser zu sehen,
sondern um besser zu glänzen.
Friedrich Nietzsche
Ich habe ein kleines Simulationsprogramm geschrieben, womit man experimentieren kann. Das vollständige Programm enthält ein paar Möglichkeiten mehr als die Kurzvariante hier. Es ist im Programmarchiv abgelegt.
Hier die Gegenüberstellung der Programmzeilen mit den Formeln des vorigen Teils:
Die Details bitte den Kommentaren entnehmen.
def sq(x):
return x * x
# Main
##############################
N = 100 # Anzahl der Saitenstücke
L = 1 # Saitenlänge
h = 0.005 # Auslenkung beim Zupfen
F0 = 1 # Spannkraft
mu = 1 # Massendichte, ein knappes Gramm pro Meter
timesteps = 75
# Anzahl der Zeitschritte bis zum Programm-Ende
zp = int(N / 2) # der Zupfpunkt N/2 für mittig, N/4 für 2 Stück
initKoor(0,L,0.1,-h ,h,h / 10)
x = [0 for i in range(N + 1)] # Koordinate entlang der Saite
Y = [0 for i in range(N + 1)] # die Querauslenkung der Saite (vorige Runde)
y = [0 for i in range(N + 1)] # die Querauslenkung der Saite, berechnte Runde
v = [0 for i in range(N + 1)] # Geschwindigkeit
dx = L / N
c2 = F0 / mu # c²
dt = 0.8 * sqrt(dx * dx / c2) # Faktor 1,0 ist die Stabilitätsgrenze
print("dt=%5.3f c=%5.3f" % (dt, sqrt(c2)))
for i in range(N + 1): # Koordinate der Saite setzen
x[i] = i * L / N
Y[zp - 3] = h * 0.1
Y[zp - 2] = h * 0.4
Y[zp - 1] = h * 0.8
Y[zp + 0] = h * 0.9
Y[zp + 1] = h * 0.8
Y[zp + 2] = h * 0.4
Y[zp + 3] = h * 0.1
Schleife über alle Zeitschritte, und dort für jedes Stückchen Seil:
for j in range(timesteps): # alle Zeiten
for i in range(1,N): # alle Masse-Elemente
Die Einzelkräfte der Seilstücke:
F1,y = F0 · (y1 - y0 ) / Δx ; F2,y = F0 · (y2 - y1 ) / Δx F1 = F0 * (Y[i] - Y[i - 1]) / dx
F2 = F0 * (Y[i + 1] - Y[i]) / dx
Fq = (F2 - F1)
Die Beschleunigung benötigt die Masse des Seilstücks FQ = µ · Δx · a . Aus dieser Beschleunigung lassen sich die Orts- und die Geschwindigkeitsänderung bestimmen:
m = mu * dx # die infinitesimale Masse
a = Fq / m
v[i] = v[i] + a * dt
y[i] = Y[i] + v[i] * dt
Die so berechneten Koordinaten werden grafisch ausgegeben. Vorher müssen die Punkte der vorigen Runde gelöscht werden (unplot). Zum Schluss wird auf der Konsole die simulierte Zeit ausgegeben und der Grafikbildschirm aufgefrischt.
for i in range(N): # berechnete Ergebnisse setzen
unplot(x[i],Y[i])
plot(x[i],y[i])
Y[i] = y[i]
print(j * dt)
update() # Bildschirm anzeigen
Das war's.
--> 1 Einführung
--> 2 Funktionsweise des Programms
--> 3 Rechenergebnisse
--> 4 Sinuswelle
--> 5 Das Soliton
--> 6 Die Kugelwelle
--> 7 Interferenz zweier Wellenzüge (Formel)
Uwe Pilz, November 2020