ReCo.jl/src/setup.jl

using Distributions: Uniform
using Dates: now

const DEFAULT_PACKING_RATIO = 0.5
const DEFAULT_δt = 1e-5
const DEFAULT_SKIN_TO_INTERACTION_R_RATIO = 1.5
const DEFAULT_EXPORTS_DIR = "exports"
const DEFAULT_PARENT_DIR = ""
const DEFAULT_COMMENT = ""

function initial_particle_grid_pos(
    i::Int64, j::Int64, grid_box_width::Float64, half_box_len::Float64
)
    term = -0.5 * grid_box_width - half_box_len
    return SVector(i * grid_box_width + term, j * grid_box_width + term)
end

function gen_particles(grid_n::Int64, grid_box_width::Float64, half_box_len::Float64)
    particles = Vector{Particle}(undef, grid_n^2)

    id = 1

    for i in 1:grid_n
        for j in 1:grid_n
            particles[id] = Particle(
                id,
                initial_particle_grid_pos(i, j, grid_box_width, half_box_len),
                rand(Uniform(-π, π)),
            )

            id += 1
        end
    end

    return particles
end

function gen_particles(bundle::Bundle, n_particles::Int64)
    particles = Vector{Particle}(undef, n_particles)

    @simd for id in 1:n_particles
        particles[id] = Particle(id, bundle.c[id, end], bundle.φ[id, end])
    end

    return particles
end

function gen_sim_consts(
    n_particles::Int64,
    v₀::Float64;
    δt::Float64=DEFAULT_δt,
    packing_ratio::Float64=DEFAULT_PACKING_RATIO,
    skin_to_interaction_r_ratio::Float64=DEFAULT_SKIN_TO_INTERACTION_R_RATIO,
)
    @assert n_particles > 0
    @assert v₀ >= 0
    @assert δt in 1e-7:1e-7:1e-4
    @assert packing_ratio > 0

    μ = 1.0
    D₀ = 1.0
    particle_radius = 0.5
    Dᵣ = 3 * D₀ / ((2 * particle_radius)^2)
    σ = 1.0
    ϵ = 100.0
    interaction_r = 2^(1 / 6) * σ

    if v₀ != 0.0
        buffer = 1.8
        max_approach_after_one_integration_step = buffer * (2 * v₀ * δt) / interaction_r
        @assert skin_to_interaction_r_ratio >= 1 + max_approach_after_one_integration_step

        skin_r = skin_to_interaction_r_ratio * interaction_r
        n_steps_before_verlet_list_update = round(
            Int64,
            (skin_to_interaction_r_ratio - 1) / max_approach_after_one_integration_step,
        )
    else
        skin_r = 1.5 * interaction_r
        n_steps_before_verlet_list_update = 100
    end

    grid_n = round(Int64, ceil(sqrt(n_particles)))

    n_particles = grid_n^2

    half_box_len = sqrt(n_particles * π / packing_ratio) * σ / 4
    grid_box_width = 2 * half_box_len / grid_n

    return (;
        # Input
        n_particles,
        v₀,
        δt,
        packing_ratio,
        # Calculated
        μ,
        D₀,
        particle_radius,
        Dᵣ,
        σ,
        ϵ,
        interaction_r,
        skin_r,
        n_steps_before_verlet_list_update,
        grid_n,
        half_box_len,
        grid_box_width,
    )
end

function init_sim_with_sim_consts(
    sim_consts;
    exports_dir::String=DEFAULT_EXPORTS_DIR,
    parent_dir::String=DEFAULT_PARENT_DIR,
    comment::String=DEFAULT_COMMENT,
)
    particles = gen_particles(
        sim_consts.grid_n, sim_consts.grid_box_width, sim_consts.half_box_len
    )
    bundle = Bundle(sim_consts.n_particles, 1)
    save_snapshot!(bundle, 1, 0.0, particles)

    dir = exports_dir

    if length(parent_dir) > 0
        dir *= "/$parent_dir"
    end

    start_datetime = now()
    dir *= "/$(start_datetime)_N=$(sim_consts.n_particles)_v=$(sim_consts.v₀)_#$(rand(1000:9999))"

    if length(comment) > 0
        dir *= "_$comment"
    end

    mkpath(dir)

    task = @async write_struct_to_json(sim_consts, "$dir/sim_consts")

    save_bundle(dir, bundle, 1, 0.0)

    runs_dir = "$dir/runs"
    mkpath(runs_dir)

    wait(task)

    return (dir, particles)
end

function init_sim(;
    n_particles::Int64,
    v₀::Float64,
    δt::Float64=DEFAULT_δt,
    packing_ratio::Float64=DEFAULT_PACKING_RATIO,
    skin_to_interaction_r_ratio::Float64=DEFAULT_SKIN_TO_INTERACTION_R_RATIO,
    exports_dir::String=DEFAULT_EXPORTS_DIR,
    parent_dir::String=DEFAULT_PARENT_DIR,
    comment::String=DEFAULT_COMMENT,
)
    sim_consts = gen_sim_consts(
        n_particles,
        v₀;
        δt=δt,
        packing_ratio=packing_ratio,
        skin_to_interaction_r_ratio=skin_to_interaction_r_ratio,
    )

    return init_sim_with_sim_consts(
        sim_consts; exports_dir=exports_dir, parent_dir=parent_dir, comment=comment
    )[1]
end