Make daisyworld runnable

Author: Datseris

Project.toml

@@ -1,7 +1,7 @@
 name = "InteractiveDynamics"
 uuid = "ec714cd0-5f51-11eb-0b6e-452e7367ff84"
 repo = "https://github.com/JuliaDynamics/InteractiveDynamics.jl.git"
-version = "0.21.5"
+version = "0.21.6"
 
 [deps]
 DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"

src/agents/abmplot.jl

@@ -20,7 +20,7 @@ Requires `Agents`. See also [`abmvideo`](@ref) and [`abmexploration`](@ref).
   each agent will be plotted as. They can each be either a constant or a *function*,
   which takes as an input a single agent and outputs the corresponding value.
 
-  Using constants: `ac = "#338c54", as = 10, am = :diamond`
+  Using constants: `ac = "#338c54", as = 15, am = :diamond`
 
   Using functions:
   ```julia
@@ -156,7 +156,7 @@ This is the internal recipe for creating an `_ABMPlot`.
 
         # Agent
         ac = JULIADYNAMICS_COLORS[1],
-        as = 10,
+        as = 15,
         am = :circle,
         offset = nothing,
         scatterkwargs = NamedTuple(),

src/agents/daisyworld_def.jl

@@ -1,96 +1,83 @@
 using Agents
-using Statistics: mean
 import StatsBase
 using Random
 
-
 @agent Daisy GridAgent{2} begin
     breed::Symbol
     age::Int
     albedo::Float64 # 0-1 fraction
 end
 
-DaisyWorld = ABM{<:GridSpaceSingle, Daisy};
+DaisyWorld = ABM{<:GridSpaceSingle, Daisy}
 
 function update_surface_temperature!(pos, model::DaisyWorld)
-    ids = ids_in_position(pos, model)
-    absorbed_luminosity = if isempty(ids) # no daisy
-        ## Set luminosity via surface albedo
+    absorbed_luminosity = if isempty(pos, model) # no daisy
         (1 - model.surface_albedo) * model.solar_luminosity
     else
-        ## Set luminosity via daisy albedo
-        (1 - model[ids[1]].albedo) * model.solar_luminosity
+        daisy = model[id_in_position(pos, model)]
+        (1 - daisy.albedo) * model.solar_luminosity
     end
-    ## We expect local heating to be 80 ᵒC for an absorbed luminosity of 1,
-    ## approximately 30 for 0.5 and approximately -273 for 0.01.
     local_heating = absorbed_luminosity > 0 ? 72 * log(absorbed_luminosity) + 80 : 80
-    ## Surface temperature is the average of the current temperature and local heating.
     model.temperature[pos...] = (model.temperature[pos...] + local_heating) / 2
 end
 
 function diffuse_temperature!(pos, model::DaisyWorld)
-    ratio = get(model.properties, :ratio, 0.5) # diffusion ratio
+    ratio = model.ratio # diffusion ratio
     npos = nearby_positions(pos, model)
     model.temperature[pos...] =
         (1 - ratio) * model.temperature[pos...] +
-        ## Each neighbor is giving up 1/8 of the diffused
-        ## amount to each of *its* neighbors
         sum(model.temperature[p...] for p in npos) * 0.125 * ratio
 end
 
-function propagate_daisy!(pos, model::DaisyWorld)
-    ids = ids_in_position(pos, model)
-    if !isempty(ids)
-        daisy = model[ids[1]]
-        temperature = model.temperature[pos...]
-        ## Set optimum growth rate to 22.5 ᵒC, with bounds of [5, 40]
-        seed_threshold = (0.1457 * temperature - 0.0032 * temperature^2) - 0.6443
-        if rand(model.rng) < seed_threshold
-            ## Collect all adjacent position that have no daisies
-            empty_neighbors = Tuple{Int,Int}[]
-            neighbors = nearby_positions(pos, model)
-            for n in neighbors
-                if isempty(ids_in_position(n, model))
-                    push!(empty_neighbors, n)
-                end
-            end
-            if !isempty(empty_neighbors)
-                ## Seed a new daisy in one of those position
-                seeding_place = rand(model.rng, empty_neighbors)
-                add_agent!(seeding_place, model, daisy.breed, 0, daisy.albedo)
+function propagate!(pos, model::DaisyWorld)
+    isempty(pos, model) && return
+    daisy = model[id_in_position(pos, model)]
+    temperature = model.temperature[pos...]
+    seed_threshold = (0.1457 * temperature - 0.0032 * temperature^2) - 0.6443
+    if rand(model.rng) < seed_threshold
+        empty_near_pos = Tuple{Int,Int}[]
+        for near_pos in nearby_positions(pos, model)
+            if isempty(near_pos, model)
+                push!(empty_near_pos, near_pos)
             end
         end
+        if !isempty(empty_near_pos)
+            seeding_place = rand(model.rng, empty_near_pos)
+            add_agent!(seeding_place, model, daisy.breed, 0, daisy.albedo)
+        end
     end
 end
 
 function daisy_step!(agent::Daisy, model::DaisyWorld)
     agent.age += 1
-    agent.age >= model.max_age && kill_agent!(agent, model)
+    agent.age ≥ model.max_age && kill_agent!(agent, model)
 end
 
 function daisyworld_step!(model)
     for p in positions(model)
         update_surface_temperature!(p, model)
         diffuse_temperature!(p, model)
-        propagate_daisy!(p, model)
+        propagate!(p, model)
     end
-    model.tick += 1
+    model.tick[] = model.tick[] + 1
     solar_activity!(model)
 end
 
 function solar_activity!(model::DaisyWorld)
     if model.scenario == :ramp
-        if model.tick > 200 && model.tick <= 400
+        if model.tick[] > 200 && model.tick[] ≤ 400
             model.solar_luminosity += model.solar_change
         end
-        if model.tick > 500 && model.tick <= 750
+        if model.tick[] > 500 && model.tick[] ≤ 750
             model.solar_luminosity -= model.solar_change / 2
         end
     elseif model.scenario == :change
         model.solar_luminosity += model.solar_change
     end
 end
 
+using Random
+
 function daisyworld(;
     griddims = (30, 30),
     max_age = 25,
@@ -107,14 +94,12 @@ function daisyworld(;
 
     rng = MersenneTwister(seed)
     space = GridSpaceSingle(griddims)
-    properties = (;max_age, surface_albedo, solar_luminosity, solar_change, scenario,)
-    properties = Dict(pairs(properties))
-    properties[:tick] = 0
-    properties[:temperature] = zeros(griddims)
+    properties = (;max_age, surface_albedo, solar_luminosity, solar_change, scenario,
+        tick = Ref(0), ratio = 0.5, temperature = zeros(griddims)
+    )
 
     model = ABM(Daisy, space; properties, rng)
 
-    ## Populate with daisies: each position has only one daisy (black or white)
     grid = collect(positions(model))
     num_positions = prod(griddims)
     white_positions =
@@ -131,10 +116,9 @@ function daisyworld(;
         add_agent_pos!(wd, model)
     end
 
-    ## Adjust temperature to initial daisy distribution
     for p in positions(model)
         update_surface_temperature!(p, model)
     end
 
-    return model, daisy_step!, daisyworld_step!
+    return model
 end

src/agents/interaction.jl

@@ -107,6 +107,7 @@ end
 function add_param_sliders!(fig, model, params, resetclick)
     datalayout = fig[end,:][1,2] = GridLayout(tellheight = true)
 
+    # TODO: This needs to become `SliderGrid`.
     slidervals = Dict{Symbol, Observable}()
     for (i, (k, vals)) in enumerate(params)
         startvalue = has_key(model[].properties, k) ?