https://hrubos.tech/blogy/content/images/20260902103835-vlak_raz_odraz.gif
using CairoMakie
using Printf
# ============================================================
# VLAK — RÁZOVÁ VLNA TAM A SPÄŤ
#
# 1 rušeň + 10 vozňov
#
# Model:
#
# RUŠEŇ | V1 | V2 | V3 | ... | V10
# →
# →
# →
# ←
# ←
#
# Spriahadlá majú mechanickú vôľu.
#
# Výsledkom je:
#
# vlak_raz_odraz.gif
#
# ============================================================
function main()
println()
println("==============================================")
println(" VLAK — RÁZOVÁ VLNA TAM A SPÄŤ")
println("==============================================")
println()
# ========================================================
# 1. VLAK
# ========================================================
N = 11
# 1 = rušeň
# 2..11 = vozne
loco_mass = 120_000.0
wagon_mass = 80_000.0
masses = fill(wagon_mass, N)
masses[1] = loco_mass
# ========================================================
# 2. GEOMETRIA
# ========================================================
wagon_length = 20.0
gap = 0.5
spacing =
wagon_length + gap
positions =
[-(i - 1) * spacing for i in 1:N]
velocities =
zeros(Float64, N)
# ========================================================
# 3. SPRIAHADLÁ
# ========================================================
# Tuhosť
K = 2.0e7
# Tlmenie
C = 3.0e4
# Mechanická vôľa
SLACK = 0.15
# ========================================================
# 4. NÁRAZ RUŠŇA
# ========================================================
IMPULSE_FORCE = 2.5e6
IMPULSE_TIME = 0.12
# ========================================================
# 5. ČAS
# ========================================================
DT = 0.0001
TOTAL_TIME = 8.0
STEPS =
Int(round(TOTAL_TIME / DT))
SAVE_EVERY = 10
# ========================================================
# 6. HISTÓRIA
# ========================================================
estimated_states =
fld(STEPS, SAVE_EVERY) + 2
times =
Vector{Float64}(
undef,
estimated_states
)
position_history =
Vector{Vector{Float64}}(
undef,
estimated_states
)
force_history =
Vector{Vector{Float64}}(
undef,
estimated_states
)
# ========================================================
# 7. SILA SPRIAHADLA
# ========================================================
function coupler_force(
x_left,
x_right,
v_left,
v_right
)
# Aktuálna vzdialenosť
distance =
x_left - x_right
# Rovnovážna vzdialenosť
equilibrium =
spacing
# Odchýlka
delta =
distance - equilibrium
# Relatívna rýchlosť
relative_velocity =
v_left - v_right
# ----------------------------------------------------
# VÔĽA
# ----------------------------------------------------
if delta > SLACK
effective_delta =
delta - SLACK
return (
K * effective_delta +
C * relative_velocity
)
elseif delta < -SLACK
effective_delta =
delta + SLACK
return (
K * effective_delta +
C * relative_velocity
)
else
return 0.0
end
end
# ========================================================
# 8. ULOŽENIE STAVU
# ========================================================
save_index = 1
times[save_index] = 0.0
position_history[save_index] =
copy(positions)
force_history[save_index] =
zeros(Float64, N - 1)
save_index += 1
# ========================================================
# 9. FYZIKÁLNA SIMULÁCIA
# ========================================================
for step in 1:STEPS
t =
step * DT
# ----------------------------------------------------
# Sily na vozidlách
# ----------------------------------------------------
forces =
zeros(Float64, N)
# ----------------------------------------------------
# Sily v spriahadlách
# ----------------------------------------------------
coupler_forces =
zeros(Float64, N - 1)
# ----------------------------------------------------
# Vypočítame všetkých 10 spriahadiel
# ----------------------------------------------------
for i in 1:(N - 1)
F =
coupler_force(
positions[i],
positions[i + 1],
velocities[i],
velocities[i + 1]
)
coupler_forces[i] =
F
# Newton III.
forces[i] -= F
forces[i + 1] += F
end
# ----------------------------------------------------
# Nárazová sila do rušňa
# ----------------------------------------------------
if t <= IMPULSE_TIME
forces[1] +=
IMPULSE_FORCE
end
# ----------------------------------------------------
# Zrýchlenie
# ----------------------------------------------------
acceleration =
forces ./ masses
# ----------------------------------------------------
# Semi-implicit Euler
# ----------------------------------------------------
velocities .+=
acceleration .* DT
positions .+=
velocities .* DT
# ----------------------------------------------------
# Uloženie
# ----------------------------------------------------
if step % SAVE_EVERY == 0
if save_index <= length(times)
times[save_index] =
t
position_history[save_index] =
copy(positions)
force_history[save_index] =
copy(coupler_forces)
save_index += 1
end
end
# ----------------------------------------------------
# Progress
# ----------------------------------------------------
if step % 5000 == 0
percent =
100.0 *
step /
STEPS
@printf(
"Simulácia: %6.1f %%\n",
percent
)
end
end
# ========================================================
# 10. OREZ HISTÓRIE
# ========================================================
actual_count =
save_index - 1
times =
times[1:actual_count]
position_history =
position_history[1:actual_count]
force_history =
force_history[1:actual_count]
println()
println(
"Simulácia dokončená: ",
length(times),
" stavov."
)
println()
# ========================================================
# 11. OBSERVABLES
# ========================================================
position_obs =
Observable(
copy(position_history[1])
)
force_obs =
Observable(
copy(force_history[1]) ./ 1.0e6
)
time_text =
Observable(
"t = 0.000 s"
)
wave_text =
Observable(
"Vlna: →"
)
peak_text =
Observable(
"Maximum: —"
)
# ========================================================
# 12. FIGURE
# ========================================================
fig =
Figure(
size = (1500, 900)
)
# ========================================================
# 13. HORNÝ PANEL — VLAK
# ========================================================
ax_train =
Axis(
fig[1, 1],
title = "RÁZOVÁ VLNA VO VLAKU",
xlabel = "poloha [m]",
ylabel = "",
limits = (
-230,
30,
0,
10
)
)
# ========================================================
# 14. STREDNÝ PANEL — SILA
# ========================================================
ax_force =
Axis(
fig[2, 1],
title = "AKTUÁLNA SILA V SPRIAHADLÁCH",
xlabel = "číslo spriahadla",
ylabel = "sila [MN]",
limits = (
0.5,
10.5,
-10,
10
)
)
# ========================================================
# 15. RUŠEŇ
# ========================================================
loco_x =
Observable(
[positions[1]]
)
scatter!(
ax_train,
loco_x,
[5.0],
markersize = 45
)
# ========================================================
# 16. VOZNE
# ========================================================
wagon_x =
Observable(
positions[2:end]
)
scatter!(
ax_train,
wagon_x,
fill(5.0, N - 1),
markersize = 35
)
# ========================================================
# 17. ČÍSLA VOZŇOV
# ========================================================
label_x =
Observable(
copy(positions)
)
label_text =
[
i == 1 ?
"RUŠEŇ" :
"V$i"
for i in 1:N
]
text!(
ax_train,
label_x,
fill(7.0, N),
text = label_text,
align = (:center, :center),
fontsize = 15
)
# ========================================================
# 18. SPRIAHADLÁ
# ========================================================
coupler_x =
Vector{Observable{Vector{Float64}}}(
undef,
N - 1
)
for i in 1:(N - 1)
coupler_x[i] =
Observable(
[
positions[i],
positions[i + 1]
]
)
lines!(
ax_train,
coupler_x[i],
[5.0, 5.0],
linewidth = 5
)
end
# ========================================================
# 19. TEXT ČASU
# ========================================================
text!(
ax_train,
-225.0,
9.3,
text = time_text,
fontsize = 25,
align = (:left, :center)
)
# ========================================================
# 20. TEXT SMERU VLNY
# ========================================================
text!(
ax_train,
-225.0,
8.5,
text = wave_text,
fontsize = 25,
align = (:left, :center)
)
# ========================================================
# 21. TEXT MAXIMA
# ========================================================
text!(
ax_train,
-225.0,
7.7,
text = peak_text,
fontsize = 18,
align = (:left, :center)
)
# ========================================================
# 22. GRAF SILY
# ========================================================
coupler_numbers =
collect(1:(N - 1))
lines!(
ax_force,
coupler_numbers,
force_obs,
linewidth = 4
)
scatter!(
ax_force,
coupler_numbers,
force_obs,
markersize = 12
)
# ========================================================
# 23. FRAME INDEXY
# ========================================================
FPS = 30
total_states =
length(times)
frame_count =
min(
total_states,
max(
1,
Int(round(TOTAL_TIME * FPS))
)
)
# Dôležité:
#
# výsledkom je Vector{Int}
#
frame_indices =
collect(
round.(
Int,
range(
1,
total_states,
length = frame_count
)
)
)
println(
"Frameov GIF: ",
length(frame_indices)
)
# ========================================================
# 24. SMER VLNY
# ========================================================
previous_peak =
1
direction =
"→"
# ========================================================
# 25. VÝSTUP
# ========================================================
OUTPUT_FILE =
"vlak_raz_odraz.gif"
println()
println(
"Vytváram: ",
OUTPUT_FILE
)
println()
# ========================================================
# 26. RECORD
# ========================================================
record(
fig,
OUTPUT_FILE,
eachindex(frame_indices);
framerate = FPS
) do frame_number
# ----------------------------------------------------
# INDEX JE VŽDY INT
# ----------------------------------------------------
idx =
frame_indices[frame_number]
# ----------------------------------------------------
# Aktuálny stav
# ----------------------------------------------------
current_positions =
position_history[idx]
current_forces =
force_history[idx]
# ----------------------------------------------------
# Vozidlá
# ----------------------------------------------------
loco_x[] =
[current_positions[1]]
wagon_x[] =
copy(current_positions[2:end])
label_x[] =
copy(current_positions)
# ----------------------------------------------------
# Spriahadlá
# ----------------------------------------------------
for i in 1:(N - 1)
coupler_x[i][] =
[
current_positions[i],
current_positions[i + 1]
]
end
# ----------------------------------------------------
# Sily
# ----------------------------------------------------
force_obs[] =
current_forces ./ 1.0e6
# ====================================================
# NÁJDEME NAJSILNEJŠIE SPRIAHADLO
# ====================================================
absolute_forces =
abs.(current_forces)
peak_index =
argmax(absolute_forces)
peak_force =
current_forces[peak_index]
# ====================================================
# URČENIE SMERU
# ====================================================
if peak_index > previous_peak
direction =
"→"
elseif peak_index < previous_peak
direction =
"←"
end
previous_peak =
peak_index
# ====================================================
# TEXT
# ====================================================
time_text[] =
@sprintf(
"t = %.3f s",
times[idx]
)
wave_text[] =
"Vlna: " *
direction
peak_text[] =
@sprintf(
"Maximum: spriahadlo %d %.2f MN",
peak_index,
peak_force / 1.0e6
)
# ====================================================
# PROGRESS
# ====================================================
if frame_number == 1 ||
frame_number % 10 == 0 ||
frame_number == length(frame_indices)
percent =
100.0 *
frame_number /
length(frame_indices)
@printf(
"GIF: %6.1f %%\n",
percent
)
end
end
# ========================================================
# HOTOVO
# ========================================================
println()
println("==============================================")
println(" HOTOVO")
println("==============================================")
println()
println(
"Výsledok: ",
OUTPUT_FILE
)
println()
end
# ============================================================
# SPUSTENIE
# ============================================================
main()
Comments “Ako sa PRAVDEPODOBNE roztrhol nákladný vlak:”