Редактирование: Задача падающей цепочки

'''Курсовой проект по [[Введение в механику дискретных сред]]'''''

'''Исполнитель:''' [[Борисенков Богдан]]

'''Группа:''' 5030103/90101

'''Семестр:''' осень 2022

===Постановка задачи===
В рамках проекта необходимо смоделировать движение одномерной цепочки: начальное положение (провисание) цепочки и дальнейшее ее падение при отпускании одного из концов под действием силы тяжести, а также исследовать зависимость ускорения крайней свободной частицы от времени.

===Математическая модель ===
Изначально запишем закон движения:
<math>
m\underline{\ddot{r}}_i(t)=\underline{F}_{i-1}+\underline{F}_{i+1} + \underline{F}_{g}\\
\underline{r}_i(0)=\underline{r}_i^0,~\underline{v}_i(0)=0~~~i=1,\ldots,n
</math>

где
<math>
\underline{F}_{i-1}, \underline{F}_{i+1}\\
</math> - силы упругости действующие на <math>i</math>-ую частицу со стороны <math>i-1</math> и <math>i+1</math> соответственно, а <math> \underline{F}_{g}=-mg\underline{k} \\ </math> - сила тяжести.

Далее распишем силу упругости как произведение модуля на соответствующий орт:
<math>
\underline{F}_{i+1}= c(|\underline{r}_{i+1}-\underline{r}_{i}| - l_0)\frac{(\underline{r}_{i+1}-\underline{r}_{i})}{|\underline{r}_{i+1}-\underline{r}_{i}|}
</math>, где <math>c</math> - коэффициент жесткости пружины.
Аналогично записывается сила <math>\underline{F}_{i-1}</math>.

Далее подставляя все силы в уравнение движения, получим:

<math>
m\underline{\ddot{r}}_i(t)= c(||\underline{r}_{i+1}-\underline{r}_i|| -l_0)\frac{(\underline{r}_{i+1}-\underline{r}_i)}{||\underline{r}_{i+1}-\underline{r}_i||} + c(||\underline{r}_{i-1}-\underline{r}_i|| - l_0)\frac{(\underline{r}_{i-1}-\underline{r}_i)}{||\underline{r}_{i-1}-\underline{r}_i||} - mg\underline{k}\\
</math>

Дальнейшее интегрирование уравнения производится с помощью явного симплектического метода Верле c нулевыми начальными условиями и условиями закрепления на концах.

