Сравнение методов Рунге-Кутта и Липфрога — различия между версиями

Материал из Department of Theoretical and Applied Mechanics
Перейти к: навигация, поиск
Строка 1: Строка 1:
 
[[Виртуальная лаборатория]]>[[Сравнение методов Рунге-Кутта и Липфрога]] <HR>
 
[[Виртуальная лаборатория]]>[[Сравнение методов Рунге-Кутта и Липфрога]] <HR>
−
Сравнение методов численного интегрирования: Рунге-Кутта и Липфрога
+
 +
== '''Постановка задачи ''' ==
 +
Дано простейшее уравнение движения грузика на пружине: <math>\ddot{x} = -cx</math>. Необходимо интегрировать его с помощью двух методов: Липфрога и Рунге-Кутта. Сравнить результаты.
 +
<br/>
  
−
В данном примере интегрируется уравнение грузика на пружинке: <math>\ddot{x} = -cx</math> методом Рунге-Кутта и Липфрога. Строится фазовая плоскость для каждого из методов.
+
== '''Метод численного интегрирования Рунге-Кутты''' ==
 +
Рассмотри задачу Коши для системы обыкновенных дифференциальных уравнений первого порядка.
 +
<math>\dot{y} = f(x, y), y(x_0)=y_0</math>
 +
Тогда приближенное значение в последующих точках вычисляется по итерационной формуле:
 +
<math>y_n+1 = y_n + \frac{h}{6}(k_1 + 2k_2 + 2k_3 + k_4)</math>
 +
=='''Основные уравнения'''==
 +
<div style="font-size: 16px;">1)''Расчет статических координат:''</div>
  
−
                        '''Метод Рунге-Кутта'''                                                                  '''Метод Липфрога'''
