1
0
Fork 0
mirror of https://gitlab.rlp.net/mobitar/ReCo.jl.git synced 2024-11-08 22:21:08 +00:00
ReCo.jl/graphics/radial_distribution.jl

87 lines
1.9 KiB
Julia
Raw Normal View History

2022-01-25 04:12:17 +00:00
using Luxor
using Random: Random
using StaticArrays: SVector
using ReCo: ReCo
function gen_rdf_graphics()
2022-01-25 04:23:26 +00:00
Random.seed!(1)
2022-01-25 04:12:17 +00:00
2022-01-26 22:37:51 +00:00
box_length = 100
2022-01-25 04:12:17 +00:00
graphics_export_dir = "exports/graphics"
mkpath(graphics_export_dir)
2022-01-26 22:37:51 +00:00
Drawing(box_length, box_length, "$graphics_export_dir/rdf_shells.pdf")
2022-01-25 04:12:17 +00:00
origin()
2022-01-25 04:23:26 +00:00
particle_radius = 6
Δr = 2.4 * particle_radius
2022-01-25 04:12:17 +00:00
setcolor("blue")
circle(Point(0, 0), particle_radius, :fill)
N = 50
2022-01-25 04:23:26 +00:00
selected_shell_ind = 3
selected_lower_radius = (selected_shell_ind - 1) * Δr
2022-01-25 04:12:17 +00:00
particle_cs = Vector{SVector{2,Float64}}(undef, N)
2022-01-27 02:17:51 +00:00
distance_limit = 0.99 * 2 * particle_radius
2022-01-25 04:12:17 +00:00
x = y = distance = 0.0
for p1_ind in 1:N
while true
2022-01-26 22:37:51 +00:00
x = (rand() - 0.5) * box_length
y = (rand() - 0.5) * box_length
2022-01-25 04:12:17 +00:00
p1_c = SVector(x, y)
distance = ReCo.norm2d(p1_c)
no_collision = true
2022-01-27 02:17:51 +00:00
if distance >= distance_limit
2022-01-25 04:12:17 +00:00
for p2_ind in 1:(p1_ind - 1)
2022-01-27 02:17:51 +00:00
if ReCo.norm2d(p1_c - particle_cs[p2_ind]) < distance_limit
2022-01-25 04:12:17 +00:00
no_collision = false
break
end
end
if no_collision
particle_cs[p1_ind] = p1_c
break
end
end
end
if selected_lower_radius <= distance < selected_lower_radius + Δr
setcolor("green")
else
setcolor("orange")
end
circle(x, y, particle_radius, :fill)
end
setcolor("black")
setline(0.2)
for shell_ind in 1:ceil(Int64, (box_length - 1) / 2 / Δr)
circle(Point(0, 0), shell_ind * Δr, :stroke)
end
setcolor("red")
setline(0.3)
line(Point(0, 0), Point(Δr, 0), :stroke)
2022-01-25 04:23:26 +00:00
fontsize(4.5)
text("Δr", Point(0.44 * Δr, -0.05 * Δr))
2022-01-25 04:12:17 +00:00
finish()
return nothing
end