<math> \begin{cases}
V_{i+1} = V_i+A_i\Delta{t}\\
X_{i+1} = X_i+V_{i+1}\Delta{t},
\end{cases} </math>
===Численное моделирование===
{{#widget:Iframe |url = http://tm.spbstu.ru/htmlets/js2020/Borisenkov/Falling_Chain/Chain2.html| width= 1200 | height = 600}}
===Выводы===

В рамках решения задачи смоделировано движение цепочки под действием силы тяжести и проиллюстрирован тот факт, что ускорение крайней массы цепочки больше, чем ускорение свободно падающего тела. Данный эффект объясняется начальным преднатяжением цепочки.
==Ссылки==
http://tm.spbstu.ru/Курсовые_работы_по_ВМДС:_2022-2023 - курсовые работы студентов 4-го курса 2022-2023 года по курсу дискретной механики
http://tm.spbstu.ru/Введение_в_механику_дискретных_сред - курс механики дискретных сред
==Код программы==
<div class="mw-collapsible mw-collapsed">
'''Код программы на языке JavaScript (разработчик Богдан Борисенков):''' <div class="mw-collapsible-content">
<syntaxhighlight lang="javascript" line start="1" enclose="div">
window.addEventListener('load',main,false);
function main() { 
	
	var N = document.getElementById('Num').value;
	var a = document.getElementById('DIST').value;
	var fps = document.getElementById('FPS').value;
	var COEF = document.getElementById('COEFF').value;
	var hop,g,b,dt,cm;
	// считаем, что у нас массы всех частиц и жесткости одинаковы
	
	PRT = []; // частицы
	var l_akt1, l_akt2, vx_dot, vy_dot;
	var dt2;
	var sp = 1;
	var flag = 0;
	
	var time = 0; // нужна для подсчета шагов по времени
	var time2 = 0;
	var M = [];
	M.m = 10;
	

	// канвасы
	var ctx = cnv.getContext('2d');
	var h = cnv.height;
	var w = cnv.width;
	var ctxG = cnv_graf.getContext('2d');
	var hG = cnv_graf.height;
	var wG = cnv_graf.width;
	var ctxG2 = cnv_graf2.getContext('2d');
	var hG2 = cnv_graf2.height;
	var wG2 = cnv_graf2.width;
	
	// для рисования цепочки по центру канваса
	var SHIFT; 
	var global_len;
	
	// массив модуля ускорения
	var VDOT = [];
	dot_speed2 = dot_speed = 0;
	RAZN = RAZN2 = 0;
	x_step = 0;
	max_value = 1;
	max_value2 = 1;
	y_tick(5);
	y_tick2(5);
	MASSIVE = [];
	// кнопки старт, Обновить, пауза
	
	// пауза
	function switchdt() {
		switch (flag) {
			case 0: {
				flag = 1;
				sp = 0;
				break;
			}
			case 1: {
				flag = 0;
				sp = 1;

			}
		}
	}	
	COEFF.oninput = function () {UpdCoeff();}
	function UpdCoeff() { COEF = parseFloat(document.getElementById('COEFF').value); SHIFT = (w - 0.65*(N-1)*COEF*parseFloat(a))/2;}
	Pause.onclick = function () { switchdt();}
	
	// Обновить 
	New.onclick = function() {
		clearInt(interv);
		Update();
	}
	
	LET.onclick = function () { 
		clearInt(interv);
		interv = setInterval(control2,1000/fps);	
	}
	
	Add.onclick = function () { 
		// вносит изменения в данные
		UpdateCoeffs();
	}
	// Функция, которая обновляет данные и запускает расчет
	function Update() { 
		// берем значение для N
		N = parseInt(document.getElementById('Num').value);
		cm = parseFloat(document.getElementById('CM').value);
		a = parseFloat(document.getElementById('DIST').value);
		dt0 = parseFloat(document.getElementById('DT').value);
		Betta = parseFloat(document.getElementById('B').value);
		g = parseFloat(document.getElementById('G').value);
		fps = parseInt(document.getElementById('FPS').value);
		hop = parseInt(document.getElementById('HOP').value);
		// вспомогательная константа
		dt2 = dt0*dt0/2;
		// теперь задаем параболу
		count(N,a);
		time = 0;
		time2 = 0;
		MASSIVE = [];
		a = Math.pow(Math.pow(PRT[2].y - PRT[1].y,2) + Math.pow(PRT[2].x - PRT[1].x,2),1/2)/2;
		interv = setInterval(control1,1000/fps);
	}
	
	function count(NUM, dist) {
		var len = 0.65*(NUM-1)*dist;
		global_len = len;
		var distt = len/(NUM-1);
		var len_2 = len/2;
		var w_2 = w/2;
		var constt2 = 2*Math.pow(2,1/2);
		var constt = constt2/len;
		
		for (var i = 0; i< NUM; i++) { 
			b = [];
			b.x = distt*i; 
			b.y = -constt*Math.pow(b.x-len,2) - constt2*(b.x-len);
			b.vx = 0;
			b.vy = 0;
			b.vx_dot = 0;
			b.vy_dot = 0;
			PRT[i] = b;
		}
		// на всякий
		PRT[0].y = 0;
		PRT[NUM-1].y = 0;
		// console.log(PRT);
		M.x = PRT[N-1].x;
		M.y = PRT[N-1].y;
		M.vy = 0;
		
		// сместим цепочку в центр канваса
		UpdCoeff();
	}
		
	
	function phys_1() {  // когда последняя частица закреплена
		// первая частица у нас всегда закреплена и последняя!
		dt = dt0*sp;
		dot_speed = dot_speed2;
		for (var i = 1; i<N-1; i++) { 
			
			l_akt1 = Math.pow(Math.pow(PRT[i+1].x - PRT[i].x,2) + Math.pow(PRT[i+1].y - PRT[i].y,2),1/2);
			l_akt2 = Math.pow(Math.pow(PRT[i].x - PRT[i-1].x,2) + Math.pow(PRT[i].y - PRT[i-1].y,2),1/2);
			// занулим силы, если расстояние меньше а
			if (l_akt1 < a) {
				FR = 0;
			} else { 
				FR = (l_akt1 - a);
			}
			if (l_akt2 < a) {
				FL = 0;
			} else { 
				FL = (l_akt2 - a);
			}
			vx_dot = cm*(FR*(PRT[i+1].x - PRT[i].x)/l_akt1 - FL*(PRT[i].x - PRT[i-1].x)/l_akt2) - Betta*PRT[i].vx;
			vy_dot = (cm*(FR*(PRT[i+1].y - PRT[i].y)/l_akt1 - FL*(PRT[i].y - PRT[i-1].y)/l_akt2) + g) - Betta*PRT[i].vy;
			PRT[i].vx_dot = vx_dot;
			PRT[i].vy_dot = vy_dot;
			PRT[i].vx += vx_dot*dt ;
			PRT[i].vy += vy_dot*dt;

		}
		for (var i = 1; i<N-1; i++) { 
			PRT[i].x += PRT[i].vx*dt;
			PRT[i].y += PRT[i].vy*dt;
		}
		
	}
	
	function phys_2() { // когда последнюю частицу отпустили
		dt = dt0*sp;
		dot_speed = dot_speed2;
		RAZN = RAZN2;
		MASSIVE.push(RAZN);
		for (var i = 1; i<N-1; i++) { 
			l_akt1 = Math.pow(Math.pow(PRT[i+1].x - PRT[i].x,2) + Math.pow(PRT[i+1].y - PRT[i].y,2),1/2);
			l_akt2 = Math.pow(Math.pow(PRT[i].x - PRT[i-1].x,2) + Math.pow(PRT[i].y - PRT[i-1].y,2),1/2);
			if (l_akt1 < a) {
				FR = 0;
			} else { 
				FR = (l_akt1 - a);
			}
			if (l_akt2 < a) {
				FL = 0;
			} else { 
				FL = (l_akt2 - a);
			}
			vx_dot = cm*(FR*(PRT[i+1].x - PRT[i].x)/l_akt1 - FL*(PRT[i].x - PRT[i-1].x)/l_akt2) - Betta*PRT[i].vx;
			vy_dot = (cm*(FR*(PRT[i+1].y - PRT[i].y)/l_akt1 - FL*(PRT[i].y - PRT[i-1].y)/l_akt2) + g) - Betta*PRT[i].vy;
			PRT[i].vx_dot = vx_dot;
			PRT[i].vy_dot = vy_dot;
			PRT[i].vx += vx_dot*dt ;
			PRT[i].vy += vy_dot*dt;

		}
		// теперь для последней частицы
		l_akt2 = Math.pow(Math.pow(PRT[N-1].x - PRT[N-2].x,2) + Math.pow(PRT[N-1].y - PRT[N-2].y,2),1/2);
		if (l_akt2 < a) {
				FL = 0;
			} else { 
				FL = (l_akt2 - a);
		}
		vx_dot =  -cm*FL*(PRT[N-1].x - PRT[N-2].x)/l_akt2 - Betta*PRT[N-1].vx;
		vy_dot =  -cm*FL*(PRT[N-1].y - PRT[N-2].y)/l_akt2 - Betta*PRT[N-1].vy;
		PRT[N-1].vx_dot = vx_dot;
		PRT[N-1].vy_dot = vy_dot;
		PRT[N-1].vx += vx_dot*dt ;
		PRT[N-1].vy += vy_dot*dt;
		
		// теперь считаем новые координаты

		for (var i = 1; i<N; i++) { 
			PRT[i].x += PRT[i].vx*dt;
			PRT[i].y += PRT[i].vy*dt;
		}
		M.y = g*Math.pow(time2*dt0,2)/2;
	}

	
	function draw1() { 
		ctx.clearRect(0,0,w,h);
		for (var i = 0; i<N; i++) { 
			ctx.beginPath();
			ctx.arc(PRT[i].x*COEF + SHIFT,PRT[i].y*COEF + 50, 3, 0, 2*Math.PI);
			ctx.stroke();
		}
		ctx.beginPath();
		ctx.moveTo(0,50);
		ctx.lineTo(w,50);
		ctx.stroke();
		
	}
	
	function draw2() { 
		ctx.clearRect(0,0,w,h);
		for (var i = 0; i<N; i++) { 
			ctx.beginPath();
			ctx.arc(PRT[i].x*COEF + SHIFT,PRT[i].y*COEF + 50, 3, 0, 2*Math.PI);
			ctx.stroke();
		}
		ctx.beginPath();
		ctx.moveTo(0,50);
		ctx.lineTo(w,50);
		ctx.stroke();
		// свободно падающее тело 
		ctx.beginPath();
		ctx.arc(M.x*COEF + SHIFT,M.y*COEF + 50, 5, 0, 2*Math.PI);
		ctx.fill();
		
	}
	
	function draw_graf() { 
		// строим график ускорений 
		dot_speed2 = Math.pow(Math.pow(PRT[N-1].vx_dot,2) + Math.pow(PRT[N-1].vy_dot,2),1/2)/g;
		RAZN2 = PRT[N-1].y - M.y;
		// надо понять масштаб графика по y
		if (dot_speed2/max_value > 0.9) { max_value = dot_speed2*1.5; y_tick(5);}
		if (RAZN2/max_value2 > 0.9) { max_value2 = RAZN2*1.5; y_tick2(5);}
		if (x_step == wG) { x_step = 0; ctxG.clearRect(0,0,wG,hG); ctxG2.clearRect(0,0,wG2,hG2); y_tick(5); y_tick2(5);}
		ctxG.beginPath();
		ctxG.moveTo(x_step, hG*(1 - dot_speed/max_value));
		x_step += 0.5;
		ctxG.lineTo(x_step,hG*(1 - dot_speed2/max_value));
		ctxG.stroke();
		// рисуем верх 
		ctxG2.clearRect(0,0,wG2,hG2);
		y_tick2(5);
		xx_step = wG2/2/MASSIVE.length;
		for (var j = 1; j < MASSIVE.length; j++) { 
			ctxG2.beginPath();
			ctxG2.moveTo(xx_step*(j-1), hG2*(1 - MASSIVE[j-1]/max_value2));
			ctxG2.lineTo(xx_step*j,hG2*(1 - MASSIVE[j]/max_value2));
			ctxG2.stroke();
		}
		

		
	}
	function y_tick(num) { 
		ctxG.clearRect(0,0,50,hG-10);
		ctxG.fillText(max_value.toFixed(2),5,0.05*hG);
		var step = hG/num;
		for (var i = 1; i< num; i++) { 
			// рисуем засечку
			val = step*i;
			ctxG.beginPath();
			ctxG.moveTo(0,hG - val); 
			ctxG.lineTo(3,hG - val);
			ctxG.fillText((val/hG*max_value).toFixed(2), 5, hG - val);
			ctxG.stroke();
		}
		ctxG.fillText('x"/G',40,0.05*hG);
		

	}
	function y_tick2(num) { 
		ctxG2.clearRect(0,0,50,hG2-10);
		ctxG2.clearRect(0,0,80,20);
		ctxG2.fillText(max_value2.toFixed(2),5,0.05*hG2);
		var step = hG/num;
		for (var i = 1; i< num; i++) { 
			// рисуем засечку
			val = step*i;
			ctxG2.beginPath();
			ctxG2.moveTo(0,hG - val); 
			ctxG2.lineTo(3,hG - val);
			ctxG2.fillText((val/hG2*max_value2).toFixed(2), 5, hG2 - val);
			ctxG2.stroke();
		}
		ctxG2.fillText('y - yM',40,0.05*hG2);
	}
	function control1() { 
		phys_1();
		if (time % hop == 0) { 
			draw1();
		}
		draw_graf();
		time++;
	}
	function control2() { 
		phys_2();
		if (time % hop == 0) { 
			draw2();
		}
		draw_graf();
		time2 += sp;
	}
	
	function clearInt(intrv) { 
		clearInterval(intrv);
	}
	FPS.oninput = function () { 
		fps = parseInt(document.getElementById('FPS').value);
		clearInt(interv);
		interv = setInterval(control2,1000/fps);
	}
	
	function UpdateCoeffs() {
		cm = parseFloat(document.getElementById('CM').value);
		dt0 = parseFloat(document.getElementById('DT').value);
		Betta = parseFloat(document.getElementById('B').value);
		g = parseFloat(document.getElementById('G').value);
		hop = parseInt(document.getElementById('HOP').value);
		dt2 = dt0*dt0/2;
	}
	UpdateCoeffs();
	Update();

}
</syntaxhighlight>
</div>
<div class="mw-collapsible mw-collapsed">
'''html - файл:''' <div class="mw-collapsible-content">
<syntaxhighlight lang="html" line start="1" enclose="div">
<!DOCTYPE html>
<html>
	<head>
		<title> Chain_2022 </title>
		<script src = 'Chain2.js'>
		</script>
	</head>
	<body>
		<div id = 'Canv'>
		<canvas id = 'cnv' width = 500 height = 500 style='border: 1px solid black;'></canvas>
		<canvas id = 'cnv_graf' width = 500 height = 250 style='border: 1px solid black;'></canvas>
		<canvas id = 'cnv_graf2' width = 500 height = 250 style='border: 1px solid black;'></canvas>
		<b>
		<br><label>C/M <input type = 'text' id = 'CM' value = '2000'></label>
		<label>DT <input type = 'text' id = 'DT' value = '0.01'></label>
		<label>FPS <input type = 'range' id = 'FPS' min = 10 step = 10 value = 120 max = 120></label>
		 <input type = "button" id = "New" value = "Обновить"/>
		 <input type="button" id = "Pause" value = "▶/||"/> <input type="button" id = "Add" value = "Внести правки"/>
		 <label> Расстояние между частицами <input type = 'text' id = 'DIST' value = '1'>
		<br>
		<br><label>B <input type = 'text' id = 'B' value = '1'></label>
		<label>N <input type = 'text' id = 'Num' value = '100'></label>
		<label>G <input type = 'text' id = 'G' value = '50'></label>
		<label>Рисовать каждый <input type = 'text' id = 'HOP' value = 1> шаг</label> 
		<input type = 'button' id = 'LET' value = 'Отпустить правый конец'>
		<label>Масштаб<input type = 'range' id = 'COEFF' min = 0.1 step = 0.01 value = 3 max = 20></label>
		<br> 
		</b>
		</div>
		
		<style type="text/css">
			#CM,#DT,#B,#G,#Num,#HOP {width:70px}
			<!-- #Canv{  -->
				<!-- width: 50%; -->
				<!-- margin: 0 auto; -->
				<!-- } -->
			#cnv_graf2 {position:absolute; left:514px;}
			#cnv_graf {position:absolute; left:514px; top:258px;}
		</style>
		
	</body>
</html>
</syntaxhighlight>
</div>