+
<math>{y}_{i}=i{l}_{0}+\sum^{i}_{k=1} {\Delta l_k}</math>
−
{{#widget:Iframe |url=http://tm.spbstu.ru/htmlets/Pogodina/Integrate/Integrate.html |width=2000 |height=600 |border=0 }}
+
<br/>
−
 
+
<br/>
−
Скачать [[Медиа:Integrate1.rar|Integrate1.rar]].
+
<math>\Delta{l_k} = \frac{(n-k+1)mg}{C} </math>
 +
<br/>
 +
<br/>
 +
где  <math>{n}</math> — количество грузиков ; <math>{i}</math> — текущий грузик ; <math>{l_0}</math> — длина не растянутой пружины; <math>\Delta{l_k}</math> — статическое удлинение пружины.
  
 +
<div style="font-size: 16px;">2) ''Уравнение движения грузиков при отсутствии верхней пружины:''</div>
 +
<math>\ddot y_i = - \frac{C}{m} ({y}_{i+1} - 2{y}_{i}+{y}_{i-1}) + g </math> 
 +
<br/>
 +
<br/>
 +
при <math>i = n </math> :
 +
<br/>
 +
<math> \ddot y_n = - \frac{C}{m}(y_n - {y}_{n-1} - l_0) + g</math>
 +
<br/>
 +
<br/>'''Программа (лучше всего смотреть ее в Mozilla Firefox)'''
 +
<br/>
 +
{{#widget:Iframe |url=http://tm.spbstu.ru/htmlets/Alexandrov/Spring_illusion.html |width=800|height=800 |border=0 }}
 +
<br />
 +
<br />
 +
Скачать программу [[Медиа:Spring_illusion.rar|Spring_illusion.rar]].
 +
<br />
 +
<br />
 
<div class="mw-collapsible mw-collapsed" style="width:100%" >
 
<div class="mw-collapsible mw-collapsed" style="width:100%" >
−
'''Текст программы на языке JavaScript (разработчик [[Погодина Валерия]], код основан на программе [[Цветков Денис]]):''' <div class="mw-collapsible-content">
+
'''Текст программы на языке JavaScript (разработчик [[Александров Сергей]]):''' <div class="mw-collapsible-content">
−
Файл '''"Integral3"'''
+
Файл '''"Spring_illusin.js"'''
 
<syntaxhighlight lang="javascript" line start="1" enclose="div">
 
<syntaxhighlight lang="javascript" line start="1" enclose="div">
−
window.addEventListener("load", Main_Integrate, true);
+
window.addEventListener("load", Main_Spring, true);
−
function Main_Integrate() {
+
function Main_Spring() {
−
+
 
−
    // Предварительные установки
+
// *** Некие исходные данные ***
−
+
 
−
    var context = Integrate_canvas.getContext("2d"); // на context происходит рисование
+
var canvas = spring_canvas;
−
var context1 = Integrate_canvas1.getContext("2d"); // на context происходит рисование
+
canvas.onselectstart = function () {return false;}; // запрет выделения canvas
−
var T0 = 1;                           // масштаб времени (период колебаний исходной системы)
+
var ctx = canvas.getContext("2d"); // на ctx происходит рисование
−
var time = 0;
+
var w = canvas.width; // ширина окна в расчетных координатах
−
+
var h = canvas.height; // высота окна в расчетных координатах
 +
 
 +
var Pi = 3.1415926; // число "пи"
 +
var g = 9.8; // гравитационная постоянная
 +
 
 +
var T0 = 0.01; // масштаб времени (период колебаний исходной системы)
 +
 
 +
var m0 = 0.15; // масштаб массы маятника
 +
var C0 = 1; // масштаб жесткости пружины
 +
 
 +
var count = true; // проводить ли расчет системы
 +
var v = 0; // скорость тела
 +
var t = 0;
 +
 
 +
// параметры полученные из размеров холста
 +
var rw = canvas.width / 4;
 +
var rh = canvas.height / 100;
 +
 
 +
// параметры грузиков
 +
var x0 = 0;
 +
var y0 = 0;
 +
var vy0 = 0;
 +
var rad0 = 7;
 +
var sTime = 0;
 +
// параметры пружины
 +
var Lp = 30; //длина пружины
 +
 
 +
 
 
// *** Задание вычислительных параметров ***
 
// *** Задание вычислительных параметров ***
−
+
 
−
    //const fps = 5500;                         // frames per second - число кадров в секунду (качеcтво отображения)
+
var fps = 60; // frames per second - число кадров в секунду (качечтво отображения)
−
    //const spf = 260;                       // steps per frame  - число шагов интегрирования между кадрами (скорость расчета)
+
var spf = 5; // steps per frame  - число шагов интегрирования между кадрами
−
    //const dt = 0.1 * T0 / fps;           // шаг интегрирования  
+
var dt = 50 * T0 / fps; // шаг интегрирования (качество расчета)
−
+
var steps = 0; // количество шагов интегрирования
−
const fps = 550;                         // frames per second - число кадров в секунду (качеcтво отображения)
+
 
−
    const spf = 5;                       // steps per frame  - число шагов интегрирования между кадрами (скорость расчета)
+
 
−
    const dt  = 1 * T0;           // шаг интегрирования
+
// *** Задание физических параметров ***
−
var n = 10;
+
var m = 1 * m0; // масса грузика
−
+
var C = 1 * C0; // жесткость пружины
−
Text_n.value = n;
+
var L0 = 0; // изначальное удлинение пружины
−
Slider_n.min = 1;                        
+
var n = 5; // количество тел
−
    Slider_n.max = 100;
+
var circle; // круглые тела
−
    Slider_n.step = 1;
+
 
−
    Slider_n.value = Text_n.value;
+
// *** Установка слайдеров для переменных величин ***
−
+
slider_m.value = (m / m0).toFixed(1);
−
    // Добавление шара
+
number_m.value = (m / m0).toFixed(1);
−
var particles = [];
+
slider_spf.value = (spf).toFixed(1);
−
var particles1 = [];
+
number_spf.value = (spf).toFixed(1);
−
for (var i=0; i<n; i++) {
+
slider_n.value = parseInt(n);
−
particles[i] = [];
+
number_n.value = parseInt(n);
−
particles1[i] = [];
+
 
 +
function setM(new_m) {
 +
m = new_m * m0;
 +
}
 +
 
 +
function setSpf(new_spf) {
 +
spf = new_spf;
 +
}
 +
function setN(new_n) {
 +
n = new_n;
 +
StaticStart();
 
}
 
}
−
        for (var i = 0; i < n; i++) {
+
 
−
for (var j = 0; j < n; j++) {
+
slider_m.oninput = function () {
−
var b = [];
+
number_m.value = slider_m.value;
−
var a = [];
+
setM(slider_m.value);
−
b.x = 0.3 / n * i;
+
};
−
b.vx =j * 0.3 / n;
+
number_m.oninput = function () {
−
a.x = 0.3 / n * i;
+
slider_m.value = number_m.value;
−
a.vx =j * 0.3 / n;
+
setM(number_m.value);
−
particles[i][j] = b;
+
};
−
particles1[i][j] = a;
+
 
 +
 
 +
slider_spf.oninput = function () {
 +
number_spf.value = slider_spf.value;
 +
setSpf(slider_spf.value);
 +
};
 +
number_spf.oninput = function () {
 +
slider_spf.value = number_spf.value;
 +
setSpf(number_spf.value);
 +
};
 +
 
 +
slider_n.oninput = function () {
 +
number_n.value = slider_n.value;
 +
setN(slider_n.value);
 +
};
 +
number_n.oninput = function () {
 +
slider_n.value = number_n.value;
 +
setN(number_n.value);
 +
};
 +
 
 +
// *** Условие падения тела ***
 +
var spring_dropped = false;
 +
drop_spring.onclick = function () {
 +
spring_dropped = true;
 +
};
 +
 
 +
// *** Пересчет статических координат
 +
new_static.onclick = function () {
 +
StaticStart()
 +
};
 +
 
 +
// *** Создание массива элементов ***
 +
function StaticStart() {
 +
spring_dropped = false;
 +
circle = [];
 +
for (var i = 0; i < n; i++) {
 +
var circ = {};
 +
circ.x = x0;
 +
circ.y = y0;
 +
circ.vy = vy0;
 +
circ.rad = rad0;
 +
circ.L = L0;
 +
 
 +
var rgb = HSVtoRGB(i / n * 2, 1, 1);
 +
 
 +
circ.fill = "rgba(" + rgb.r + ", " + rgb.g + ", " + rgb.b + ", 1)";
 +
circle[i] = circ;
 +
}
 +
for (var i = 0; i < n; i++) {
 +
if (i == 0) {
 +
circle[i].x = w / 2;
 +
circle[i].y0 = Lp + (n - i) * m * g / C;
 +
circle[i].y = circle[i].y0;
 +
 
 +
} else {
 +
circle[i].x = w / 2;
 +
circle[i].y0 = circle[i - 1].y0 + Lp + (n - i) * m * g / C;
 +
circle[i].y = circle[i].y0;
 +
console.log(circle[i].y0);
 +
 
 
}
 
}
−
        }
+
}
−
+
}
−
        this.set_n = function(input)
 
−
        {
 
−
                n = Number(input);
 
−
                Text_n.value  = n;
 
−
                Slider_n.value = Text_n.value; //записываем значение начального смещения
 
−
for (var i=0; i<n; i++) {
 
−
particles[i] = [];
 
−
particles1[i] = [];
 
−
}
 
−
for (var i = 0; i < n; i++) {
 
−
for (var j = 0; j < n; j++) {
 
−
var b = [];
 
−
var a = [];
 
−
b.x = 0.3 / n * i;
 
−
b.vx = j * 0.3 / n;
 
−
a.x = 0.3 / n * i;
 
−
a.vx = j * 0.3 / n;
 
−
particles[i][j] = b;
 
−
particles1[i][j] = a;
 
−
}
 
−
} ; 
 
−
                context.clearRect(0, 0, 300, 300);     
 
−
                context1.clearRect(0, 0, 300, 300);
 
−
        }
 
 
 
 +
// *** функция разукрашивания объектов
 +
function HSVtoRGB(h, s, v) {             
 +
    var r, g, b, i, f, p, q, t;
 +
    i = Math.floor(h * 6);
 +
    f = h * 6 - i;
 +
    p = v * (1 - s);
 +
    q = v * (1 - f * s);
 +
    t = v * (1 - (1 - f) * s);
 +
    switch (i % 6) {
 +
        case 0: r = v, g = t, b = p; break;
 +
        case 1: r = q, g = v, b = p; break;
 +
        case 2: r = p, g = v, b = t; break;
 +
        case 3: r = p, g = q, b = v; break;
 +
        case 4: r = t, g = p, b = v; break;
 +
        case 5: r = v, g = p, b = q; break;
 +
    }
 +
    return {
 +
        r: Math.floor(r * 255),
 +
        g: Math.floor(g * 255),
 +
        b: Math.floor(b * 255)
 +
    };
 +
}
 +
 +
    // график
 +
    var vGraph1 = new TM_graph(                  // определить график
 +
        "#vGraph1",                              // на html-элементе #vGraph
 +
        250,                                    // сколько шагов по оси "x" отображается
 +
        0, 100, 5);                            // мин. значение оси Y, макс. значение оси Y, шаг по оси Y
  
−
setInterval(control1, 1500 / fps);
+
// *** Функция обеспечивающая "жизнь" пружин ***
 +
function control() {
 +
calculate();
 +
draw();
 +
requestAnimationFrame(control);
 +
}
 +
StaticStart();
 +
control();
 +
 
 +
 
 +
 
 +
// *** Функция расчетов координат ***
 +
function calculate() {
 +
 
 +
for (var s = 1; s <= spf; s++) {
 +
var k = n;
 +
for (var i = 1; i < n; i++) {
 +
if (spring_dropped)
 +
circle[0].L = 0;
 +
else
 +
circle[0].L = circle[0].y - Lp;
 +
circle[i].L = circle[i].y - circle[i - 1].y - Lp;
 +
}
 +
 
 +
for (var i = 0; i < n; i++) {
 +
if (i == n - 1)
 +
circle[i].vy += ((-1) * C * (circle[i].L) / m + g) * dt;
 +
else
 +
circle[i].vy += ((-1) * C * (circle[i].L - circle[i + 1].L) / m + g) * dt;
 +
 +
 +
circle[i].y += circle[i].vy * dt;
 +
 
 +
}
 +
steps++;
 +
 
 +
if (spring_dropped)
 +
{ sTime ++;
 +
if (sTime % 50 == 0) vGraph1.graphIter(sTime/50, (circle[n-1].y-circle[n-1].y0)) ;
 +
};
 +
}
 +
}
 +
 
 +
 
 +
// *** Функция рисования объектов ***
 +
function draw() {
 +
ctx.clearRect(0, 0, w, h);
 +
 
 +
// Рисование закрепления
 +
ctx.lineWidth = 6;
 +
ctx.strokeStyle = "#7394cb";
 +
ctx.beginPath();
 +
ctx.moveTo(200, 0);
 +
ctx.lineTo(300, 0);
 +
ctx.stroke();
 +
 
 +
for (var i = 0; i < n; i++) {
 +
 
 +
// Рисование пружинок
 +
if (i == 0) {
 +
if (!spring_dropped)
 +
draw_spring(0, circle[i].y, w / 2, 10, 40);
 +
} else
 +
draw_spring(circle[i - 1].y, circle[i].y, w / 2, 10, 40);
 +
 
 +
 
 +
}
 +
for (var i = 0; i < n; i++) {
 
 
−
function control1() {
+
// Рисование грузиков
−
physics1();
+
ctx.lineWidth = 6;
−
draw();
+
ctx.strokeStyle = "#8B008B";
 +
ctx.beginPath();
 +
ctx.moveTo(circle[i].x-20, circle[i].y);
 +
ctx.lineTo(circle[i].x+20, circle[i].y);
 +
ctx.stroke();
 
}
 
}
  
−
function physics1() {                   // то, что происходит каждый шаг времени
+
}
−
        for (var s = 1; s <= spf; s++) {
+
 
−
for (i = 0; i < n; i++) {
+
// *** Функция рисования пружины ***
−
for (j = 0; j < n; j++) {
+
function draw_spring(x_start, x_end, y, n, h) {
−
runge();
+
ctx.lineWidth = 1;
−
leapfrog();
+
 
−
}
+
var L = x_end - x_start;
−
}
+
for (var i = 0; i < n; i++) {
−
        }
+
var rgb = HSVtoRGB(i / n * 0.5, 1, 1);
−
    }
+
ctx.strokeStyle = "rgba(" + rgb.r + ", " + rgb.g + ", " + rgb.b + ", 1)";
−
 
−
function leapfrog(){
 
−
particles1[i][j].vx = particles1[i][j].vx - particles1[i][j].x * dt;         
 
−
particles1[i][j].x = particles1[i][j].x + particles1[i][j].vx * dt;
 
−
console.log(particles1[0][0].x);
 
−
        }
 
−
 
−
function runge() {
 
−
var k1, k2, k3, k4;
 
−
k1 = - (particles[i][j].x);
 
−
k2 = - (particles[i][j].x + k1 * dt / 2);
 
−
k3 = - (particles[i][j].x + k2 * dt / 2);
 
−
k4 = - (particles[i][j].x + k3 * dt);
 
−
particles[i][j].vx = particles[i][j].vx + dt / 6 * (k1 + 2 * k2 + 2 * k3 + k4);
 
−
k1 = (particles[i][j].vx);
 
−
k2 = (particles[i][j].vx + k1 * dt / 2);
 
−
k3 = (particles[i][j].vx + k2 * dt / 2);
 
−
k4 = (particles[i][j].vx + k3 * dt);
 
−
particles[i][j].x = particles[i][j].x + dt / 6 * (k1 + 2 * k2 + 2 * k3 + k4);
 
−
console.log(particles[0][0].x);
 
−
}
 
  
−
function draw() {
+
var x_st = x_start + L / n * i;
−
context.clearRect(0, 0, 600, 600);
+
var x_end = x_start + L / n * (i + 1);
−
context.fillStyle="#000000";
+
var l = x_end - x_st;
−
context.beginPath();
+
ctx.beginPath();
−
context.moveTo(0, 300);
+
ctx.bezierCurveTo(y, x_st, y + h, x_st + l / 4, y, x_st + l / 2);
−
context.lineTo(600, 300);
+
ctx.bezierCurveTo(y, x_st + l / 2, y - h, x_st + 3 * l / 4, y, x_st + l);
−
context.stroke();
+
ctx.stroke();
−
context.beginPath();
 
−
context.moveTo(300, 0);
 
−
context.lineTo(300, 600);
 
−
context.stroke();
 
−
context1.clearRect(0, 0, 600, 600);
 
−
context1.fillStyle="#000000";
 
−
context1.beginPath();
 
−
context1.moveTo(0, 300);
 
−
context1.lineTo(600, 300);
 
−
context1.stroke();
 
−
context1.beginPath();
 
−
context1.moveTo(300, 0);
 
−
context1.lineTo(300, 600);
 
−
context1.stroke();
 
−
for (var i = 0; i < n; i++) {
 
−
for (var j = 0; j < n; j++) {
 
−
context.beginPath();
 
−
context.arc(particles[i][j].x * 300 + 300, 300 - particles[i][j].vx * 300, 1, 0, 2 * Math.PI, false);                  //рисуем шар
 
−
context.fill();
 
−
context.closePath();
 
−
 
−
context1.beginPath();
 
−
context1.arc(particles1[i][j].x * 300 + 300, 300 - particles1[i][j].vx * 300, 1, 0, 2 * Math.PI, false);                 //рисуем шар
 
−
context1.fill();
 
−
context1.closePath();
 
−
}
 
−
}
 
−
}  
 
 
}
 
}
−
</syntaxhighlight>
+
}
−
Файл '''"Chain_v3_release.html"'''
+
}
−
<syntaxhighlight lang="html5" line start="1" enclose="div">
 
−
<!DOCTYPE html>
 
−
    <html>
 
−
    <head>
 
−
        <meta charset="UTF-8"/>
 
−
        <title> Integrate </title>
 
−
        <script src="Integral3.js"></script>
 
−
    </head>
 
−
    <body>
 
−
        <canvas id="Integrate_canvas" width="600" height="600" style="border:1px solid #000000;"></canvas>
 
−
<canvas id="Integrate_canvas1" width="600" height="600" style="border:1px solid #000000;"></canvas>
 
−
            <div>
 
−
                            <!-- n-->
 
−
            <div>
 
−
                n =
 
−
            <input id="Text_n" style="width: 4.2ex;" required pattern="[-+]?([0-9]*\.[0-9]+|[0-9]+)" oninput="
 
−
              if (!this.checkValidity()) return;
 
−
              app.set_n(this.value);
 
−
              document.getElementById('Slider_n').value = this.value;
 
−
          ">
 
−
                    <input type = "range"  id="Slider_n" style="width: 100px;" oninput="app.set_n(this.value); document.getElementById('Text_n').value = this.value;">
 
−
            </I></font>
 
−
        </div>
 
−
        </div>   
 
−
            <script type="text/javascript"> app = new Main_Integrate(document.getElementById('Integrate_canvas'),document.getElementById('Integrate_canvas1'));</script>   
 
−
    </body>
 
−
    </html>
 
 
</syntaxhighlight>
 
</syntaxhighlight>
 
</div>
 
</div>
 
</div>
 
</div>
  
 +
=='''Физическое объяснение'''==
 +
Мы наблюдаем иллюзию "зависания" нижнего грузика потому, что до него должна дойти волна возмущений верхней пружины, прежде чем он придет в движение.
 +
 +
=='''Ссылки'''==
 +
* [[Виртуальная лаборатория]]
 +
<br/>
 
[[Category: Виртуальная лаборатория]]
 
[[Category: Виртуальная лаборатория]]
 
[[Category: Программирование]]
 
[[Category: Программирование]]
 
[[Category: JavaScript]]
 
[[Category: JavaScript]]

Версия 20:35, 12 января 2016

Виртуальная лаборатория>Сравнение методов Рунге-Кутта и Липфрога

Постановка задачи

Дано простейшее уравнение движения грузика на пружине: [math]\ddot{x} = -cx[/math]. Необходимо интегрировать его с помощью двух методов: Липфрога и Рунге-Кутта. Сравнить результаты.

Метод численного интегрирования Рунге-Кутты

Рассмотри задачу Коши для системы обыкновенных дифференциальных уравнений первого порядка. [math]\dot{y} = f(x, y), y(x_0)=y_0[/math] Тогда приближенное значение в последующих точках вычисляется по итерационной формуле: [math]y_n+1 = y_n + \frac{h}{6}(k_1 + 2k_2 + 2k_3 + k_4)[/math]

Основные уравнения

1)Расчет статических координат:

[math]{y}_{i}=i{l}_{0}+\sum^{i}_{k=1} {\Delta l_k}[/math]

[math]\Delta{l_k} = \frac{(n-k+1)mg}{C} [/math]

где [math]{n}[/math] — количество грузиков ; [math]{i}[/math] — текущий грузик ; [math]{l_0}[/math] — длина не растянутой пружины; [math]\Delta{l_k}[/math] — статическое удлинение пружины.

2) Уравнение движения грузиков при отсутствии верхней пружины:

[math]\ddot y_i = - \frac{C}{m} ({y}_{i+1} - 2{y}_{i}+{y}_{i-1}) + g [/math]

при [math]i = n [/math] :
[math] \ddot y_n = - \frac{C}{m}(y_n - {y}_{n-1} - l_0) + g[/math]

Программа (лучше всего смотреть ее в Mozilla Firefox)


Скачать программу Spring_illusion.rar.

Текст программы на языке JavaScript (разработчик Александров Сергей):

Файл "Spring_illusin.js" 

  1 window.addEventListener("load", Main_Spring, true);
  2 function Main_Spring() {
  3 
  4 	// *** Некие исходные данные ***
  5 
  6 	var canvas = spring_canvas;
  7 	canvas.onselectstart = function () {return false;}; // запрет выделения canvas
  8 	var ctx = canvas.getContext("2d"); // на ctx происходит рисование
  9 	var w = canvas.width; // ширина окна в расчетных координатах
 10 	var h = canvas.height; // высота окна в расчетных координатах
 11 
 12 	var Pi = 3.1415926; // число "пи"
 13 	var g = 9.8; // гравитационная постоянная
 14 
 15 	var T0 = 0.01; // масштаб времени (период колебаний исходной системы)
 16 
 17 	var m0 = 0.15; // масштаб массы маятника
 18 	var C0 = 1; // масштаб жесткости пружины
 19 
 20 	var count = true; // проводить ли расчет системы
 21 	var v = 0; // скорость тела
 22 	var t = 0;
 23 
 24 	// параметры полученные из размеров холста
 25 	var rw = canvas.width / 4;
 26 	var rh = canvas.height / 100;
 27 
 28 	// параметры грузиков
 29 	var x0 = 0;
 30 	var y0 = 0;
 31 	var vy0 = 0;
 32 	var rad0 = 7;
 33 	var sTime = 0;
 34 	// параметры пружины
 35 	var Lp = 30; //длина пружины
 36 
 37 
 38 	// *** Задание вычислительных параметров ***
 39 
 40 	var fps = 60; // frames per second - число кадров в секунду (качечтво отображения)
 41 	var spf = 5; // steps per frame   - число шагов интегрирования между кадрами
 42 	var dt = 50 * T0 / fps; // шаг интегрирования (качество расчета)
 43 	var steps = 0; // количество шагов интегрирования
 44 
 45 
 46 	// *** Задание физических параметров ***
 47 	var m = 1 * m0; // масса грузика
 48 	var C = 1 * C0; // жесткость пружины
 49 	var L0 = 0; // изначальное удлинение пружины
 50 	var n = 5; // количество тел
 51 	var circle; // круглые тела
 52 
 53 	// *** Установка слайдеров для переменных величин ***
 54 	slider_m.value = (m / m0).toFixed(1);
 55 	number_m.value = (m / m0).toFixed(1);
 56 	slider_spf.value = (spf).toFixed(1);
 57 	number_spf.value = (spf).toFixed(1);
 58 	slider_n.value = parseInt(n);
 59 	number_n.value = parseInt(n);
 60 
 61 	function setM(new_m) {
 62 		m = new_m * m0;
 63 	}
 64 
 65 	function setSpf(new_spf) {
 66 		spf = new_spf;
 67 	}
 68 	function setN(new_n) {
 69 		n = new_n;
 70 		StaticStart();
 71 	}
 72 
 73 	slider_m.oninput = function () {
 74 		number_m.value = slider_m.value;
 75 		setM(slider_m.value);
 76 	};
 77 	number_m.oninput = function () {
 78 		slider_m.value = number_m.value;
 79 		setM(number_m.value);
 80 	};
 81 
 82 
 83 	slider_spf.oninput = function () {
 84 		number_spf.value = slider_spf.value;
 85 		setSpf(slider_spf.value);
 86 	};
 87 	number_spf.oninput = function () {
 88 		slider_spf.value = number_spf.value;
 89 		setSpf(number_spf.value);
 90 	};
 91 
 92 	slider_n.oninput = function () {
 93 		number_n.value = slider_n.value;
 94 		setN(slider_n.value);
 95 	};
 96 	number_n.oninput = function () {
 97 		slider_n.value = number_n.value;
 98 		setN(number_n.value);
 99 	};
100 
101 	// *** Условие падения тела ***
102 	var spring_dropped = false;
103 	drop_spring.onclick = function () {
104 		spring_dropped = true;
105 	};
106 
107 	// *** Пересчет статических координат
108 	new_static.onclick = function () {
109 		StaticStart()
110 	};
111 
112 	// *** Создание массива элементов ***
113 	function StaticStart() {
114 		spring_dropped = false;
115 		circle = [];
116 		for (var i = 0; i < n; i++) {
117 			var circ = {};
118 			circ.x = x0;
119 			circ.y = y0;
120 			circ.vy = vy0;
121 			circ.rad = rad0;
122 			circ.L = L0;
123 
124 			var rgb = HSVtoRGB(i / n * 2, 1, 1);
125 
126 			circ.fill = "rgba(" + rgb.r + ", " + rgb.g + ", " + rgb.b + ", 1)";
127 			circle[i] = circ;
128 		}
129 		for (var i = 0; i < n; i++) {
130 			if (i == 0) {
131 				circle[i].x = w / 2;
132 				circle[i].y0 = Lp + (n - i) * m * g / C;
133 				circle[i].y = circle[i].y0;
134 
135 			} else {
136 				circle[i].x = w / 2;
137 				circle[i].y0 = circle[i - 1].y0 + Lp + (n - i) * m * g / C;
138 				circle[i].y = circle[i].y0;
139 				console.log(circle[i].y0);
140 
141 			}
142 		}
143 	}
144 	
145 // *** функция разукрашивания объектов
146 function HSVtoRGB(h, s, v) {              
147     var r, g, b, i, f, p, q, t;
148     i = Math.floor(h * 6);
149     f = h * 6 - i;
150     p = v * (1 - s);
151     q = v * (1 - f * s);
152     t = v * (1 - (1 - f) * s);
153     switch (i % 6) {
154         case 0: r = v, g = t, b = p; break;
155         case 1: r = q, g = v, b = p; break;
156         case 2: r = p, g = v, b = t; break;
157         case 3: r = p, g = q, b = v; break;
158         case 4: r = t, g = p, b = v; break;
159         case 5: r = v, g = p, b = q; break;
160     }
161     return {
162         r: Math.floor(r * 255),
163         g: Math.floor(g * 255),
164         b: Math.floor(b * 255)
165     };
166 }
167 
168 	    // график
169     var vGraph1 = new TM_graph(                  // определить график
170         "#vGraph1",                              // на html-элементе #vGraph
171         250,                                    // сколько шагов по оси "x" отображается
172         0, 100, 5);                            // мин. значение оси Y, макс. значение оси Y, шаг по оси Y
173 
174 // *** Функция обеспечивающая "жизнь" пружин ***
175 function control() {
176 	calculate();
177 	draw();
178 	requestAnimationFrame(control);
179 }
180 StaticStart();
181 control();
182 
183 
184 
185 // *** Функция расчетов координат ***
186 function calculate() {
187 
188 	for (var s = 1; s <= spf; s++) {
189 	var k = n;
190 		for (var i = 1; i < n; i++) {
191 			if (spring_dropped)
192 				circle[0].L = 0;
193 			else
194 				circle[0].L = circle[0].y - Lp;
195 			circle[i].L = circle[i].y - circle[i - 1].y - Lp;
196 		}
197 
198 		for (var i = 0; i < n; i++) {
199 			if (i == n - 1)
200 				circle[i].vy += ((-1) * C * (circle[i].L) / m + g) * dt;
201 			else
202 				circle[i].vy += ((-1) * C * (circle[i].L - circle[i + 1].L) / m + g) * dt;
203 				
204 			
205 			circle[i].y += circle[i].vy * dt;
206 
207 		}
208 		steps++;
209 
210 		if (spring_dropped)
211 		{ sTime ++;
212 		if (sTime % 50 == 0) vGraph1.graphIter(sTime/50, (circle[n-1].y-circle[n-1].y0)) ;
213 		};
214 	}
215 }
216 
217 
218 // *** Функция рисования объектов ***
219 function draw() {
220 	ctx.clearRect(0, 0, w, h);
221 
222 	// Рисование закрепления
223 	ctx.lineWidth = 6;
224 	ctx.strokeStyle = "#7394cb";
225 	ctx.beginPath();
226 	ctx.moveTo(200, 0);
227 	ctx.lineTo(300, 0);
228 	ctx.stroke();
229 
230 	for (var i = 0; i < n; i++) {
231 
232 		// Рисование пружинок
233 		if (i == 0) {
234 			if (!spring_dropped)
235 				draw_spring(0, circle[i].y, w / 2, 10, 40);
236 		} else
237 			draw_spring(circle[i - 1].y, circle[i].y, w / 2, 10, 40);
238 
239 
240 	}
241 	for (var i = 0; i < n; i++) {
242 	
243 		// Рисование грузиков	
244 		ctx.lineWidth = 6;
245 		ctx.strokeStyle = "#8B008B";
246 		ctx.beginPath();
247 		ctx.moveTo(circle[i].x-20, circle[i].y);
248 		ctx.lineTo(circle[i].x+20, circle[i].y);
249 		ctx.stroke();
250 	}
251 
252 }
253 
254 // *** Функция рисования пружины ***
255 function draw_spring(x_start, x_end, y, n, h) {
256 	ctx.lineWidth = 1;
257 
258 	var L = x_end - x_start;
259 	for (var i = 0; i < n; i++) {
260 		var rgb = HSVtoRGB(i / n * 0.5, 1, 1);
261 		ctx.strokeStyle = "rgba(" + rgb.r + ", " + rgb.g + ", " + rgb.b + ", 1)";
262 
263 		var x_st = x_start + L / n * i;
264 		var x_end = x_start + L / n * (i + 1);
265 		var l = x_end - x_st;
266 		ctx.beginPath();
267 		ctx.bezierCurveTo(y, x_st, y + h, x_st + l / 4, y, x_st + l / 2);
268 		ctx.bezierCurveTo(y, x_st + l / 2, y - h, x_st + 3 * l / 4, y, x_st + l);
269 		ctx.stroke();
270 	}
271 }
272 }

Физическое объяснение

Мы наблюдаем иллюзию "зависания" нижнего грузика потому, что до него должна дойти волна возмущений верхней пружины, прежде чем он придет в движение.

Ссылки


Источник — «http://tm.spbstu.ru/?title=Сравнение_методов_Рунге-Кутта_и_Липфрога&oldid=38741»