diff --git a/src/core/simulation.jl b/src/core/simulation.jl index 8ed33f1db1923e160a4bae9455ddb12141d84e93..0a18828b688d6af0eaca99f03215e84bc4c46b9b 100644 --- a/src/core/simulation.jl +++ b/src/core/simulation.jl @@ -12,7 +12,7 @@ This is the heart of the model - a struct that holds all data and state for one simulation run. It is created by [`initialise`](@ref) and passed as input to most model functions. """ -mutable struct AgricultureModel <: SimulationModel +mutable struct AgricultureModel{Tcroptype,Tcropstate} <: SimulationModel settings::Dict{String,Any} rng::AbstractRNG logger::AbstractLogger @@ -21,9 +21,9 @@ mutable struct AgricultureModel <: SimulationModel date::Date landscape::Matrix{Pixel} weather::Dict{Date,Weather} - crops::Dict{String,AbstractCropType} + crops::Dict{String,Tcroptype} farmers::Vector{Farmer} - farmplots::Vector{AbstractFarmPlot} + farmplots::Vector{FarmPlot{Tcropstate}} animals::Vector{Union{Animal,Nothing}} migrants::Vector{Pair{Animal,Date}} events::Vector{FarmEvent} @@ -120,35 +120,38 @@ function initmodel(settings::Dict{String, Any}) # TODO: make this switching on cropmodel simpler if settings["crop.cropmodel"] == "almass" + Tcroptype = ALMaSS.CropType + Tcropstate = ALMaSS.CropState crops = ALMaSS.readcropparameters(settings["crop.cropfile"], settings["crop.growthfile"]) - Tfarmplot = ALMaSS.FarmPlot elseif settings["crop.cropmodel"] == "simple" - # TODO: generalise ALMaSS.CropType? - # crops = Dict{String, ALMaSS.CropType}() - crops = ALMaSS.readcropparameters(settings["crop.cropfile"], - settings["crop.growthfile"]) - Tfarmplot = SimpleCrop.FarmPlot + Tcroptype = SimpleCrop.CropType + Tcropstate = SimpleCrop.CropState + crops_almass = ALMaSS.readcropparameters(settings["crop.cropfile"], + settings["crop.growthfile"]) + crops = Dict(name => SimpleCrop.CropType(ct.name) for (name, ct) in crops_almass) else error("Init for crop model \"$(settings["crop.cropmodel"])\" not implemented") end farmers = Vector{Farmer}() - farmplots = Vector{Tfarmplot}() - model = AgricultureModel(settings, - StableRNG(settings["core.seed"]), - logger, - Vector{DataOutput}(), - Dict{String, DataFrame}(), - settings["core.startdate"], - landscape, - weather, - crops, - farmers, - farmplots, - Vector{Union{Animal,Nothing}}(), - Vector{Pair{Animal, Date}}(), - Vector{FarmEvent}()) + farmplots = Vector{FarmPlot{Tcropstate}}() + model = AgricultureModel{Tcroptype,Tcropstate}( + settings, + StableRNG(settings["core.seed"]), + logger, + Vector{DataOutput}(), + Dict{String, DataFrame}(), + settings["core.startdate"], + landscape, + weather, + crops, + farmers, + farmplots, + Vector{Union{Animal,Nothing}}(), + Vector{Pair{Animal, Date}}(), + Vector{FarmEvent}() + ) saveinputfiles(model) if settings["crop.cropmodel"] == "almass" diff --git a/src/crop/almass.jl b/src/crop/almass.jl index 86175c9a5964e05fbc788591c21645b7f5ececd8..116cd50ca2b17f1084e57b8b442170e7d840d8f2 100644 --- a/src/crop/almass.jl +++ b/src/crop/almass.jl @@ -8,9 +8,8 @@ module ALMaSS using Persefone: - AbstractCropType, - AbstractCropState, EventType, + FarmPlot, Length, m, SimulationModel, @@ -19,7 +18,6 @@ using Persefone: maxtemp, mintemp import Persefone: - Persefone, stepagent!, croptype, cropname, @@ -56,9 +54,7 @@ end The type struct for all crops. Currently follows the crop growth model as implemented in ALMaSS. """ -struct CropType <: AbstractCropType - #TODO make this into an abstract type and create subtypes for different - # crop submodels (#70) +struct CropType name::String minsowdate::Union{Missing,Date} maxsowdate::Union{Missing,Date} @@ -78,7 +74,7 @@ cropname(ct::CropType) = ct.name The state data for an ALMaSS vegetation point calculation as used in FarmPlot. """ -mutable struct CropState <: AbstractCropState +mutable struct CropState #TODO add Unitful croptype::CropType phase::GrowthPhase @@ -94,8 +90,6 @@ croptype(cs::CropState) = cs.croptype cropname(cs::CropState) = cropname(croptype(cs)) cropheight(cs::CropState) = cs.height -const FarmPlot = Persefone.FarmPlot{CropState} - """ Base.tryparse(type, str) @@ -189,7 +183,7 @@ end Update a farm plot by one day. """ -function stepagent!(farmplot::FarmPlot, model::SimulationModel) +function stepagent!(farmplot::FarmPlot{CropState}, model::SimulationModel) # update growing degree days # if no crop-specific base temperature is given, default to 5°C # (https://www.eea.europa.eu/publications/europes-changing-climate-hazards-1/heat-and-cold/heat-and-cold-2014-mean) @@ -213,7 +207,7 @@ end Sow the specified crop on this farmplot. """ -function sow!(cropname::String, farmplot::FarmPlot, model::SimulationModel) +function sow!(cropname::String, farmplot::FarmPlot{CropState}, model::SimulationModel) cs = farmplot.crop_state createevent!(model, farmplot.pixels, sowing) cs.croptype = model.crops[cropname] @@ -226,7 +220,7 @@ end Harvest the crop on this farmplot. """ -function harvest!(farmplot::FarmPlot, model::SimulationModel) +function harvest!(farmplot::FarmPlot{CropState}, model::SimulationModel) cs = farmplot.crop_state createevent!(model, farmplot.pixels, harvesting) cs.phase in [ALMaSS.harvest1, ALMaSS.harvest2] ? @@ -246,7 +240,7 @@ end Apply the relevant crop growth model to update the plants on this farm plot. Currently only supports the ALMaSS crop growth model by Topping et al. """ -function growcrop!(farmplot::FarmPlot, model::SimulationModel) +function growcrop!(farmplot::FarmPlot{CropState}, model::SimulationModel) cs = farmplot.crop_state fertiliser in cs.events ? curve = cs.croptype.lownutrientgrowth : diff --git a/src/crop/croptypes.jl b/src/crop/croptypes.jl index 92bc4baa4d768af80f22e342c37fb50e7a1d4413..37020ce60bd4719f108c5f23ff75083364f1579a 100644 --- a/src/crop/croptypes.jl +++ b/src/crop/croptypes.jl @@ -3,6 +3,6 @@ ### This file contains abstract types that crop models must derive from. ### -abstract type AbstractCropType end +# abstract type AbstractCropType end -abstract type AbstractCropState end +# abstract type AbstractCropState end diff --git a/src/crop/farmplot.jl b/src/crop/farmplot.jl index b570df0c3271a7c7d746561cfa550939289767aa..a3c283ee9b15b5acc048d8ecb34d990b0010956d 100644 --- a/src/crop/farmplot.jl +++ b/src/crop/farmplot.jl @@ -3,218 +3,10 @@ ### This file contains code for the fields that farmers manage. ### -#XXX not sure whether it makes sense to have this as an agent type, -# or perhaps better a grid property? -""" - FarmPlot - -This represents one field, i.e. a collection of pixels with the same management. -This is the spatial unit with which the crop growth model and the farm model work. -""" -mutable struct FarmPlot <: ModelAgent - #TODO add Unitful - const id::Int64 - pixels::Vector{Tuple{Int64, Int64}} - croptype::ALMaSS.CropType - phase::ALMaSS.GrowthPhase - growingdegreedays::Float64 - height::Float64 - LAItotal::Float64 - LAIgreen::Float64 - #biomass::Float64 #XXX I need to figure out how to calculate this - events::Vector{Management} -end - -""" - initfields!(model) - -Initialise the model with its farm plots. -""" -function initfields!(model::SimulationModel) - n = 0 - convertid = Dict{Int64,Int64}() - width, height = size(model.landscape) - for x in 1:width - for y in 1:height - # for each pixel, we need to extract the field ID given by the map input - # file, and convert it into the internal object ID used by Agents.jl, - # creating a new agent object if necessary - rawid = model.landscape[x,y].fieldid - (ismissing(rawid)) && continue - if rawid in keys(convertid) - objectid = convertid[rawid] - model.landscape[x,y].fieldid = objectid - push!(model.farmplots[objectid].pixels, (x,y)) - else - #XXX does this phase calculation work out? - month(model.date) < 3 ? phase = ALMaSS.janfirst : phase = ALMaSS.marchfirst - fp = FarmPlot(length(model.farmplots)+1, [(x,y)], - model.crops["natural grass"], phase, - 0.0, 0.0, 0.0, 0.0, Vector{Management}()) - push!(model.farmplots, fp) - model.landscape[x,y].fieldid = fp.id - convertid[rawid] = fp.id - n += 1 - end - end - end - @info "Initialised $n farm plots." -end - -""" - stepagent!(farmplot, model) - -Update a farm plot by one day. -""" -function stepagent!(farmplot::FarmPlot, model::SimulationModel) - # update growing degree days - # if no crop-specific base temperature is given, default to 5°C - # (https://www.eea.europa.eu/publications/europes-changing-climate-hazards-1/heat-and-cold/heat-and-cold-2014-mean) - basetemp = farmplot.croptype.mingrowthtemp - ismissing(basetemp) && (basetemp = 5.0) - gdd = (maxtemp(model)+mintemp(model))/2 - basetemp - gdd > 0 && (farmplot.growingdegreedays += gdd) - # update the phase on key dates - monthday(model.date) == (1,1) && (farmplot.phase = ALMaSS.janfirst) - monthday(model.date) == (3,1) && (farmplot.phase = ALMaSS.marchfirst) - # update crop growth - growcrop!(farmplot, model) -end - - -## CROP MANAGEMENT AND GROWTH FUNCTIONS - -""" - sow!(cropname, farmplot, model) - -Sow the specified crop on this farmplot. -""" -function sow!(cropname::String, farmplot::FarmPlot, model::SimulationModel) - createevent!(model, farmplot.pixels, sowing) - farmplot.croptype = model.crops[cropname] - farmplot.phase = ALMaSS.sow - #XXX test if the crop is sowable? -end - -""" - harvest!(farmplot, model) - -Harvest the crop on this farmplot. -""" -function harvest!(farmplot::FarmPlot, model::SimulationModel) - createevent!(model, farmplot.pixels, harvesting) - farmplot.phase in [ALMaSS.harvest1, ALMaSS.harvest2] ? - farmplot.phase = ALMaSS.harvest2 : - farmplot.phase = ALMaSS.harvest1 - # height & LAI will be automatically adjusted by the growth function - #TODO calculate and return yield -end - -#TODO fertilise!() -#TODO spray!() -#TODO till!() - -""" - growcrop!(farmplot, model) - -Apply the relevant crop growth model to update the plants on this farm plot. -Currently only supports the ALMaSS crop growth model by Topping et al. -""" -function growcrop!(farmplot::FarmPlot, model::SimulationModel) - fertiliser in farmplot.events ? - curve = farmplot.croptype.lownutrientgrowth : - curve = farmplot.croptype.highnutrientgrowth - points = curve.GDD[farmplot.phase] - for p in 1:length(points) - if points[p] == 99999 - return # the marker that there is no further growth this phase - elseif points[p] == -1 # the marker to set all variables to specified values - farmplot.height = curve.height[farmplot.phase][p] - farmplot.LAItotal = curve.LAItotal[farmplot.phase][p] - farmplot.LAIgreen = curve.LAIgreen[farmplot.phase][p] - return - else - gdd = farmplot.growingdegreedays - # figure out which is the correct slope value to use for growth - if p == length(points) || gdd < points[p+1] - farmplot.height += curve.height[farmplot.phase][p] - farmplot.LAItotal += curve.LAItotal[farmplot.phase][p] - farmplot.LAIgreen += curve.LAIgreen[farmplot.phase][p] - return - end - #XXX To be precise, we ought to figure out if one or more inflection - # points have been passed between yesterday and today, and calculate the - # growth exactly up to the inflection point with the old slope, and onward - # with the new slope. Not doing so will introduce a small amount of error, - # although I think this is acceptable. - end - end -end - - -## UTILITY FUNCTIONS - -""" - averagefieldsize(model) - -Calculate the average field size in hectares for the model landscape. -""" -function averagefieldsize(model::SimulationModel) - conversionfactor = 100 #our pixels are currently 10x10m, so 100 pixels per hectare - sizes::Vector{Float64} = [] - for fp in model.farmplots - push!(sizes, size(fp.pixels)[1]/conversionfactor) - end - round(sum(sizes)/size(sizes)[1], digits=2) -end - -""" - croptype(model, position) - -Return the crop at this position, or nothing if there is no crop here (utility wrapper). -""" -function croptype(pos::Tuple{Int64,Int64}, model::SimulationModel) - ismissing(model.landscape[pos...].fieldid) ? nothing : - model.farmplots[model.landscape[pos...].fieldid].croptype -end - -""" - cropname(model, position) - -Return the name of the crop at this position, or an empty string if there is no crop here -(utility wrapper). -""" -function cropname(pos::Tuple{Int64,Int64}, model::SimulationModel) - field = model.landscape[pos...].fieldid - ismissing(field) ? "" : model.farmplots[field].croptype.name -end - -""" - cropheight(model, position) - -Return the height of the crop at this position, or 0 if there is no crop here -(utility wrapper). -""" -function cropheight(pos::Tuple{Int64,Int64}, model::SimulationModel) - ismissing(model.landscape[pos...].fieldid) ? 0cm : #FIXME can I return something better than 0? - model.farmplots[model.landscape[pos...].fieldid].height*1cm #FIXME units -end - -""" - cropcover(model, position) - -Return the percentage ground cover of the crop at this position, or 0 if there is no crop -here (utility wrapper). -""" -function cropcover(pos::Tuple{Int64,Int64}, model::SimulationModel) - #FIXME LAItotal != ground cover? - ismissing(model.landscape[pos...].fieldid) ? 0 : #FIXME can I return something better than 0? - model.farmplots[model.landscape[pos...].fieldid].LAItotal -end - - +# abstract type AbstractFarmPlot <: ModelAgent end -mutable struct FarmPlot{T <: AbstractCropState} <: AbstractFarmPlot +# mutable struct FarmPlot{T <: AbstractCropState} <: AbstractFarmPlot +mutable struct FarmPlot{T} <: ModelAgent const id::Int64 pixels::Vector{Tuple{Int64, Int64}} crop_state :: T diff --git a/src/crop/simplecrop.jl b/src/crop/simplecrop.jl index 8e074b34c2b9a4aae05c0228f1f8d21c6d41a68d..7c622e7ea7c8571adac8f5fa05b0da21303973f8 100644 --- a/src/crop/simplecrop.jl +++ b/src/crop/simplecrop.jl @@ -1,26 +1,24 @@ module SimpleCrop using Persefone: - AbstractCropType, - AbstractCropState, + FarmPlot, Length, m, SimulationModel, initfields_fill_with! import Persefone: - Persefone, stepagent!, croptype, cropname, cropheight -struct CropType <: AbstractCropType +struct CropType name::String end cropname(ct::CropType) = ct.name -mutable struct CropState <: AbstractCropState +mutable struct CropState croptype::CropType cropheight::Length{Float64} end @@ -29,14 +27,12 @@ croptype(cs::CropState) = cs.croptype cropname(cs::CropState) = cropname(croptype(cs)) cropheight(cs::CropState) = cs.height -const FarmPlot = Persefone.FarmPlot{CropState} - """ stepagent!(farmplot, model) Update a farm plot by one day. """ -function stepagent!(farmplot::FarmPlot, model::SimulationModel) +function stepagent!(farmplot::FarmPlot{CropState}, model::SimulationModel) # TODO: do something simple end @@ -47,7 +43,7 @@ Initialise the model with its farm plots. """ function initfields!(model::SimulationModel) initfields_fill_with!(model) do model, x, y - FarmPlot(length(model.farmplots) + 1, [(x,y)], CropState(CropType("nothing"), 0.0m)) + FarmPlot(length(model.farmplots) + 1, [(x,y)], CropState(model.crops["natural grass"], 0.0m)) end end diff --git a/test/runtests.jl b/test/runtests.jl index b3180102cbbeabba079c448d675672d50e7a4595..e826c76d1fef94a54ed1fec13784faa13678369b 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -43,22 +43,27 @@ function inittestmodel(smallmap=true) TESTSETTINGS["world.weatherfile"]), TESTSETTINGS["core.startdate"], TESTSETTINGS["core.enddate"]) + # TODO: allow init of other crop models besides ALMaSS crops = Ps.ALMaSS.readcropparameters(TESTSETTINGS["crop.cropfile"], TESTSETTINGS["crop.growthfile"]) - model = AgricultureModel(TESTSETTINGS, - StableRNG(TESTSETTINGS["core.seed"]), - global_logger(), - Vector{DataOutput}(), - Dict{String, DataFrame}(), - TESTSETTINGS["core.startdate"], - landscape, - weather, - crops, - Vector{Farmer}(), - Vector{Ps.ALMaSS.FarmPlot}(), - Vector{Union{Animal,Nothing}}(), - Vector{Pair{Animal, Date}}(), - Vector{FarmEvent}()) + farmers = Vector{Farmer}() + farmplots = Vector{Ps.FarmPlot{Ps.ALMaSS.CropState}}() + model = AgricultureModel{Ps.ALMaSS.CropType,Ps.ALMaSS.CropState}( + TESTSETTINGS, + StableRNG(TESTSETTINGS["core.seed"]), + global_logger(), + Vector{DataOutput}(), + Dict{String, DataFrame}(), + TESTSETTINGS["core.startdate"], + landscape, + weather, + crops, + farmers, + farmplots, + Vector{Union{Animal,Nothing}}(), + Vector{Pair{Animal, Date}}(), + Vector{FarmEvent}() + ) model end