Вынужденные колебания цепочки в вязкой среде — различия между версиями

Материал из Department of Theoretical and Applied Mechanics
Перейти к: навигация, поиск
(→‎Теоретическая сводка)
(→‎Теоретическая сводка)
 
(не показано 14 промежуточных версий 2 участников)
Строка 12: Строка 12:
  
 
==Построение модели==
 
==Построение модели==
−
Процесс моделируется как одномерные колебания цепочки частиц.  
+
Процесс моделируется как одномерные колебания цепочки частиц.
−
Уравнение взаимодействия : (строчка с F)
 
−
 
 
 
==Теоретическая сводка==
 
==Теоретическая сводка==
  
 
Парное взаимодействие определяется формулой:
 
Парное взаимодействие определяется формулой:
  
−
:<math>F = c(x_{k+1}-x_{k}-1)</math>
+
:<math>F = c(x_{k+1}-x_{k})</math>
  
 
где <math>F</math> — сила взаимодействия, <math>c</math> — жесткость связи, <math>x_k</math> — перемещение частицы, <math>k</math> — номер частицы.
 
где <math>F</math> — сила взаимодействия, <math>c</math> — жесткость связи, <math>x_k</math> — перемещение частицы, <math>k</math> — номер частицы.
Строка 29: Строка 27:
 
где <math>F</math> — сила взаимодействия, <math>m</math> — масса частицы, <math>μ</math> — коэффициент вязкости, <math>V_{k}</math> — скорость частицы.
 
где <math>F</math> — сила взаимодействия, <math>m</math> — масса частицы, <math>μ</math> — коэффициент вязкости, <math>V_{k}</math> — скорость частицы.
  
−
==Решение==
+
===Граничные и начальные условия===
−
При построении модели были приняты следующие значения параметров:
+
Начальные условия нулевые
 +
 
 +
Периодические граничные условия:
 +
:<math>V_{N}=c(x_{1}+aN-x_{N}-1)\frac{dt}{m}-2μV_{N}+Qcos(ωt)</math>
 +
:<math>V_{1}=-c(x_{1}+aN-x_{N}-1)\frac{dt}{m}-2μV_{1}</math>
 +
где <math>Q</math> — амплитуда возмущающей силы, <math>ω</math> — частота возмущающей силы, <math>N</math> — количество частиц.
 +
 
 +
==Визуализация==
 +
 
 +
 
 +
[[File:vksvt.gif]]
 +
<div class="mw-collapsible mw-collapsed" style="width:100%" >
 +
'''Текст программы на языке Python (среда разработки PyCharm):''' <div class="mw-collapsible-content">
 +
<syntaxhighlight lang="python" line start="1" enclose="div">
 +
import numpy as np
 +
import math
 +
import matplotlib.pyplot as plt
 +
from matplotlib.patches import Circle
 +
from matplotlib.collections import PatchCollection
 +
from celluloid import Camera
 +
 
 +
 
 +
def Plotting(X, rad, u):
 +
    n = np.shape(X)[1]
 +
    patches = []
 +
    for it in range(n):
 +
        circle = Circle((X[0, it], 0), rad)
 +
        patches.append(circle)
 +
    p = PatchCollection(patches)
 +
    ax_1 = plt.subplot(2, 1, 1)
 +
    plt.plot()
 +
    ax_1.add_collection(p)
 +
 
 +
    plt.xlim([-1, N])
 +
    plt.ylim([-4, 4])
 +
    plt.grid(True)
 +
 
 +
 
 +
    ax_2 = plt.subplot(2, 1, 2)
 +
    plt.plot(np.arange(0, N), u.T, '*-k')
 +
    plt.xlim([-3, (n + 2) * a])
 +
    plt.ylim([-1, 1])
 +
    plt.grid(True)
 +
    plt.xlabel("x_i")
 +
    plt.ylabel("Amplitude")
 +
 
 +
    plt.pause(0.01)
 +
    plt.delaxes(ax_1)
 +
    plt.delaxes(ax_2)
 +
 
 +
 
 +
 
 +
#Constants
 +
N = 20
 +
mu = 0.05
 +
c = 1
 +
a = 1
 +
dt = 0.2
 +
m = 1
 +
rad = 0.4
 +
q = -0.1
 +
ST = 1000
  
−
<math>a=1</math>
 
  
−
<math>D=1</math>
+
#Creating arrays
 +
X = np.zeros((1, N))
 +
U = np.zeros((N, ST))
 +
X0 = np.arange(N)
 +
X[0, :] = [it*a for it in range(N)]
  
−
<math>m=1</math>
+
V = np.zeros((1, N))
  
−
<math>dt=0.001</math> -шаг по времени
+
#fig = plt.figure()
 +
#fig, axes = plt.subplots(2)
 +
#camera = Camera(fig)
  
−
[[File:Пуля3.gif|500px|]]
+
for st in range(ST):
 +
    #print(st)
 +
    for sp in range(0, N-1): #loop for springs
  
 +
        F = c*(X[0, sp+1] - X[0, sp] - 1)
 +
        V[0, sp] += F / m * dt - 2*mu*V[0, sp]
 +
        V[0, sp+1] += -F / m * dt - 2*mu*V[0, sp+1]
  
 +
