/* Variablendeklaration */ var ODE_DP; // ODE-Solver 2. Ordnung für 2 verkoppelte DGl ... var t = 0; // Zeitvariable var dt; // Increment der Zeitv. var frmRate; // Screen-Refreshrate var M; // dyn. Maßstab var y, vy; // Bewegungsvariable var l0 = 1; // Ruhelängen der Feder in m var d; // Däpfungskonstante Feder var m = 0.1; // Masse des Objektes = 0.1 kg var n0; // Federkonstante var kappa; // Nichtlinearität var Exponent = 5; // Potenzreihe, Exponent var yExcite; // Anregungsamplitude var f0 = 3, omega0; // Eigenfrequenz in Hz var yExcAmpl; // max. Anregungsamplitude function setup() { ... xi0 = 0.2*width; // int. Nullpunkt für kart. Koordinatensystem rel. zum internen K. yi0 = 0.9*height; M = height/(10*l0); // Berechnung des dyn. Maßstabs frmRate = 60; frameRate(frmRate); dt = 1.0/frmRate; ... } function draw() { ... //****************************** Berechnung der Bewegung ****************************** if (START) { ... } else { // Berechnung der Zeitfunktionen if (j++ == frmRate) j = 0; // Anregungsfunktion if (j > 20 & j < 30) yExcite = yExcAmpl; else yExcite = 0; result = ODE_DP.DPsO(grad, y, vy, dt); y = result[0]; vy = result[1]; t = t + dt; // Zeitachse für Simulation } //****************************************** Administration ******************************************** // Hier wird in Pixeln gerechnet! prepareScreen("Schwingungen in Feder/Masse-Systemen", "Feder mit nichtlinearer Federkonstante", 0, 50, 5); // Titel des Programms } function grad(s, vs) { return -(2*d*vs + springForce(s, yExcite, l0, kappa, Exponent)); } function springConst(s, l0, kappa, Imax) // Federkonstante { var n = 1; // relative Federkonstante n/n0 var delta_s = (s - l0); for (var i = 2; i < Imax; i++) // Anzahl der Reihenglieder { n = n + pow(kappa*abs(delta_s)/l0, 2*i-1)/(kappa*(2*i-1)); } return(n); } function springForce(s, sExc, l0, kappa, Imax) // Federkraft { var s_ = 0; // Init Summe var delta_s = (s - l0 - sExc); for (var i = 2; i < Imax; i++) { s_ = s_ + pow(kappa*(delta_s)/l0, 2*i-1)/(2*i*(2*i-1)); // rechte Hälfte der Kennlinie } return(sq(omega0)*(delta_s)*(1 + abs(s_/kappa))); // komplette Kennlinie }