cheetah conditions wip
2211567
parent
ece09a29f5
commit
f5d3ebd0c2
|
|
@ -5,7 +5,7 @@ include("../src/visualization.jl")
|
|||
using .Visualization
|
||||
|
||||
N = 256
|
||||
tspan = (0.0, 1500.0)
|
||||
tspan = (1.0, 3000.0)
|
||||
|
||||
sol = AnimalFurFHN.run_simulation(tspan, N)
|
||||
|
||||
|
|
|
|||
|
|
@ -3,6 +3,7 @@ module AnimalFurFHN
|
|||
|
||||
include("constants.jl")
|
||||
include("laplacian.jl")
|
||||
include("utils.jl")
|
||||
include("solver.jl")
|
||||
|
||||
export run_simulation # Make sure this is here!
|
||||
|
|
|
|||
174
src/solver.jl
174
src/solver.jl
|
|
@ -30,162 +30,6 @@ function fhn!(du, u, p::FHNParams, t=0)
|
|||
du .= vcat(vec(fu), vec(fv))
|
||||
end
|
||||
|
||||
# helper functions for filling cells in specific places of the matrix
|
||||
function blocks_ic(N)
|
||||
u = fill(1.0, N, N)
|
||||
v = fill(0.0, N, N)
|
||||
p = div(N, 2)
|
||||
|
||||
function safe_block!(u, row_center, col_center)
|
||||
row_start = max(row_center - 8, 1)
|
||||
row_end = min(row_center + 7, N)
|
||||
col_start = max(col_center - 8, 1)
|
||||
col_end = min(col_center + 7, N)
|
||||
u[row_start:row_end, col_start:col_end] .= -0.01
|
||||
end
|
||||
|
||||
safe_block!(u, p, p)
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function column_ic(N)
|
||||
u = fill(0.01, N, N)
|
||||
v = fill(0.99, N, N)
|
||||
|
||||
col_center = div(N, 2)
|
||||
col_width = 8 # You can adjust this
|
||||
|
||||
col_start = max(col_center - div(col_width, 2), 1)
|
||||
col_end = min(col_center + div(col_width, 2) - 1, N)
|
||||
|
||||
u[col_start:col_end, :] .= -0.01
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function two_rows_edge_distance_ic(N)
|
||||
row_width = 8
|
||||
distance_from_edge = 50
|
||||
u = fill(0.01, N, N)
|
||||
v = fill(0.99, N, N)
|
||||
|
||||
# --- Input Validation ---
|
||||
if row_width <= 0 || distance_from_edge < 0
|
||||
error("row_width must be positive and distance_from_edge must be non-negative.")
|
||||
end
|
||||
|
||||
# Calculate column 1 (from the left edge)
|
||||
col1_start = distance_from_edge + 1
|
||||
col1_end = col1_start + row_width - 1
|
||||
|
||||
# Calculate column 2 (from the right edge)
|
||||
col2_end = N - distance_from_edge
|
||||
col2_start = col2_end - row_width + 1
|
||||
|
||||
# --- Further Validation for placement ---
|
||||
if col1_end > N || col2_start < 1
|
||||
error("Columns go out of bounds. Adjust N, row_width, or distance_from_edge.")
|
||||
end
|
||||
if col1_end >= col2_start # Check for overlap or touching
|
||||
error("Columns overlap or touch. Adjust N, row_width, or distance_from_edge such that 2 * (distance_from_edge + row_width) <= N.")
|
||||
end
|
||||
|
||||
# Apply the first column
|
||||
u[:, col1_start:col1_end] .= -0.01
|
||||
|
||||
# Apply the second column
|
||||
u[:, col2_start:col2_end] .= -0.01
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function center_band_ic(N)
|
||||
u = fill(0.0, N, N)
|
||||
v = fill(0.0, N, N)
|
||||
|
||||
band_width = div(N, 8)
|
||||
row_start = div(N, 2) - div(band_width, 2)
|
||||
row_end = div(N, 2) + div(band_width, 2)
|
||||
|
||||
u[row_start:row_end, :] .= 0.1 .+ 0.01 .* randn(band_width + 1, N)
|
||||
v[row_start:row_end, :] .= 0.1 .+ 0.01 .* randn(band_width + 1, N)
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function circle_ic(N)
|
||||
u = fill(0.01, N, N)
|
||||
v = fill(0.99, N, N)
|
||||
cx, cy = div(N, 2), div(N, 2) # center of matrix
|
||||
radius = 0.125 * N # circle radius = 3/4 of N divided by 2
|
||||
|
||||
for i in 1:N, j in 1:N
|
||||
if (i - cx)^2 + (j - cy)^2 ≤ radius^2
|
||||
u[i, j] = -0.01
|
||||
end
|
||||
end
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function three_circles_random_ic(N)
|
||||
u = fill(0.01, N, N)
|
||||
v = fill(0.99, N, N)
|
||||
radius = 0.125 * N
|
||||
|
||||
# Define the bounds for random centers to ensure the circle stays within the matrix
|
||||
min_coord = ceil(Int, radius) + 1
|
||||
max_coord = floor(Int, N - radius)
|
||||
|
||||
if min_coord > max_coord
|
||||
error("Matrix size N is too small to place circles of this radius without overlap or going out of bounds.")
|
||||
end
|
||||
|
||||
for _ in 1:5 # Place 3 circles
|
||||
# Generate random center coordinates
|
||||
cx = rand(min_coord:max_coord)
|
||||
cy = rand(min_coord:max_coord)
|
||||
|
||||
# Apply the circle to the matrix
|
||||
for i in 1:N, j in 1:N
|
||||
if (i - cx)^2 + (j - cy)^2 ≤ radius^2
|
||||
u[i, j] = -0.01
|
||||
end
|
||||
end
|
||||
end
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function squiggle_ic(N, Lx=400.0, Ly=400.0)
|
||||
uplus = 0.01
|
||||
vplus = 0.99
|
||||
uminus = -uplus
|
||||
|
||||
# Create coordinate grids
|
||||
x = LinRange(0, Lx, N)
|
||||
y = LinRange(0, Ly, N)
|
||||
X = repeat(x', N, 1) # Transposed to align with meshgrid X
|
||||
Y = repeat(y, 1, N) # Broadcasted to align with meshgrid Y
|
||||
|
||||
# Squiggle pattern
|
||||
cos_term = 0.05 * Lx .* sin.(10 * 2π .* Y ./ Ly .+ π * 0.3)
|
||||
exp_term = exp.(-((Y .- Ly / 2) ./ (0.1 * Ly)) .^ 2)
|
||||
width = 0.05 * Lx
|
||||
Z = exp.(-((X .- Lx / 2 .+ cos_term .* exp_term) ./ width) .^ 2)
|
||||
|
||||
|
||||
u = fill(uplus, N, N)
|
||||
v = fill(vplus, N, N)
|
||||
|
||||
# Apply squiggle
|
||||
u[Z .> 0.8] .= uminus
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
|
||||
"""
|
||||
run_simulation(tspan::Tuple{Float64,Float64}, N::Int)
|
||||
|
||||
|
|
@ -199,21 +43,11 @@ end
|
|||
- `sol`: solved differential equation (ODE)
|
||||
"""
|
||||
function run_simulation(tspan::Tuple{Float64,Float64}, N::Int)
|
||||
# Turing-spot parameters
|
||||
p = FHNParams(N=N, dx=1.0, Du=0.016, Dv=0.1, ϵ=0.1, a=0.5, b=0.9)
|
||||
# Initial conditions
|
||||
# params, y0 = zebra_conditions(N) # tspan of ~3500 enough
|
||||
params, y0 = cheetah_conditions(N)
|
||||
|
||||
# Initial conditions (random noise)
|
||||
#Random.seed!(4321)
|
||||
# Create two vectors with length N*N with numbers between 0.1 and 0.11
|
||||
#u0 = vec(0.4 .+ 0.01 .* rand(N, N))
|
||||
#v0 = vec(0.4 .+ 0.01 .* rand(N, N))
|
||||
|
||||
# Or use this
|
||||
u0, v0 = two_rows_edge_distance_ic(N)
|
||||
|
||||
y0 = vcat(vcat(u0, v0))
|
||||
|
||||
prob = ODEProblem(fhn!, y0, tspan, p)
|
||||
prob = ODEProblem(fhn!, y0, tspan, params)
|
||||
sol = solve(prob, BS3(), saveat=10.0) # You can try `Rosenbrock23()` too
|
||||
|
||||
return sol
|
||||
|
|
|
|||
|
|
@ -0,0 +1,184 @@
|
|||
|
||||
# initial conditions for different patterns
|
||||
function zebra_conditions(N)
|
||||
# Turing-spot parameters
|
||||
params = FHNParams(N=N, dx=1.0, Du=0.016, Dv=0.1, ϵ=0.1, a=0.5, b=0.9)
|
||||
|
||||
# Or use this
|
||||
u0, v0 = two_rows_edge_distance_ic(N)
|
||||
|
||||
return params, vcat(u0, v0)
|
||||
end
|
||||
|
||||
function cheetah_conditions(N)
|
||||
# Turing-spot parameters
|
||||
params = FHNParams(N=N, dx=1.0, Du=0.0025, Dv=0.6, ϵ=0.05, a=0.7, b=0.8)
|
||||
|
||||
# Approximate Homogenous Steady State (HSS) for a=0.7, b=0.8
|
||||
u_hss = -1.2
|
||||
v_hss = -0.625
|
||||
|
||||
# Small perturbations around the HSS
|
||||
# rand(N,N) generates values between 0 and 1.
|
||||
# (2 .* rand(N,N) .- 1) generates values between -1 and 1 symmetrically around 0.
|
||||
perturbation_amplitude = 0.01
|
||||
u0 = vec(u_hss .+ perturbation_amplitude .* (2 .* rand(N, N) .- 1))
|
||||
v0 = vec(v_hss .+ perturbation_amplitude .* (2 .* rand(N, N) .- 1))
|
||||
|
||||
return params, vcat(u0, v0)
|
||||
end
|
||||
|
||||
# helper functions for filling cells in specific places of the matrix
|
||||
function blocks_ic(N)
|
||||
u = fill(1.0, N, N)
|
||||
v = fill(0.0, N, N)
|
||||
p = div(N, 2)
|
||||
|
||||
function safe_block!(u, row_center, col_center)
|
||||
row_start = max(row_center - 8, 1)
|
||||
row_end = min(row_center + 7, N)
|
||||
col_start = max(col_center - 8, 1)
|
||||
col_end = min(col_center + 7, N)
|
||||
u[row_start:row_end, col_start:col_end] .= -0.01
|
||||
end
|
||||
|
||||
safe_block!(u, p, p)
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function column_ic(N)
|
||||
u = fill(0.01, N, N)
|
||||
v = fill(0.99, N, N)
|
||||
|
||||
col_center = div(N, 2)
|
||||
col_width = 8 # You can adjust this
|
||||
|
||||
col_start = max(col_center - div(col_width, 2), 1)
|
||||
col_end = min(col_center + div(col_width, 2) - 1, N)
|
||||
|
||||
u[col_start:col_end, :] .= -0.01
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function two_rows_edge_distance_ic(N)
|
||||
row_width = 8
|
||||
distance_from_edge = 50
|
||||
u = fill(0.01, N, N)
|
||||
v = fill(0.99, N, N)
|
||||
|
||||
# --- Input Validation ---
|
||||
if row_width <= 0 || distance_from_edge < 0
|
||||
error("row_width must be positive and distance_from_edge must be non-negative.")
|
||||
end
|
||||
|
||||
# Calculate column 1 (from the left edge)
|
||||
col1_start = distance_from_edge + 1
|
||||
col1_end = col1_start + row_width - 1
|
||||
|
||||
# Calculate column 2 (from the right edge)
|
||||
col2_end = N - distance_from_edge
|
||||
col2_start = col2_end - row_width + 1
|
||||
|
||||
# --- Further Validation for placement ---
|
||||
if col1_end > N || col2_start < 1
|
||||
error("Columns go out of bounds. Adjust N, row_width, or distance_from_edge.")
|
||||
end
|
||||
if col1_end >= col2_start # Check for overlap or touching
|
||||
error("Columns overlap or touch. Adjust N, row_width, or distance_from_edge such that 2 * (distance_from_edge + row_width) <= N.")
|
||||
end
|
||||
|
||||
# Apply the first column
|
||||
u[:, col1_start:col1_end] .= -0.01
|
||||
|
||||
# Apply the second column
|
||||
u[:, col2_start:col2_end] .= -0.01
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function center_band_ic(N)
|
||||
u = fill(0.0, N, N)
|
||||
v = fill(0.0, N, N)
|
||||
|
||||
band_width = div(N, 8)
|
||||
row_start = div(N, 2) - div(band_width, 2)
|
||||
row_end = div(N, 2) + div(band_width, 2)
|
||||
|
||||
u[row_start:row_end, :] .= 0.1 .+ 0.01 .* randn(band_width + 1, N)
|
||||
v[row_start:row_end, :] .= 0.1 .+ 0.01 .* randn(band_width + 1, N)
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function circle_ic(N)
|
||||
u = fill(0.01, N, N)
|
||||
v = fill(0.99, N, N)
|
||||
cx, cy = div(N, 2), div(N, 2) # center of matrix
|
||||
radius = 0.125 * N # circle radius = 3/4 of N divided by 2
|
||||
|
||||
for i in 1:N, j in 1:N
|
||||
if (i - cx)^2 + (j - cy)^2 ≤ radius^2
|
||||
u[i, j] = -0.01
|
||||
end
|
||||
end
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function three_circles_random_ic(N)
|
||||
u = fill(0.01, N, N)
|
||||
v = fill(0.99, N, N)
|
||||
radius = 0.125 * N
|
||||
|
||||
# Define the bounds for random centers to ensure the circle stays within the matrix
|
||||
min_coord = ceil(Int, radius) + 1
|
||||
max_coord = floor(Int, N - radius)
|
||||
|
||||
if min_coord > max_coord
|
||||
error("Matrix size N is too small to place circles of this radius without overlap or going out of bounds.")
|
||||
end
|
||||
|
||||
for _ in 1:5 # Place 3 circles
|
||||
# Generate random center coordinates
|
||||
cx = rand(min_coord:max_coord)
|
||||
cy = rand(min_coord:max_coord)
|
||||
|
||||
# Apply the circle to the matrix
|
||||
for i in 1:N, j in 1:N
|
||||
if (i - cx)^2 + (j - cy)^2 ≤ radius^2
|
||||
u[i, j] = -0.01
|
||||
end
|
||||
end
|
||||
end
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
||||
function squiggle_ic(N, Lx=400.0, Ly=400.0)
|
||||
uplus = 0.01
|
||||
vplus = 0.99
|
||||
uminus = -uplus
|
||||
|
||||
# Create coordinate grids
|
||||
x = LinRange(0, Lx, N)
|
||||
y = LinRange(0, Ly, N)
|
||||
X = repeat(x', N, 1) # Transposed to align with meshgrid X
|
||||
Y = repeat(y, 1, N) # Broadcasted to align with meshgrid Y
|
||||
|
||||
# Squiggle pattern
|
||||
cos_term = 0.05 * Lx .* sin.(10 * 2π .* Y ./ Ly .+ π * 0.3)
|
||||
exp_term = exp.(-((Y .- Ly / 2) ./ (0.1 * Ly)) .^ 2)
|
||||
width = 0.05 * Lx
|
||||
Z = exp.(-((X .- Lx / 2 .+ cos_term .* exp_term) ./ width) .^ 2)
|
||||
|
||||
|
||||
u = fill(uplus, N, N)
|
||||
v = fill(vplus, N, N)
|
||||
|
||||
# Apply squiggle
|
||||
u[Z .> 0.8] .= uminus
|
||||
|
||||
return vec(u), vec(v)
|
||||
end
|
||||
|
|
@ -21,7 +21,7 @@ function step_through_solution(sol, N::Int)
|
|||
# Initialize heatmap with first time step
|
||||
u0 = reshape(sol[1][1:N^2], N, N)
|
||||
heat_obs = Observable(u0)
|
||||
hmap = heatmap!(ax, heat_obs, colormap=:berlin)
|
||||
hmap = heatmap!(ax, heat_obs, colormap=:RdGy)
|
||||
|
||||
# Update heatmap on slider movement
|
||||
on(slider.value) do i
|
||||
|
|
|
|||
Loading…
Reference in New Issue