#Periodical Boundary cond
 +
    V[0, N-1] += c * (X[0, 0] + a * N - X[0, N-1] - 1) * dt / m - 2*mu*V[0, N-1] + q*math.cos((st)*2*math.pi/100)
 +
    V[0, 0] += -c * (X[0, 0] + a * N - X[0, N-1] - 1) * dt / m - 2*mu*V[0, 0]
 +
#Positions
 +
    X[0, :] += V[0, :]*dt
  
 +
    Plotting(X, rad, X-X0)
 +
    #camera.snap()
 +
    #anim = camera.animate(blit=True)
 +
    #anim.save('MDS.gif')
 +
</syntaxhighlight>
 +
</div>
 +
</div>
 
=== См. также ===
 
=== См. также ===
  

Текущая версия на 11:16, 24 января 2020

Курсовой проект по Механике дискретных сред

Исполнитель: Шпетный Даниил

Группа: 3630103/60101

Семестр: осень 2019


Постановка задачи[править]

Исследовать вынужденные колебания в вязкой среде. Построить график перемещений частиц от времени.

Построение модели[править]

Процесс моделируется как одномерные колебания цепочки частиц.

Теоретическая сводка[править]

Парное взаимодействие определяется формулой:

[math]F = c(x_{k+1}-x_{k})[/math]

где [math]F[/math] — сила взаимодействия, [math]c[/math] — жесткость связи, [math]x_k[/math] — перемещение частицы, [math]k[/math] — номер частицы.

Уравнение скорости частиц [math]k[/math] и [math]k+1[/math]:

[math]V_{k}=\frac{F}{m}dt-2μV_{k}[/math]
[math]V_{k+1}=-\frac{F}{m}dt-2μV_{k+1}[/math]

где [math]F[/math] — сила взаимодействия, [math]m[/math] — масса частицы, [math]μ[/math] — коэффициент вязкости, [math]V_{k}[/math] — скорость частицы.

Граничные и начальные условия[править]

Начальные условия нулевые

Периодические граничные условия:

[math]V_{N}=c(x_{1}+aN-x_{N}-1)\frac{dt}{m}-2μV_{N}+Qcos(ωt)[/math]
[math]V_{1}=-c(x_{1}+aN-x_{N}-1)\frac{dt}{m}-2μV_{1}[/math]

где [math]Q[/math] — амплитуда возмущающей силы, [math]ω[/math] — частота возмущающей силы, [math]N[/math] — количество частиц.

Визуализация[править]

Vksvt.gif

Текст программы на языке Python (среда разработки PyCharm): 
 1 import numpy as np
 2 import math
 3 import matplotlib.pyplot as plt
 4 from matplotlib.patches import Circle
 5 from matplotlib.collections import PatchCollection
 6 from celluloid import Camera
 7 
 8 
 9 def Plotting(X, rad, u):
10     n = np.shape(X)[1]
11     patches = []
12     for it in range(n):
13         circle = Circle((X[0, it], 0), rad)
14         patches.append(circle)
15     p = PatchCollection(patches)
16     ax_1 = plt.subplot(2, 1, 1)
17     plt.plot()
18     ax_1.add_collection(p)
19 
20     plt.xlim([-1, N])
21     plt.ylim([-4, 4])
22     plt.grid(True)
23 
24 
25     ax_2 = plt.subplot(2, 1, 2)
26     plt.plot(np.arange(0, N), u.T, '*-k')
27     plt.xlim([-3, (n + 2) * a])
28     plt.ylim([-1, 1])
29     plt.grid(True)
30     plt.xlabel("x_i")
31     plt.ylabel("Amplitude")
32 
33     plt.pause(0.01)
34     plt.delaxes(ax_1)
35     plt.delaxes(ax_2)
36 
37 
38 
39 #Constants
40 N = 20
41 mu = 0.05
42 c = 1
43 a = 1
44 dt = 0.2
45 m = 1
46 rad = 0.4
47 q = -0.1
48 ST = 1000
49 
50 
51 #Creating arrays
52 X = np.zeros((1, N))
53 U = np.zeros((N, ST))
54 X0 = np.arange(N)
55 X[0, :] = [it*a for it in range(N)]
56 
57 V = np.zeros((1, N))
58 
59 #fig = plt.figure()
60 #fig, axes = plt.subplots(2)
61 #camera = Camera(fig)
62 
63 for st in range(ST):
64     #print(st)
65     for sp in range(0, N-1): #loop for springs
66 
67         F = c*(X[0, sp+1] - X[0, sp] - 1)
68         V[0, sp] += F / m * dt - 2*mu*V[0, sp]
69         V[0, sp+1] += -F / m * dt - 2*mu*V[0, sp+1]
70 
71 #Periodical Boundary cond
72     V[0, N-1] += c * (X[0, 0] + a * N - X[0, N-1] - 1) * dt / m - 2*mu*V[0, N-1] + q*math.cos((st)*2*math.pi/100)
73     V[0, 0] += -c * (X[0, 0] + a * N - X[0, N-1] - 1) * dt / m - 2*mu*V[0, 0]
74 #Positions
75     X[0, :] += V[0, :]*dt
76 
77     Plotting(X, rad, X-X0)
78     #camera.snap()
79     #anim = camera.animate(blit=True)
80     #anim.save('MDS.gif')

См. также[править]

Источник — «http://tm.spbstu.ru/?title=Вынужденные_колебания_цепочки_в_вязкой_среде&oldid=55074»