From c5c68cab4d43e9f51b9caf04907e958c2da4971a Mon Sep 17 00:00:00 2001 From: Samuel-AMAP Date: Thu, 16 Oct 2025 11:11:02 +0200 Subject: [PATCH 1/8] Simple single-variable hacky ad hoc connector test ; as a model it's awkward --- src/Abstract_model_structs.jl | 45 +++++++++- src/mtg/GraphSimulation.jl | 8 +- src/run.jl | 20 ++++- test/test-simulation.jl | 151 +++++++++++++++++++++++++++++++++- 4 files changed, 217 insertions(+), 7 deletions(-) diff --git a/src/Abstract_model_structs.jl b/src/Abstract_model_structs.jl index 9cef7b08..59405e93 100644 --- a/src/Abstract_model_structs.jl +++ b/src/Abstract_model_structs.jl @@ -25,4 +25,47 @@ model_(m::AbstractModel) = m get_models(m::AbstractModel) = [model_(m)] # Get the models of an AbstractModel # Note: it is returning a vector of models, because in this case the user provided a single model instead of a vector of. get_status(m::AbstractModel) = nothing -get_mapped_variables(m::AbstractModel) = Pair{Symbol,String}[] \ No newline at end of file +get_mapped_variables(m::AbstractModel) = Pair{Symbol,String}[] + + +#using Dates +struct TimestepRange + lower_bound::Period + upper_bound::Period +end + +# Default, no specified range, meaning the model either doesn't depend on time or uses the simulation's default (eg smallest) timestep +TimestepRange() = TimestepRange(Second(0), Second(0)) +# Only a single timestep type possible +TimestepRange(p::Period) = TimestepRange(p, p) + +""" + timestep_range_(tsr::TimestepRange) + +Return the model's valid range for timesteps (which corresponds to the simulation base timestep in the default case). +""" +function timestep_range_(model::AbstractModel) + return TimestepRange() +end + +""" + timestep_valid(tsr::TimestepRange) + +Checks whether a TimestepRange +""" +timestep_valid(tsr::TimestepRange) = tsr.lower_bound <= tsr.upper_bound + +function model_timestep_range_compatible_with_timestep(tsr::TimestepRange, p::Period) + if !timestep_valid(tsr) + return false + end + + # 0 means any timestep is valid, no timestep constraints + if tsr.upper_bound == Seconds(0) + return true + end + + return p >= tsr.lower_bound && p <= tsr.lower_bound +end + +# TODO should i set all timestep ranges to default and hope the modeler gets it right or should i force them to write something ? \ No newline at end of file diff --git a/src/mtg/GraphSimulation.jl b/src/mtg/GraphSimulation.jl index 93376240..3aa491c9 100644 --- a/src/mtg/GraphSimulation.jl +++ b/src/mtg/GraphSimulation.jl @@ -16,7 +16,7 @@ A type that holds all information for a simulation over a graph. - `models`: a dictionary of models - `outputs`: a dictionary of outputs """ -struct GraphSimulation{T,S,U,O,V} +struct GraphSimulation{T,S,U,O,V,W} graph::T statuses::S status_templates::Dict{String,Dict{Symbol,Any}} @@ -26,10 +26,12 @@ struct GraphSimulation{T,S,U,O,V} models::Dict{String,U} outputs::Dict{String,O} outputs_index::Dict{String, Int} + default_timestep::Int # TODO make it a period ? + model_timesteps::Dict{W, Int} #where {W <: AbstractModel} end -function GraphSimulation(graph, mapping; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false) - GraphSimulation(init_simulation(graph, mapping; nsteps=nsteps, outputs=outputs, type_promotion=type_promotion, check=check, verbose=verbose)...) +function GraphSimulation(graph, mapping; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false, default_timestep=1, model_timesteps=Dict{String, Int}()) + GraphSimulation(init_simulation(graph, mapping; nsteps=nsteps, outputs=outputs, type_promotion=type_promotion, check=check, verbose=verbose)..., default_timestep, model_timesteps) end dep(g::GraphSimulation) = g.dependency_graph diff --git a/src/run.jl b/src/run.jl index d3312473..03974654 100644 --- a/src/run.jl +++ b/src/run.jl @@ -365,17 +365,21 @@ function run!( meteo=nothing, constants=PlantMeteo.Constants(), extra=nothing; + orchestrator::Orchestrator=nothing, nsteps=nothing, tracked_outputs=nothing, check=true, - executor=ThreadedEx() + executor=ThreadedEx(), + default_timestep::Int, + model_timesteps::Dict{T, Int} where {T} + ) isnothing(nsteps) && (nsteps = get_nsteps(meteo)) meteo_adjusted = adjust_weather_timesteps_to_given_length(nsteps, meteo) # NOTE : replace_mapping_status_vectors_with_generated_models is assumed to have already run if used # otherwise there might be vector length conflicts with timesteps - sim = GraphSimulation(object, mapping, nsteps=nsteps, check=check, outputs=tracked_outputs) + sim = GraphSimulation(object, mapping, nsteps=nsteps, check=check, outputs=tracked_outputs, default_timestep=default_timestep, model_timesteps=model_timesteps) run!( sim, meteo_adjusted, @@ -455,6 +459,18 @@ function run_node_multiscale!( return nothing end + model_timestep = object.model_timesteps[typeof(node.value)] + + if model_timestep != object.default_timestep + # do accumulation + + + # run if necessary + if i % model_timestep != 0 + return nothing + end + end + node_statuses = status(object)[node.scale] # Get the status of the nodes at the current scale models_at_scale = models[node.scale] diff --git a/test/test-simulation.jl b/test/test-simulation.jl index 4ff83a64..c4439308 100644 --- a/test/test-simulation.jl +++ b/test/test-simulation.jl @@ -292,4 +292,153 @@ end; end end end -end \ No newline at end of file +end + + + + +using PlantSimEngine +# Include the example dummy processes: +using PlantSimEngine.Examples +using Test, Aqua +using Tables, DataFrames, CSV +using MultiScaleTreeGraph +using PlantMeteo, Statistics +using Documenter # for doctests + +using PlantMeteo.Dates +include("helper-functions.jl") + + + +# These models might be worth exposing in the future ? +PlantSimEngine.@process "basic_current_timestep" verbose = false + +struct HelperCurrentTimestepModel <: AbstractBasic_Current_TimestepModel +end + +PlantSimEngine.inputs_(::HelperCurrentTimestepModel) = (next_timestep=1,) +PlantSimEngine.outputs_(m::HelperCurrentTimestepModel) = (current_timestep=1,) + +function PlantSimEngine.run!(m::HelperCurrentTimestepModel, models, status, meteo, constants=nothing, extra=nothing) + status.current_timestep = status.next_timestep + end + + PlantSimEngine.ObjectDependencyTrait(::Type{<:HelperCurrentTimestepModel}) = PlantSimEngine.IsObjectDependent() + PlantSimEngine.TimeStepDependencyTrait(::Type{<:HelperCurrentTimestepModel}) = PlantSimEngine.IsTimeStepDependent() + +PlantSimEngine.timestep_range_(m::HelperCurrentTimestepModel) = Day(1) + + + PlantSimEngine.@process "basic_next_timestep" verbose = false + struct HelperNextTimestepModel <: AbstractBasic_Next_TimestepModel + end + + PlantSimEngine.inputs_(::HelperNextTimestepModel) = (current_timestep=1,) + PlantSimEngine.outputs_(m::HelperNextTimestepModel) = (next_timestep=1,) + + function PlantSimEngine.run!(m::HelperNextTimestepModel, models, status, meteo, constants=nothing, extra=nothing) + status.next_timestep = status.current_timestep + 1 + end + +PlantSimEngine.timestep_range_(m::HelperNextTimestepModel) = Day(1) + + + + + +PlantSimEngine.@process "ToyDay" verbose = false + +struct MyToyDayModel <: AbstractToydayModel end + +PlantSimEngine.inputs_(m::MyToyDayModel) = (a=1,) +PlantSimEngine.outputs_(m::MyToyDayModel) = (daily_temperature=-Inf,) + +function PlantSimEngine.run!(m::MyToyDayModel, models, status, meteo, constants=nothing, extra=nothing) + status.daily_temperature = meteo.T +end + +PlantSimEngine.@process "ToyWeek" verbose = false + +struct MyToyWeekModel <: AbstractToyweekModel + temperature_threshold::Float64 +end + +MyToyWeekModel() = MyToyWeekModel(30.0) +function PlantSimEngine.inputs_(::MyToyWeekModel) + (weekly_max_temperature=-Inf,) +end +PlantSimEngine.outputs_(m::MyToyWeekModel) = (hot = false,) + +function PlantSimEngine.run!(m::MyToyWeekModel, models, status, meteo, constants=nothing, extra=nothing) + status.hot = status.weekly_max_temperature > m.temperature_threshold +end + +PlantSimEngine.timestep_range_(m::MyToyWeekModel) = Week(1) + + + +PlantSimEngine.@process "DWConnector" verbose = false + +struct MyDwconnectorModel <: AbstractDwconnectorModel + T_daily::Array{Float64} +end + +MyDwconnectorModel() = MyDwconnectorModel(Array{Float64}(undef, 7)) + +function PlantSimEngine.inputs_(::MyDwconnectorModel) + (daily_temperature=-Inf, current_timestep=1,) +end +PlantSimEngine.outputs_(m::MyDwconnectorModel) = (weekly_max_temperature = 0.0,) + +function PlantSimEngine.run!(m::MyDwconnectorModel, models, status, meteo, constants=nothing, extra=nothing) + m.T_daily[1 + (status.current_timestep % 7)] = status.daily_temperature + + if(status.current_timestep % 7 == 1) + status.weekly_max_temperature = sum(m.T_daily)/7.0 + else + status.weekly_max_temperature = 0 + end +end + + PlantSimEngine.timestep_range_(m::MyDwconnectorModel) = Day(1) + + + + + +meteo_day = read_weather(joinpath(pkgdir(PlantSimEngine), "examples/meteo_day.csv"), duration=Day) + +m = Dict("Default" => ( + MyToyDayModel(), + MyToyWeekModel(), + MyDwconnectorModel(), + HelperNextTimestepModel(), + MultiScaleModel( + model=HelperCurrentTimestepModel(), + mapped_variables=[PreviousTimeStep(:next_timestep),], + ), + Status(a=1,))) + +to_initialize(m) + +models_timestep = Dict(MyToyDayModel=>1, MyDwconnectorModel => 1, MyToyWeekModel =>7, HelperNextTimestepModel => 1, HelperCurrentTimestepModel => 1) + +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) + +out = run!(mtg, m, meteo_day, default_timestep=1, model_timesteps=models_timestep) + +@testset "Test varying timestep" begin + + + @test + @test + +end + + + # NOTE : replace_mapping_status_vectors_with_generated_models is assumed to have already run if used + # otherwise there might be vector length conflicts with timesteps + sim = @enter PlantSimEngine.GraphSimulation(mtg, m, nsteps=nothing, check=true, outputs=nothing, default_timestep=1, model_timesteps=models_timestep) + +using PlantSimEngine From a3c89735b6d6c193de6d76b19c2573f511b3729a Mon Sep 17 00:00:00 2001 From: Samuel-AMAP Date: Thu, 16 Oct 2025 11:15:38 +0200 Subject: [PATCH 2/8] Start to implement an orchestrator struct that centralizes timestep issues and connector data --- src/PlantSimEngine.jl | 5 ++++ src/mtg/GraphSimulation.jl | 11 ++++----- src/mtg/initialisation.jl | 2 +- src/timestep/timestep_mapping.jl | 42 ++++++++++++++++++++++++++++++++ 4 files changed, 53 insertions(+), 7 deletions(-) create mode 100644 src/timestep/timestep_mapping.jl diff --git a/src/PlantSimEngine.jl b/src/PlantSimEngine.jl index ebf59dcb..b31d9847 100644 --- a/src/PlantSimEngine.jl +++ b/src/PlantSimEngine.jl @@ -26,6 +26,7 @@ import Statistics import SHA: sha1 using PlantMeteo +using PlantMeteo.Dates # UninitializedVar + PreviousTimeStep: include("variables_wrappers.jl") @@ -65,6 +66,9 @@ include("dependencies/printing.jl") include("dependencies/dependencies.jl") include("dependencies/get_model_in_dependency_graph.jl") +# Timesteps. : +include("timestep/timestep_mapping.jl") + # MTG compatibility: include("mtg/GraphSimulation.jl") include("mtg/mapping/getters.jl") @@ -103,6 +107,7 @@ include("examples_import.jl") export PreviousTimeStep export AbstractModel export ModelList, MultiScaleModel +export Orchestrator export RMSE, NRMSE, EF, dr export Status, TimeStepTable, status export init_status! diff --git a/src/mtg/GraphSimulation.jl b/src/mtg/GraphSimulation.jl index 3aa491c9..001a7f3f 100644 --- a/src/mtg/GraphSimulation.jl +++ b/src/mtg/GraphSimulation.jl @@ -14,9 +14,10 @@ A type that holds all information for a simulation over a graph. - `var_need_init`: a dictionary indicating if a variable needs to be initialized - `dependency_graph`: the dependency graph of the models applied to the graph - `models`: a dictionary of models +- `Orchestrator : the structure that handles timestep peculiarities - `outputs`: a dictionary of outputs """ -struct GraphSimulation{T,S,U,O,V,W} +struct GraphSimulation{T,S,U,O,V} graph::T statuses::S status_templates::Dict{String,Dict{Symbol,Any}} @@ -24,14 +25,12 @@ struct GraphSimulation{T,S,U,O,V,W} var_need_init::Dict{String,V} dependency_graph::DependencyGraph models::Dict{String,U} + orchestrator::Orchestrator outputs::Dict{String,O} - outputs_index::Dict{String, Int} - default_timestep::Int # TODO make it a period ? - model_timesteps::Dict{W, Int} #where {W <: AbstractModel} end -function GraphSimulation(graph, mapping; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false, default_timestep=1, model_timesteps=Dict{String, Int}()) - GraphSimulation(init_simulation(graph, mapping; nsteps=nsteps, outputs=outputs, type_promotion=type_promotion, check=check, verbose=verbose)..., default_timestep, model_timesteps) +function GraphSimulation(graph, mapping; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false, orchestrator=Orchestrator()) + GraphSimulation(init_simulation(graph, mapping; nsteps=nsteps, outputs=outputs, type_promotion=type_promotion, check=check, verbose=verbose)..., orchestrator) end dep(g::GraphSimulation) = g.dependency_graph diff --git a/src/mtg/initialisation.jl b/src/mtg/initialisation.jl index 57aef81c..95fb6e5c 100644 --- a/src/mtg/initialisation.jl +++ b/src/mtg/initialisation.jl @@ -305,7 +305,7 @@ The value is not a reference to the one in the attribute of the MTG, but a copy a value in a Dict. If you need a reference, you can use a `Ref` for your variable in the MTG directly, and it will be automatically passed as is. """ -function init_simulation(mtg, mapping; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false) +function init_simulation(mtg, mapping; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false, orchestrator=Orchestrator()) # Ensure the user called the model generation function to handle vectors passed into a status # before we keep going diff --git a/src/timestep/timestep_mapping.jl b/src/timestep/timestep_mapping.jl new file mode 100644 index 00000000..e30fc97c --- /dev/null +++ b/src/timestep/timestep_mapping.jl @@ -0,0 +1,42 @@ + +# Those names all suck, need to change them +# Some of them are probably not ideal for new users, too + +# Some types can also be constrained a lot more, probably + +struct TimestepMapper#{V} + variable_from#::V + timestep_from::Int + mapping_function +end + +struct SimulationTimestepHandler#{W,V} + model_timesteps::Dict{Any, Int} # where {W <: AbstractModel} # if a model isn't in there, then it follows the default, todo check if the given timestep respects the model's range + timestep_variable_mapping::Dict{Any, TimestepMapper} #where {V} +end + +SimulationTimestepHandler() = SimulationTimestepHandler(Dict(), Dict()) #Dict{W, Int}(), Dict{V, TimestepMapper}()) where {W, V} + +mutable struct Orchestrator + # This is actually a general simulation parameter, not-scale specific + # todo change to Period + default_timestep::Int64 + + # This needs to be per-scale : if a model is used at two different scales, + # and the same variable of that model maps to some other timestep to two *different* variables + # then I believe we can only rely on the different scale to disambiguate + non_default_timestep_data_per_scale::Dict{String, SimulationTimestepHandler} + + function Orchestrator(default::Int64, per_scale::Dict{String, SimulationTimestepHandler}) + @assert default >= 0 "The default_timestep should be greater than or equal to 0." + return new(default, per_scale) + end +end + +# TODO have a default constructor take in a meteo or something, and set up the default timestep automagically to be the finest weather timestep +# Other options are possible +Orchestrator() = Orchestrator(1, Dict{String, SimulationTimestepHandler}()) + + +#o = Orchestrator() +#oo = Orchestrator(1, Dict{String, SimulationTimestepHandler}()) \ No newline at end of file From cb77476e7288f3d0cffcae946c04e7c18c41d3f5 Mon Sep 17 00:00:00 2001 From: Samuel-AMAP Date: Tue, 28 Oct 2025 15:03:51 +0100 Subject: [PATCH 3/8] Many changes : new orchestrator attempt, many more fixes. Multiple timesteps handled, but bad meshing with refvalue/refvectors cause overwrites if we use multiscale mappings, and unexplored issues if we avoid them. Also, many-node to many-node issues if multiscale mapping. Some TODO comments, structs, commented out code and notes are outdated. More exploration and cleanup required, but current state is much more interesting than the previous commit. --- src/Abstract_model_structs.jl | 6 +- src/PlantSimEngine.jl | 8 +- src/dependencies/dependencies.jl | 6 +- src/dependencies/dependency_graph.jl | 12 + src/dependencies/hard_dependencies.jl | 42 ++- src/dependencies/soft_dependencies.jl | 101 +++++- src/mtg/GraphSimulation.jl | 8 +- src/mtg/add_organ.jl | 5 + src/mtg/initialisation.jl | 120 ++++++- src/mtg/mapping/compute_mapping.jl | 4 +- src/mtg/mapping/reverse_mapping.jl | 2 +- src/processes/model_initialisation.jl | 2 +- src/run.jl | 85 +++-- src/timestep/timestep_mapping.jl | 93 +++++- test/test_multitimestep.jl | 457 ++++++++++++++++++++++++++ 15 files changed, 887 insertions(+), 64 deletions(-) create mode 100644 test/test_multitimestep.jl diff --git a/src/Abstract_model_structs.jl b/src/Abstract_model_structs.jl index 59405e93..f8a079e0 100644 --- a/src/Abstract_model_structs.jl +++ b/src/Abstract_model_structs.jl @@ -50,18 +50,16 @@ end """ timestep_valid(tsr::TimestepRange) - -Checks whether a TimestepRange """ timestep_valid(tsr::TimestepRange) = tsr.lower_bound <= tsr.upper_bound -function model_timestep_range_compatible_with_timestep(tsr::TimestepRange, p::Period) +function is_timestep_in_range(tsr::TimestepRange, p::Period) if !timestep_valid(tsr) return false end # 0 means any timestep is valid, no timestep constraints - if tsr.upper_bound == Seconds(0) + if tsr.upper_bound == Second(0) return true end diff --git a/src/PlantSimEngine.jl b/src/PlantSimEngine.jl index b31d9847..5ffaea10 100644 --- a/src/PlantSimEngine.jl +++ b/src/PlantSimEngine.jl @@ -57,6 +57,9 @@ include("component_models/get_status.jl") # Transform into a dataframe: include("dataframe.jl") +# Timesteps. : +include("timestep/timestep_mapping.jl") + # Computing model dependencies: include("dependencies/soft_dependencies.jl") include("dependencies/hard_dependencies.jl") @@ -66,9 +69,6 @@ include("dependencies/printing.jl") include("dependencies/dependencies.jl") include("dependencies/get_model_in_dependency_graph.jl") -# Timesteps. : -include("timestep/timestep_mapping.jl") - # MTG compatibility: include("mtg/GraphSimulation.jl") include("mtg/mapping/getters.jl") @@ -107,7 +107,7 @@ include("examples_import.jl") export PreviousTimeStep export AbstractModel export ModelList, MultiScaleModel -export Orchestrator +export Orchestrator, Orchestrator2, TimestepRange, Var_to, Var_from, ModelTimestepMapping export RMSE, NRMSE, EF, dr export Status, TimeStepTable, status export init_status! diff --git a/src/dependencies/dependencies.jl b/src/dependencies/dependencies.jl index 9252e7b0..27b9fcde 100644 --- a/src/dependencies/dependencies.jl +++ b/src/dependencies/dependencies.jl @@ -100,16 +100,16 @@ function dep(m::NamedTuple, nsteps=1; verbose::Bool=true) dep(nsteps; verbose=verbose, m...) end -function dep(mapping::Dict{String,T}; verbose::Bool=true) where {T} +function dep(mapping::Dict{String,T}; verbose::Bool=true, orchestrator=Orchestrator2()) where {T} # First step, get the hard-dependency graph and create SoftDependencyNodes for each hard-dependency root. In other word, we want # only the nodes that are not hard-dependency of other nodes. These nodes are taken as roots for the soft-dependency graph because they # are independant. - soft_dep_graphs_roots, hard_dep_dict = hard_dependencies(mapping; verbose=verbose) + soft_dep_graphs_roots, hard_dep_dict = hard_dependencies(mapping; verbose=verbose, orchestrator=Orchestrator2()) # Second step, compute the soft-dependency graph between SoftDependencyNodes computed in the first step. To do so, we search the # inputs of each process into the outputs of the other processes, at the same scale, but also between scales. Then we keep only the # nodes that have no soft-dependencies, and we set them as root nodes of the soft-dependency graph. The other nodes are set as children # of the nodes that they depend on. - dep_graph = soft_dependencies_multiscale(soft_dep_graphs_roots, mapping, hard_dep_dict) + dep_graph = soft_dependencies_multiscale(soft_dep_graphs_roots, mapping, hard_dep_dict, orchestrator=orchestrator) # During the building of the soft-dependency graph, we identified the inputs and outputs of each dependency node, # and also defined **inputs** as MappedVar if they are multiscale, i.e. if they take their values from another scale. # What we are missing is that we need to also define **outputs** as multiscale if they are needed by another scale. diff --git a/src/dependencies/dependency_graph.jl b/src/dependencies/dependency_graph.jl index 2be06b64..33b5fc62 100644 --- a/src/dependencies/dependency_graph.jl +++ b/src/dependencies/dependency_graph.jl @@ -12,6 +12,16 @@ mutable struct HardDependencyNode{T} <: AbstractDependencyNode children::Vector{HardDependencyNode} end +mutable struct TimestepMapping + variable_from::Symbol + variable_to::Symbol + node_to # SoftDependencyNode causes a circular reference # TODO could it be a harddependencynode... ? + mapping_function::Function + mapping_data_template + mapping_data::Dict{Int, Any} # TODO Any's type is the variable's type, also, is Int good here ? Prob not +end + +# can hard dependency nodes also handle timestep mapped variables... ? mutable struct SoftDependencyNode{T} <: AbstractDependencyNode value::T process::Symbol @@ -23,6 +33,8 @@ mutable struct SoftDependencyNode{T} <: AbstractDependencyNode parent_vars::Union{Nothing,NamedTuple} children::Vector{SoftDependencyNode} simulation_id::Vector{Int} # id of the simulation + timestep::Period + timestep_mapping_data::Union{Nothing, Vector{TimestepMapping}} # TODO : this approach might not play too well with parallelisation over MTG nodes end # Add methods to check if a node is parallelizable: diff --git a/src/dependencies/hard_dependencies.jl b/src/dependencies/hard_dependencies.jl index 0d1fa5b3..eda0573f 100644 --- a/src/dependencies/hard_dependencies.jl +++ b/src/dependencies/hard_dependencies.jl @@ -112,7 +112,7 @@ end # When we use a mapping (multiscale), we return the set of soft-dependencies (we put the hard-dependencies as their children): -function hard_dependencies(mapping::Dict{String,T}; verbose::Bool=true) where {T} +function hard_dependencies(mapping::Dict{String,T}; verbose::Bool=true, orchestrator::Orchestrator2=Orchestrator2()) where {T} full_vars_mapping = Dict(first(mod) => Dict(get_mapped_variables(last(mod))) for mod in mapping) soft_dep_graphs = Dict{String,Any}() not_found = Dict{Symbol,DataType}() @@ -226,6 +226,42 @@ function hard_dependencies(mapping::Dict{String,T}; verbose::Bool=true) where {T end end +#= + # TODO check whether this is a bit late in the game + # maybe the timestep mapping should be done before we enter this function + if length(orchestrator.non_default_timestep_data_per_scale) > 0 + if haskey(orchestrator.non_default_timestep_data_per_scale, symbol(node)) + tvm = orchestrator.non_default_timestep_data_per_scale[symbol(node)].timestep_variable_mapping + if haskey(twm, var) + + end + end + end + error("Variable `$(var)` is not computed by any model, not mapped from a different scale or timestep not initialised by the user in the status, and not found in the MTG at scale $(symbol(node)) (checked for MTG node $(node_id(node))).") +=# + + + # Once multiscale mapping has been dealt with, check if any variable has a timestep mapping + # Which will add potential new dependencies + #=if !isempty(orchestrator.non_default_timestep_data_per_scale) + # TODO the user can get away with not declaring the model, only the scale if necessary + # a prepass that recomputes everything might simplify code here and make the simulation require less variable digging + for (scale, tsh) in non_default_timestep_data_per_scale + # TODO find which model the variable is pulled from + # TODO check the variable exists + for (model, timestep) in tsh.model_timesteps + # TODO check the timestep is within the model's accepted timestep range + # TODO recover the right variables + end + + for (variable, tvm) in tsh.timestep_variable_mapping + # TODO check the variable isn't already mapped + # If it is, ensure there are no name conflicts + # and the model of the variable it is taken from has the expected timestep + # If it isn't, create a new link + end + end + end=# for (organ, model) in mapping soft_dep_graph = Dict( process_ => SoftDependencyNode( @@ -238,7 +274,9 @@ function hard_dependencies(mapping::Dict{String,T}; verbose::Bool=true) where {T nothing, nothing, SoftDependencyNode[], - [0] # Vector of zeros of length = number of time-steps + [0], # Vector of zeros of length = number of time-steps + orchestrator.default_timestep, + nothing ) for (process_, soft_dep_vars) in hard_deps[organ].roots # proc_ = :carbon_assimilation ; soft_dep_vars = hard_deps.roots[proc_] ) diff --git a/src/dependencies/soft_dependencies.jl b/src/dependencies/soft_dependencies.jl index b772a207..241a63bf 100644 --- a/src/dependencies/soft_dependencies.jl +++ b/src/dependencies/soft_dependencies.jl @@ -70,7 +70,9 @@ function soft_dependencies(d::DependencyGraph{Dict{Symbol,HardDependencyNode}}, nothing, nothing, SoftDependencyNode[], - fill(0, nsteps) + fill(0, nsteps), + Day(1), # TODO + nothing ) for (process_, soft_dep_vars) in d.roots ) @@ -138,8 +140,18 @@ function soft_dependencies(d::DependencyGraph{Dict{Symbol,HardDependencyNode}}, return DependencyGraph(independant_process_root, d.not_found) end +function timestep_mapped_variables(orchestrator) + +#=struct SimulationTimestepHandler#{W,V} + model_timesteps::Dict{Any, Period} # where {W <: AbstractModel} # if a model isn't in there, then it follows the default, todo check if the given timestep respects the model's range + timestep_variable_mapping::Dict{Any, TimestepMapper} #where {V} +end + non_default_timestep_data_per_scale::Dict{String, SimulationTimestepHandler} +=# +end + # For multiscale mapping: -function soft_dependencies_multiscale(soft_dep_graphs_roots::DependencyGraph{Dict{String,Any}}, mapping::Dict{String,A}, hard_dep_dict::Dict{Pair{Symbol,String},HardDependencyNode}) where {A<:Any} +function soft_dependencies_multiscale(soft_dep_graphs_roots::DependencyGraph{Dict{String,Any}}, mapping::Dict{String,A}, hard_dep_dict::Dict{Pair{Symbol,String},HardDependencyNode}; orchestrator::Orchestrator2=Orchestrator2()) where {A<:Any} mapped_vars = mapped_variables(mapping, soft_dep_graphs_roots, verbose=false) rev_mapping = reverse_mapping(mapped_vars, all=false) @@ -324,10 +336,93 @@ function soft_dependencies_multiscale(soft_dep_graphs_roots::DependencyGraph{Dic end end - return DependencyGraph(independant_process_root, soft_dep_graphs_roots.not_found) + dep_graph = DependencyGraph(independant_process_root, soft_dep_graphs_roots.not_found) + traverse_dependency_graph!(x -> set_non_default_timestep_in_node(x, orchestrator), dep_graph, visit_hard_dep=false) + traverse_dependency_graph!(x -> add_timestep_data_to_node(x, orchestrator), dep_graph, visit_hard_dep=false) + + return dep_graph +end + +# # set the timestep for everyone first, else we might not use the correct timestep when looking at the parents later +function set_non_default_timestep_in_node(soft_dependency_node, orchestrator::Orchestrator2) + for mtsm in orchestrator.non_default_timestep_mapping + if mtsm.scale == soft_dependency_node.scale && (mtsm.model) == typeof(model_(soft_dependency_node.value)) + soft_dependency_node.timestep = mtsm.timestep + end + end +end + +function add_timestep_data_to_node(soft_dependency_node, orchestrator::Orchestrator2) + + # now we can create the mapping + for mtsm in orchestrator.non_default_timestep_mapping + if mtsm.scale == soft_dependency_node.scale && (mtsm.model) == typeof(model_(soft_dependency_node.value)) + for (var_to, var_from) in mtsm.var_to_var + if !isnothing(soft_dependency_node.parent) + parent = nothing + variable_mapping = nothing + for parent_node in soft_dependency_node.parent + if typeof(parent_node.value) == var_from.model && parent_node.scale == var_from.scale + parent = parent_node + variable_mapping = create_timestep_mapping(soft_dependency_node, parent, var_to, var_from) + break + end + end + if isnothing(parent) + #error + end + if isnothing(parent.timestep_mapping_data) + parent.timestep_mapping_data = Vector{TimestepMapping}() + end + push!(parent.timestep_mapping_data, variable_mapping) + else + # Error + end + end + end + end end +# TODO this is incorrect, there may be multiple variables mapped between the two nodes +function create_timestep_mapping(node::SoftDependencyNode, parent::SoftDependencyNode, var_to::Var_to, var_from::Var_from) + + @assert parent.timestep != 0 "Error : node timestep internally set to 0" + + timestep_ratio = node.timestep / parent.timestep + + # Keeping things simple for now, only integers allowed + @assert timestep_ratio == trunc(timestep_ratio) "Error : non-integer timestep ratio" + + # TODO ensure type compatibility between var_to and var_from + # Simplification probably possible by doing the check earlier + + # TODO test previoustimestep + var_type = DataType + + for (symbol, var_dump) in node.inputs + for var in var_dump + if isa(var, MappedVar) + # check the source variable, because the sink one might be a vector...? + # TODO multinode mapping + if var.source_variable == var_from.name + # This should be a fixed size array, ideally + var_type = eltype(mapped_default(var)) + break + end + else + if var.symbol == var_from.name + @assert "untested" + var_type = eltype(mapped_default(var)) + break + end + end + end + end + mapping_data_template = Vector{var_type}(undef, convert(Int64, timestep_ratio)) + # TODO : type shouldn't be Any but Vector{var_type} + return TimestepMapping(var_from.name, var_to.name, node, var_from.mapping_function, mapping_data_template, Dict{MultiScaleTreeGraph.NodeMTG, Any}()) +end """ drop_process(proc_vars, process) diff --git a/src/mtg/GraphSimulation.jl b/src/mtg/GraphSimulation.jl index 001a7f3f..efe17ec5 100644 --- a/src/mtg/GraphSimulation.jl +++ b/src/mtg/GraphSimulation.jl @@ -25,12 +25,14 @@ struct GraphSimulation{T,S,U,O,V} var_need_init::Dict{String,V} dependency_graph::DependencyGraph models::Dict{String,U} - orchestrator::Orchestrator outputs::Dict{String,O} + outputs_index::Dict{String, Int} + orchestrator::Orchestrator2 + end -function GraphSimulation(graph, mapping; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false, orchestrator=Orchestrator()) - GraphSimulation(init_simulation(graph, mapping; nsteps=nsteps, outputs=outputs, type_promotion=type_promotion, check=check, verbose=verbose)..., orchestrator) +function GraphSimulation(graph, mapping; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false, orchestrator=Orchestrator2()) + GraphSimulation(init_simulation(graph, mapping; nsteps=nsteps, outputs=outputs, type_promotion=type_promotion, check=check, verbose=verbose, orchestrator=orchestrator)...) end dep(g::GraphSimulation) = g.dependency_graph diff --git a/src/mtg/add_organ.jl b/src/mtg/add_organ.jl index 1c2925ef..58c71a9e 100644 --- a/src/mtg/add_organ.jl +++ b/src/mtg/add_organ.jl @@ -33,5 +33,10 @@ function add_organ!(node::MultiScaleTreeGraph.Node, sim_object, link, symbol, sc new_node = MultiScaleTreeGraph.Node(id, node, MultiScaleTreeGraph.NodeMTG(link, symbol, index, scale), attributes) st = init_node_status!(new_node, sim_object.statuses, sim_object.status_templates, sim_object.reverse_multiscale_mapping, sim_object.var_need_init, check=check) + # TODO add the node to the timestep mappings + # TODO initialise the MTG nodes in the timestep mappings + # NOTE : this isn't ideal, as it constrains the add_organ! function usage + init_timestep_mapping_data(new_node, sim_object.dependency_graph) + return st end \ No newline at end of file diff --git a/src/mtg/initialisation.jl b/src/mtg/initialisation.jl index 95fb6e5c..eac56d0c 100644 --- a/src/mtg/initialisation.jl +++ b/src/mtg/initialisation.jl @@ -21,7 +21,7 @@ a dictionary of variables that need to be initialised or computed by other model `(;statuses, status_templates, reverse_multiscale_mapping, vars_need_init, nodes_with_models)` """ -function init_statuses(mtg, mapping, dependency_graph=first(hard_dependencies(mapping; verbose=false)); type_promotion=nothing, verbose=false, check=true) +function init_statuses(mtg, mapping, dependency_graph=first(hard_dependencies(mapping; verbose=false, orchestrator=Orchestrator2())); type_promotion=nothing, verbose=false, check=true, orchestrator=Orchestrator2()) # We compute the variables mapping for each scale: mapped_vars = mapped_variables(mapping, dependency_graph, verbose=verbose) @@ -38,9 +38,17 @@ function init_statuses(mtg, mapping, dependency_graph=first(hard_dependencies(ma convert_reference_values!(mapped_vars) # Get the variables that are not initialised or computed by other models in the output: - vars_need_init = Dict(org => filter(x -> isa(last(x), UninitializedVar), vars) |> keys for (org, vars) in mapped_vars) |> + vars_need_init = Dict(org => filter(x -> isa(last(x), UninitializedVar), vars) |> keys |> collect for (org, vars) in mapped_vars) |> filter(x -> length(last(x)) > 0) + # Filter out variables that are timestep-mapped by the user, + # as those variables are initialized by another model, but are currently flagged as needing initialization + # At this stage, data present in the orchestrator is expected to be valid, so we can take it into account + # A model with a different timestep can still have unitialized vars found in a node, the meteo, or to be initialized by the user + # in which case it'll be absent from the timestep mapping, but this needs testing + #filter_timestep_mapped_variables!(vars_need_init, orchestrator) + + # Note: these variables may be present in the MTG attributes, we check that below when traversing the MTG. # We traverse the MTG to initialise the statuses linked to the nodes: @@ -91,7 +99,7 @@ The `check` argument is a boolean indicating if variables initialisation should in the node attributes (using the variable name). If `true`, the function returns an error if the attribute is missing, otherwise it uses the default value from the model. """ -function init_node_status!(node, statuses, mapped_vars, reverse_multiscale_mapping, vars_need_init=Dict{String,Any}(), type_promotion=nothing; check=true, attribute_name=:plantsimengine_status) +function init_node_status!(node, statuses, mapped_vars, reverse_multiscale_mapping, vars_need_init=Dict{String,Any}(), type_promotion=nothing; check=true, attribute_name=:plantsimengine_status, orchestrator=Orchestrator2()) # Check if the node has a model defined for its symbol, if not, no need to compute symbol(node) ∉ collect(keys(mapped_vars)) && return @@ -113,7 +121,6 @@ function init_node_status!(node, statuses, mapped_vars, reverse_multiscale_mappi end continue end - error("Variable `$(var)` is not computed by any model, not initialised by the user in the status, and not found in the MTG at scale $(symbol(node)) (checked for MTG node $(node_id(node))).") end # Applying the type promotion to the node attribute if needed: if isnothing(type_promotion) @@ -305,7 +312,7 @@ The value is not a reference to the one in the attribute of the MTG, but a copy a value in a Dict. If you need a reference, you can use a `Ref` for your variable in the MTG directly, and it will be automatically passed as is. """ -function init_simulation(mtg, mapping; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false, orchestrator=Orchestrator()) +function init_simulation(mtg, mapping; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false, orchestrator=Orchestrator2()) # Ensure the user called the model generation function to handle vectors passed into a status # before we keep going @@ -314,9 +321,11 @@ function init_simulation(mtg, mapping; nsteps=1, outputs=nothing, type_promotion @assert false "Error : Mapping status at $organ_with_vector level contains a vector. If this was intentional, call the function generate_models_from_status_vectors on your mapping before calling run!. And bear in mind this is not meant for production. If this wasn't intentional, then it's likely an issue on the mapping definition, or an unusual model." end + # preliminary_check_timestep_data(mapping, orchestrator) + # Get the status of each node by node type, pre-initialised considering multi-scale variables: statuses, status_templates, reverse_multiscale_mapping, vars_need_init = - init_statuses(mtg, mapping, first(hard_dependencies(mapping; verbose=false)); type_promotion=type_promotion, verbose=verbose, check=check) + init_statuses(mtg, mapping, first(hard_dependencies(mapping; verbose=false, orchestrator=orchestrator)); type_promotion=type_promotion, verbose=verbose, check=check, orchestrator=orchestrator) # Print an info if models are declared for nodes that don't exist in the MTG: if check && any(x -> length(last(x)) == 0, statuses) @@ -329,5 +338,102 @@ function init_simulation(mtg, mapping; nsteps=1, outputs=nothing, type_promotion outputs = pre_allocate_outputs(statuses, status_templates, reverse_multiscale_mapping, vars_need_init, outputs, nsteps, type_promotion=type_promotion, check=check) outputs_index = Dict{String, Int}(s => 1 for s in keys(outputs)) - return (; mtg, statuses, status_templates, reverse_multiscale_mapping, vars_need_init, dependency_graph=dep(mapping, verbose=verbose), models, outputs, outputs_index) + + dependency_graph = dep(mapping, verbose=verbose, orchestrator=orchestrator) + + # Samuel : Once the dependency graph is created, and the timestep mappings are added into it + # We need to register the existing MTG nodes to initialize their individual data + # The current implementation is heavy, it may be quite slow for MTGs that already contain many nodes + # A faster way would be to count nodes by type once, store their ids, and only traverse the dependency graph once to add them + MultiScaleTreeGraph.traverse!(mtg, init_timestep_mapping_data, dependency_graph) + + return (; mtg, statuses, status_templates, reverse_multiscale_mapping, vars_need_init, dependency_graph=dependency_graph, models, outputs, outputs_index, orchestrator) +end + +function preliminary_check_timestep_data(mapping, orchestrator) + + if isempty(orchestrator.non_default_timestep_data_per_scale) + return + end + + # First, check timesteps are within models' accepted ranges + for (organ, models_status) in mapping + models = get_models(models_status) + + for model in models + checked = false + + if haskey(orchestrator.non_default_timestep_data_per_scale, organ) + tsh = orchestrator.non_default_timestep_data_per_scale[organ] + + if typeof(model) in collect(keys(tsh.model_timesteps)) + checked = true + # Check the timestep is within the model's accepted timestep range + if !is_timestep_in_range(timestep_range_(model), tsh.model_timesteps[typeof(model)]) + # TODO return error + end + end + end + # if it wasn't found, it means it's a model set to the default timestep + if !checked + if !is_timestep_in_range(timestep_range_(model), orchestrator.default_timestep) + # TODO return error + end + end + end + end + + # Next, check timestep mapped variables : + # They should all exist as input/output somewhere (not dealing with mtg stuff for now) + # If they are already mapped in a MultiScaleModel, there should be no differences + # They should not cause name conflicts + # The scales should all be in the mapping + # They should map to different timesteps, and if the timestep isn't the default, then the downstream model should be in the list as well + for (organ, tsh) in orchestrator.non_default_timestep_data_per_scale + #if !organ in collect(keys(mapping)) + # TODO error + #end + + for (model, tvm) in tsh.model_timesteps + end + end> + + return +end + + + + +#= +function filter_timestep_mapped_variables!(vars_need_init, orchestrator) + for (org, vars) in vars_need_init + if haskey(orchestrator.non_default_timestep_data_per_scale, org) + + filter!(o -> !haskey(orchestrator.non_default_timestep_data_per_scale[org].timestep_variable_mapping, o), vars) + #for var in vars + # if haskey(orchestrator.non_default_timestep_data_per_scale[org].timestep_variable_mapping, var) + # delete!(vars, var) + # end + #end + if isempty(vars) + delete!(vars_need_init, org) + end + end + end +end=# + +function filter_timestep_mapped_variables!(vars_need_init, orchestrator) + for (org, vars) in vars_need_init + for tmst in orchestrator.non_default_timestep_mapping + if tmst.scale == org + + for (var_to, var_from) in tmst.var_to_var + filter!(o -> o == var_to.name || o == var_from.name, vars) + end + if isempty(vars) + delete!(vars_need_init, org) + end + end + end + end end \ No newline at end of file diff --git a/src/mtg/mapping/compute_mapping.jl b/src/mtg/mapping/compute_mapping.jl index 21d64008..029b9c3d 100644 --- a/src/mtg/mapping/compute_mapping.jl +++ b/src/mtg/mapping/compute_mapping.jl @@ -11,7 +11,7 @@ However, models that are identified as hard-dependencies are not given individua nodes under other models. - `verbose::Bool`: whether to print the stacktrace of the search for the default value in the mapping. """ -function mapped_variables(mapping, dependency_graph=first(hard_dependencies(mapping; verbose=false)); verbose=false) +function mapped_variables(mapping, dependency_graph=first(hard_dependencies(mapping; verbose=false, orchestrator=Orchestrator2())); verbose=false) # Initialise a dict that defines the multiscale variables for each organ type: mapped_vars = mapped_variables_no_outputs_from_other_scale(mapping, dependency_graph) @@ -54,7 +54,7 @@ This function returns a dictionary with the (multiscale-) inputs and outputs var Note that this function does not include the variables that are outputs from another scale and not computed by this scale, see `mapped_variables_with_outputs_as_inputs` for that. """ -function mapped_variables_no_outputs_from_other_scale(mapping, dependency_graph=first(hard_dependencies(mapping; verbose=false))) +function mapped_variables_no_outputs_from_other_scale(mapping, dependency_graph=first(hard_dependencies(mapping; verbose=false, orchestrator=Orchestrator2()))) nodes_insouts = Dict(organ => (inputs=ins, outputs=outs) for (organ, (soft_dep_graph, ins, outs)) in dependency_graph.roots) ins = Dict{String,NamedTuple}(organ => flatten_vars(vcat(values(ins)...)) for (organ, (ins, outs)) in nodes_insouts) outs = Dict{String,NamedTuple}(organ => flatten_vars(vcat(values(outs)...)) for (organ, (ins, outs)) in nodes_insouts) diff --git a/src/mtg/mapping/reverse_mapping.jl b/src/mtg/mapping/reverse_mapping.jl index 059b4294..71ba46df 100644 --- a/src/mtg/mapping/reverse_mapping.jl +++ b/src/mtg/mapping/reverse_mapping.jl @@ -69,7 +69,7 @@ Dict{String, Dict{String, Dict{Symbol, Any}}} with 3 entries: """ function reverse_mapping(mapping::Dict{String,T}; all=true) where {T<:Any} # Method for the reverse mapping applied directly on the mapping (not used in the code base) - mapped_vars = mapped_variables(mapping, first(hard_dependencies(mapping; verbose=false)), verbose=false) + mapped_vars = mapped_variables(mapping, first(hard_dependencies(mapping; verbose=false, orchestrator=Orchestrator2())), verbose=false) reverse_mapping(mapped_vars, all=all) end diff --git a/src/processes/model_initialisation.jl b/src/processes/model_initialisation.jl index 95a8cc76..04702072 100755 --- a/src/processes/model_initialisation.jl +++ b/src/processes/model_initialisation.jl @@ -148,7 +148,7 @@ function to_initialize(mapping::Dict{String,T}, graph=nothing) where {T} end to_init = Dict(org => Symbol[] for org in keys(mapping)) - mapped_vars = mapped_variables(mapping, first(hard_dependencies(mapping; verbose=false)), verbose=false) + mapped_vars = mapped_variables(mapping, first(hard_dependencies(mapping; verbose=false, orchestrator=Orchestrator2())), verbose=false) for (org, vars) in mapped_vars for (var, val) in vars if isa(val, UninitializedVar) && var ∉ vars_in_mtg diff --git a/src/run.jl b/src/run.jl index 03974654..8a62a200 100644 --- a/src/run.jl +++ b/src/run.jl @@ -365,21 +365,18 @@ function run!( meteo=nothing, constants=PlantMeteo.Constants(), extra=nothing; - orchestrator::Orchestrator=nothing, nsteps=nothing, tracked_outputs=nothing, check=true, - executor=ThreadedEx(), - default_timestep::Int, - model_timesteps::Dict{T, Int} where {T} - + executor=ThreadedEx(), + orchestrator::Orchestrator2=Orchestrator2(), ) isnothing(nsteps) && (nsteps = get_nsteps(meteo)) meteo_adjusted = adjust_weather_timesteps_to_given_length(nsteps, meteo) # NOTE : replace_mapping_status_vectors_with_generated_models is assumed to have already run if used # otherwise there might be vector length conflicts with timesteps - sim = GraphSimulation(object, mapping, nsteps=nsteps, check=check, outputs=tracked_outputs, default_timestep=default_timestep, model_timesteps=model_timesteps) + sim = GraphSimulation(object, mapping, nsteps=nsteps, check=check, outputs=tracked_outputs, orchestrator=orchestrator) run!( sim, meteo_adjusted, @@ -452,6 +449,15 @@ function run_node_multiscale!( executor ) where {T<:GraphSimulation} # T is the status of each node by organ type + # Hack until default timestep is inferred from the meteo + node_ratio = node.timestep / Day(1) + + # Check if the model needs to run this timestep + if (1 + (i-1) % node_ratio) != node_ratio + node.simulation_id[1] += 1 + return nothing + end + # run!(status(object), dependency_node, meteo, constants, extra) # Check if all the parents have been called before the child: if !AbstractTrees.isroot(node) && any([p.simulation_id[1] <= node.simulation_id[1] for p in node.parent]) @@ -459,33 +465,60 @@ function run_node_multiscale!( return nothing end - model_timestep = object.model_timesteps[typeof(node.value)] - - if model_timestep != object.default_timestep - # do accumulation - - - # run if necessary - if i % model_timestep != 0 - return nothing - end - end - node_statuses = status(object)[node.scale] # Get the status of the nodes at the current scale models_at_scale = models[node.scale] - for st in node_statuses # for each node status at the current scale (potentially in parallel over nodes) - # Actual call to the model: - run!(node.value, models_at_scale, st, meteo, constants, extra) - end + if isnothing(node.timestep_mapping_data) || isempty(node.timestep_mapping_data) + # Samuel : this is the happy path, no further timestep mapping checks needed - node.simulation_id[1] += 1 # increment the simulation id, to remember that the model has been called already + for st in node_statuses # for each node status at the current scale (potentially in parallel over nodes) + # Actual call to the model: + run!(node.value, models_at_scale, st, meteo, constants, extra) + end + node.simulation_id[1] += 1 # increment the simulation id, to remember that the model has been called already - # Recursively visit the children (soft dependencies only, hard dependencies are handled by the model itself): - for child in node.children + # Recursively visit the children (soft dependencies only, hard dependencies are handled by the model itself): + for child in node.children #! check if we can run this safely in a @floop loop. I would say no, #! because we are running a parallel computation above already, modifying the node.simulation_id, #! which is not thread-safe yet. - run_node_multiscale!(object, child, i, models, meteo, constants, extra, check, executor) + run_node_multiscale!(object, child, i, models, meteo, constants, extra, check, executor) + end + else + # we have a child with a different timestep than the parent, + # which requires accumulation/reduce and running only some of the time + + for st in node_statuses + run!(node.value, models_at_scale, st, meteo, constants, extra) + # TODO do the accumulation + # then write into the child's status if need be ? + for tmst in node.timestep_mapping_data + ratio = Int(tmst.node_to.timestep / node.timestep) + # TODO assert etc. + # do the accumulation for each variable + accumulation_index = 1 + ((i-1)%ratio) + tmst.mapping_data[node_id(st.node)][accumulation_index] = st[tmst.variable_from] + + # if we have completed a *full* cycle, then presumably we need to write out the value to + # the mapped model + # A full cycle isn't just the ratio to the parent, it's the ratio to the finest-grained timestep + if accumulation_index == ratio + node_statuses_to = status(object)[tmst.node_to.scale] + + # TODO : INCORRECT in a scale with multiple mtg nodes + for st_to in node_statuses_to + # TODO might be able to catch mapping_function type incompatibility errors and make them clearer + st_to[tmst.variable_to] = tmst.mapping_function(tmst.mapping_data[node_id(st.node)]) + end + end + end + end + + node.simulation_id[1] += 1 + + for child in node.children + + run_node_multiscale!(object, child, i, models, meteo, constants, extra, check, executor) + end end end \ No newline at end of file diff --git a/src/timestep/timestep_mapping.jl b/src/timestep/timestep_mapping.jl index e30fc97c..20ad038e 100644 --- a/src/timestep/timestep_mapping.jl +++ b/src/timestep/timestep_mapping.jl @@ -4,39 +4,116 @@ # Some types can also be constrained a lot more, probably +# Many shortcuts will be taken, I'll try and comment what's missing/implicit/etc. + +# TODO specify scale ? struct TimestepMapper#{V} variable_from#::V - timestep_from::Int + timestep_from::Period # ? Not sure whether this is the best bit of info... Also, to or from ? And it should be a Period, no ? mapping_function + mapping_data # TODO this should be internal end struct SimulationTimestepHandler#{W,V} - model_timesteps::Dict{Any, Int} # where {W <: AbstractModel} # if a model isn't in there, then it follows the default, todo check if the given timestep respects the model's range + model_timesteps::Dict{Any, Period} # where {W <: AbstractModel} # if a model isn't in there, then it follows the default, todo check if the given timestep respects the model's range timestep_variable_mapping::Dict{Any, TimestepMapper} #where {V} end -SimulationTimestepHandler() = SimulationTimestepHandler(Dict(), Dict()) #Dict{W, Int}(), Dict{V, TimestepMapper}()) where {W, V} +SimulationTimestepHandler() = SimulationTimestepHandler(Dict{Any, Int}(), Dict{Any, TimestepMapper}()) #Dict{W, Int}(), Dict{V, TimestepMapper}()) where {W, V} mutable struct Orchestrator # This is actually a general simulation parameter, not-scale specific # todo change to Period - default_timestep::Int64 + default_timestep::Period # This needs to be per-scale : if a model is used at two different scales, # and the same variable of that model maps to some other timestep to two *different* variables # then I believe we can only rely on the different scale to disambiguate non_default_timestep_data_per_scale::Dict{String, SimulationTimestepHandler} - function Orchestrator(default::Int64, per_scale::Dict{String, SimulationTimestepHandler}) - @assert default >= 0 "The default_timestep should be greater than or equal to 0." + function Orchestrator(default::Period, per_scale::Dict{String, SimulationTimestepHandler}) + @assert default >= Second(0) "The default_timestep should be greater than or equal to 0." return new(default, per_scale) end end # TODO have a default constructor take in a meteo or something, and set up the default timestep automagically to be the finest weather timestep # Other options are possible -Orchestrator() = Orchestrator(1, Dict{String, SimulationTimestepHandler}()) +Orchestrator() = Orchestrator(Day(1), Dict{String, SimulationTimestepHandler}()) #o = Orchestrator() -#oo = Orchestrator(1, Dict{String, SimulationTimestepHandler}()) \ No newline at end of file +#oo = Orchestrator(1, Dict{String, SimulationTimestepHandler}()) + + +# TODO issue : what if the user wants to force initialize a variable that isn't at the default timestep ? +# This requires the ability to fit data with a vector that isn't at the default timestep +# Automated model generation does not have that feature +# As is, the current workaround is for the user to write their own model, I think, which is not ideal + +# TODO check for cycles (and other errors) before timestep mapping, then do it again afterwards, as the new mapping dependencies might cause specific cycles. + + +# TODO status initialisation ? +# TODO type promotion for mapped timestep variables ? +# TODO check type if a variable is timestep mapped and scale mapped +# TODO simulation_id change consequences ? + + + +# TODO prev timestep ? Vector mapping ? +struct Var_from + model + scale::String + name::Symbol + mapping_function::Function + # mapping_data::Dict{NodeMTG, Vector{stuff}} +end + +struct Var_to + name::Symbol +end + +struct ModelTimestepMapping + model + scale::String + timestep::Period + var_to_var::Dict{Var_to, Var_from} +end + +mutable struct Orchestrator2 + default_timestep::Period + non_default_timestep_mapping::Vector{ModelTimestepMapping} + + function Orchestrator2(default::Period, non_default_timestep_mapping::Vector{ModelTimestepMapping}) + @assert default >= Second(0) "The default_timestep should be greater than or equal to 0." + return new(default, non_default_timestep_mapping) + end +end + +Orchestrator2() = Orchestrator2(Day(1), Vector{ModelTimestepMapping}()) + + +# TODO parallelization, + +function init_timestep_mapping_data(node_mtg::MultiScaleTreeGraph.Node, dependency_graph) + traverse_dependency_graph!(x -> register_mtg_node_in_timestep_mapping(x, node_mtg), dependency_graph, visit_hard_dep=false) +end + +function register_mtg_node_in_timestep_mapping(node_dep::SoftDependencyNode, node_mtg::MultiScaleTreeGraph.Node) + if isnothing(node_dep.timestep_mapping_data) + return + end + + # no need to check the current softdependencynode's scale, I think + # only the mapped downstream softdependencynodes + + # TODO this structure doesn't play well with parallelisation... ? + # TODO having an extra level of indirection, mapping the MTG node to an index into a vector + # Allows one to resize! the vector when it lacks space, saving in terms of # of memory allocations/copies + for mtsm in node_dep.timestep_mapping_data + if node_dep.scale == symbol(node_mtg) + push!(mtsm.mapping_data, node_id(node_mtg) => deepcopy(mtsm.mapping_data_template)) + end + end +end \ No newline at end of file diff --git a/test/test_multitimestep.jl b/test/test_multitimestep.jl new file mode 100644 index 00000000..9e25f1ce --- /dev/null +++ b/test/test_multitimestep.jl @@ -0,0 +1,457 @@ + +########################### +# Simple test using an ad hoc connector model +# Broken by subsequent changes, left just in case for now (TODO remove once prototyping is over) +########################### + +#= +using PlantSimEngine +# Include the example dummy processes: +using PlantSimEngine.Examples +using Test, Aqua +using Tables, DataFrames, CSV +using MultiScaleTreeGraph +using PlantMeteo, Statistics +using Documenter # for doctests + +using PlantMeteo.Dates +include("helper-functions.jl") + + + +# These models might be worth exposing in the future ? +PlantSimEngine.@process "basic_current_timestep" verbose = false + +struct HelperCurrentTimestepModel <: AbstractBasic_Current_TimestepModel +end + +PlantSimEngine.inputs_(::HelperCurrentTimestepModel) = (next_timestep=1,) +PlantSimEngine.outputs_(m::HelperCurrentTimestepModel) = (current_timestep=1,) + +function PlantSimEngine.run!(m::HelperCurrentTimestepModel, models, status, meteo, constants=nothing, extra=nothing) + status.current_timestep = status.next_timestep + end + + PlantSimEngine.ObjectDependencyTrait(::Type{<:HelperCurrentTimestepModel}) = PlantSimEngine.IsObjectDependent() + PlantSimEngine.TimeStepDependencyTrait(::Type{<:HelperCurrentTimestepModel}) = PlantSimEngine.IsTimeStepDependent() + +PlantSimEngine.timestep_range_(m::HelperCurrentTimestepModel) = Day(1) + + + PlantSimEngine.@process "basic_next_timestep" verbose = false + struct HelperNextTimestepModel <: AbstractBasic_Next_TimestepModel + end + + PlantSimEngine.inputs_(::HelperNextTimestepModel) = (current_timestep=1,) + PlantSimEngine.outputs_(m::HelperNextTimestepModel) = (next_timestep=1,) + + function PlantSimEngine.run!(m::HelperNextTimestepModel, models, status, meteo, constants=nothing, extra=nothing) + status.next_timestep = status.current_timestep + 1 + end + +PlantSimEngine.timestep_range_(m::HelperNextTimestepModel) = Day(1) + + + + + +PlantSimEngine.@process "ToyDay" verbose = false + +struct MyToyDayModel <: AbstractToydayModel end + +PlantSimEngine.inputs_(m::MyToyDayModel) = (a=1,) +PlantSimEngine.outputs_(m::MyToyDayModel) = (daily_temperature=-Inf,) + +function PlantSimEngine.run!(m::MyToyDayModel, models, status, meteo, constants=nothing, extra=nothing) + status.daily_temperature = meteo.T +end + +PlantSimEngine.@process "ToyWeek" verbose = false + +struct MyToyWeekModel <: AbstractToyweekModel + temperature_threshold::Float64 +end + +MyToyWeekModel() = MyToyWeekModel(30.0) +function PlantSimEngine.inputs_(::MyToyWeekModel) + (weekly_max_temperature=-Inf,) +end +PlantSimEngine.outputs_(m::MyToyWeekModel) = (hot = false,) + +function PlantSimEngine.run!(m::MyToyWeekModel, models, status, meteo, constants=nothing, extra=nothing) + status.hot = status.weekly_max_temperature > m.temperature_threshold +end + +PlantSimEngine.timestep_range_(m::MyToyWeekModel) = Week(1) + + + +PlantSimEngine.@process "DWConnector" verbose = false + +struct MyDwconnectorModel <: AbstractDwconnectorModel + T_daily::Array{Float64} +end + +MyDwconnectorModel() = MyDwconnectorModel(Array{Float64}(undef, 7)) + +function PlantSimEngine.inputs_(::MyDwconnectorModel) + (daily_temperature=-Inf, current_timestep=1,) +end +PlantSimEngine.outputs_(m::MyDwconnectorModel) = (weekly_max_temperature = 0.0,) + +function PlantSimEngine.run!(m::MyDwconnectorModel, models, status, meteo, constants=nothing, extra=nothing) + m.T_daily[1 + (status.current_timestep % 7)] = status.daily_temperature + + if(status.current_timestep % 7 == 1) + status.weekly_max_temperature = sum(m.T_daily)/7.0 + else + status.weekly_max_temperature = 0 + end +end + + PlantSimEngine.timestep_range_(m::MyDwconnectorModel) = Day(1) + + + + + +meteo_day = read_weather(joinpath(pkgdir(PlantSimEngine), "examples/meteo_day.csv"), duration=Day) + +m = Dict("Default" => ( + MyToyDayModel(), + MyToyWeekModel(), + MyDwconnectorModel(), + HelperNextTimestepModel(), + MultiScaleModel( + model=HelperCurrentTimestepModel(), + mapped_variables=[PreviousTimeStep(:next_timestep),], + ), + Status(a=1,))) + +to_initialize(m) + +models_timestep = Dict(MyToyDayModel=>1, MyDwconnectorModel => 1, MyToyWeekModel =>7, HelperNextTimestepModel => 1, HelperCurrentTimestepModel => 1) + +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) + +out = run!(mtg, m, meteo_day, default_timestep=1, model_timesteps=models_timestep) + + +# NOTE : replace_mapping_status_vectors_with_generated_models is assumed to have already run if used +# otherwise there might be vector length conflicts with timesteps +sim = PlantSimEngine.GraphSimulation(mtg, m, nsteps=nothing, check=true, outputs=nothing, default_timestep=1, model_timesteps=models_timestep) + +=# + + +########################### +# First attempt at an orchetrator, broken by subsequent changes +########################### +#= +using MultiScaleTreeGraph +using PlantSimEngine +using PlantMeteo +using PlantMeteo.Dates + +PlantSimEngine.@process "ToyDay" verbose = false + +struct MyToyDayModel <: AbstractToydayModel end + +PlantSimEngine.inputs_(m::MyToyDayModel) = (a=1,) +PlantSimEngine.outputs_(m::MyToyDayModel) = (daily_temperature=-Inf,) + +function PlantSimEngine.run!(m::MyToyDayModel, models, status, meteo, constants=nothing, extra=nothing) + status.daily_temperature = meteo.T +end + +PlantSimEngine.@process "ToyWeek" verbose = false + +struct MyToyWeekModel <: AbstractToyweekModel + temperature_threshold::Float64 +end + +MyToyWeekModel() = MyToyWeekModel(30.0) +function PlantSimEngine.inputs_(::MyToyWeekModel) + (weekly_max_temperature=-Inf,) +end +PlantSimEngine.outputs_(m::MyToyWeekModel) = (hot = false,) + +function PlantSimEngine.run!(m::MyToyWeekModel, models, status, meteo, constants=nothing, extra=nothing) + status.hot = status.weekly_max_temperature > m.temperature_threshold +end + +PlantSimEngine.timestep_range_(m::MyToyWeekModel) = TimestepRange(Week(1)) + + +m = Dict("Default" => ( + MyToyDayModel(), + MyToyWeekModel(), + Status(a=1,))) + + meteo_day = read_weather(joinpath(pkgdir(PlantSimEngine), "examples/meteo_day.csv"), duration=Day) +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) + + + model_timesteps_defaultscale = Dict(MyToyWeekModel =>Week(1)) + tsm = PlantSimEngine.TimestepMapper(:daily_temperature, Day(1), max, nothing) + sth = PlantSimEngine.SimulationTimestepHandler(model_timesteps_defaultscale, Dict(:weekly_max_temperature => tsm )) + + orchestrator = Orchestrator(Day(1), Dict("Default" => sth)) + + #out = @enter run!(mtg, m, meteo_day, orchestrator=orchestrator) + +# TODO could some mapping happen automatically for variables directly taken from weather data ? +# Does this happen often in a typical model ? + +#=m_multiscale = Dict("Default" => ( + MyToyDayModel(), + Status(a=1,) + ), + "Default2" => ( + MyToyWeekModel(), + ),) +=# +m_multiscale = Dict("Default" => ( + MyToyDayModel(), + Status(a=1,) + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeekModel(), + mapped_variables=[:weekly_max_temperature => "Default" => :daily_temperature], + ), + ),) + + + +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) +mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 2)) + + orchestrator_multiscale = Orchestrator(Day(1), Dict("Default2" => sth)) + +#out = @enter run!(mtg, m_multiscale, meteo_day, orchestrator=orchestrator_multiscale) +=# + + +########################### +# Simple test with a second attempt at an orchestrator +# Functional, except for the revalue/refvector overwriting issues +########################### +using MultiScaleTreeGraph +using PlantSimEngine +using PlantMeteo +using PlantMeteo.Dates + +PlantSimEngine.@process "ToyDay" verbose = false + +struct MyToyDayModel <: AbstractToydayModel end + +PlantSimEngine.inputs_(m::MyToyDayModel) = (a=1,) +PlantSimEngine.outputs_(m::MyToyDayModel) = (daily_temperature=-Inf,) + +function PlantSimEngine.run!(m::MyToyDayModel, models, status, meteo, constants=nothing, extra=nothing) + status.daily_temperature = meteo.T +end + +PlantSimEngine.@process "ToyWeek" verbose = false + +struct MyToyWeekModel <: AbstractToyweekModel + temperature_threshold::Float64 +end + +MyToyWeekModel() = MyToyWeekModel(30.0) +function PlantSimEngine.inputs_(::MyToyWeekModel) + (weekly_max_temperature=-Inf,) +end +PlantSimEngine.outputs_(m::MyToyWeekModel) = (hot = false,) + +function PlantSimEngine.run!(m::MyToyWeekModel, models, status, meteo, constants=nothing, extra=nothing) + status.hot = status.weekly_max_temperature > m.temperature_threshold +end + +PlantSimEngine.timestep_range_(m::MyToyWeekModel) = TimestepRange(Week(1)) + + + meteo_day = read_weather(joinpath(pkgdir(PlantSimEngine), "examples/meteo_day.csv"), duration=Day) + + model_timesteps_defaultscale = Dict(MyToyWeekModel =>Week(1)) + tsm = PlantSimEngine.TimestepMapper(:daily_temperature, Day(1), max, nothing) + sth = PlantSimEngine.SimulationTimestepHandler(model_timesteps_defaultscale, Dict(:weekly_max_temperature => tsm )) + +m_multiscale = Dict("Default" => ( + MyToyDayModel(), + Status(a=1,) + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeekModel(), + #mapped_variables=[:weekly_max_temperature => ["Default" => :daily_temperature]], # TODO test this + mapped_variables=[:weekly_max_temperature => "Default" => :daily_temperature], + ), + ),) + + + +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) +mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 2)) + + +to = PlantSimEngine.Var_to(:weekly_max_temperature) +from = PlantSimEngine.Var_from(MyToyDayModel, "Default", :daily_temperature, maximum) + +dict_to_from = Dict(to => from) +mtsm = PlantSimEngine.ModelTimestepMapping(MyToyWeekModel, "Default2", Week(1), dict_to_from) + +orch2 = PlantSimEngine.Orchestrator2(Day(1), [mtsm,]) + +out = @enter run!(mtg, m_multiscale, meteo_day, orchestrator=orch2) + + +########################### +# Test with three timesteps, multiscale +# Issues with data overwriting (refvector/refvalue) +########################### + +using MultiScaleTreeGraph +using PlantSimEngine +using PlantMeteo +using PlantMeteo.Dates + +PlantSimEngine.@process "ToyDay2" verbose = false + +struct MyToyDay2Model <: AbstractToyday2Model end + +PlantSimEngine.inputs_(m::MyToyDay2Model) = NamedTuple() +PlantSimEngine.outputs_(m::MyToyDay2Model) = (out_day=-Inf,) + +function PlantSimEngine.run!(m::MyToyDay2Model, models, status, meteo, constants=nothing, extra=nothing) + status.out_day = meteo.data +end + +PlantSimEngine.@process "ToyWeek2" verbose = false + +struct MyToyWeek2Model <: AbstractToyweek2Model end + +PlantSimEngine.inputs_(::MyToyWeek2Model) = (in_week=-Inf,) +PlantSimEngine.outputs_(m::MyToyWeek2Model) = (out_week=-Inf,) + +function PlantSimEngine.run!(m::MyToyWeek2Model, models, status, meteo, constants=nothing, extra=nothing) + status.out_week = status.in_week +end + +PlantSimEngine.timestep_range_(m::MyToyWeek2Model) = TimestepRange(Week(1)) + + +PlantSimEngine.@process "ToyFourWeek2" verbose = false + +struct MyToyFourWeek2Model <: AbstractToyfourweek2Model end + +PlantSimEngine.inputs_(::MyToyFourWeek2Model) = (in_four_week_from_week=-Inf, in_four_week_from_day=-Inf,) +PlantSimEngine.outputs_(m::MyToyFourWeek2Model) = (inputs_agreement=false,) + +function PlantSimEngine.run!(m::MyToyFourWeek2Model, models, status, meteo, constants=nothing, extra=nothing) + status.inputs_agreement = status.in_four_week_from_week == status.in_four_week_from_day +end + +PlantSimEngine.timestep_range_(m::MyToyFourWeek2Model) = TimestepRange(Week(4)) + + + +df = DataFrame(:data => [1 for i in 1:365], ) + + # TODO can make this optional if the timestep range is actually a single value + model_timesteps_defaultscale = Dict(MyToyWeek2Model =>Week(1), MyToyFourWeek2Model =>Week(4), ) + tsm_d = PlantSimEngine.TimestepMapper(:out_day, Day(1), sum, nothing) + sth_d = PlantSimEngine.SimulationTimestepHandler(model_timesteps_defaultscale, Dict(:in_week => tsm_d )) + + tsm_w = PlantSimEngine.TimestepMapper(:out_week, Week(1), sum, nothing) + sth_w = PlantSimEngine.SimulationTimestepHandler(model_timesteps_defaultscale, Dict(:in_four_week_from_day => tsm_d, :in_four_week_from_week => tsm_w )) + + to_w = PlantSimEngine.Var_to(:in_week) + from_d = PlantSimEngine.Var_from(MyToyDay2Model, "Default", :out_day, sum) + dict_to_from_w = Dict(to_w => from_d) + + to_w4_d = PlantSimEngine.Var_to(:in_four_week_from_day) + to_w4_w = PlantSimEngine.Var_to(:in_four_week_from_week) + from_w = PlantSimEngine.Var_from(MyToyWeek2Model, "Default2", :out_week, sum) + + dict_to_from_w4 = Dict(to_w4_d => from_d, to_w4_w => from_w) + + mtsm_w = PlantSimEngine.ModelTimestepMapping(MyToyWeek2Model, "Default2", Week(1), dict_to_from_w) + mtsm_w4 = PlantSimEngine.ModelTimestepMapping(MyToyFourWeek2Model, "Default3", Week(4), dict_to_from_w4) + + orch2 = PlantSimEngine.Orchestrator2(Day(1), [mtsm_w, mtsm_w4]) + + +m_multiscale = Dict("Default" => ( + MyToyDay2Model(), + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeek2Model(), + mapped_variables=[:in_week => "Default" => :out_day], + ), + ), + "Default3" => ( + MultiScaleModel(model=MyToyFourWeek2Model(), + mapped_variables=[ + :in_four_week_from_day => "Default" => :out_day, + :in_four_week_from_week => "Default2" => :out_week, + ], + ), + ),) + + +# TODO test with multiple nodes +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) +mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 2)) +mtg3 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default3", 1, 3)) +#out = @enter run!(mtg, m_multiscale, df, orchestrator=orch2) + + + +########################### +# Three timestep model that is single-scale, to circumvent refvector/refvalue overwriting +# and explore alternatives +# (eg filtering out timestep-mapped variables from vars_need_init and storing the values elsewhere) +########################### + + m_singlescale = Dict("Default" => ( + MyToyDay2Model(), + MyToyWeek2Model(), + MyToyFourWeek2Model(), + ),) + + + model_timesteps_defaultscale = Dict(MyToyWeek2Model =>Week(1), MyToyFourWeek2Model =>Week(4), ) + tsm_d = PlantSimEngine.TimestepMapper(:out_day, Day(1), sum, nothing) + sth_d = PlantSimEngine.SimulationTimestepHandler(model_timesteps_defaultscale, Dict(:in_week => tsm_d )) + + tsm_w = PlantSimEngine.TimestepMapper(:out_week, Week(1), sum, nothing) + sth_w = PlantSimEngine.SimulationTimestepHandler(model_timesteps_defaultscale, Dict(:in_four_week_from_day => tsm_d, :in_four_week_from_week => tsm_w )) + + to_w = PlantSimEngine.Var_to(:in_week) + from_d = PlantSimEngine.Var_from(MyToyDay2Model, "Default", :out_day, sum) + dict_to_from_w = Dict(to_w => from_d) + + to_w4_d = PlantSimEngine.Var_to(:in_four_week_from_day) + to_w4_w = PlantSimEngine.Var_to(:in_four_week_from_week) + from_w = PlantSimEngine.Var_from(MyToyWeek2Model, "Default", :out_week, sum) + + dict_to_from_w4 = Dict(to_w4_d => from_d, to_w4_w => from_w) + + mtsm_w = PlantSimEngine.ModelTimestepMapping(MyToyWeek2Model, "Default", Week(1), dict_to_from_w) + mtsm_w4 = PlantSimEngine.ModelTimestepMapping(MyToyFourWeek2Model, "Default", Week(4), dict_to_from_w4) + + orch2 = PlantSimEngine.Orchestrator2(Day(1), [mtsm_w, mtsm_w4]) + +mtg_single = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) +out = @run run!(mtg_single, m_singlescale, df, orchestrator=orch2) + + + +using Test + @test unique([out["Default3"][i].in_four_week_from_day for i in 1:length(out["Default3"])]) == [0.0, 28.0] + @test unique([out["Default3"][i].in_four_week_from_week for i in 1:length(out["Default3"])]) == [0.0, 7.0] + @test unique([out["Default3"][i].inputs_agreement for i in 1:length(out["Default3"])]) == [1] + + unique([out["Default3"][i].in_four_week_from_week for i in 1:length(out["Default3"])]) + unique([out["Default3"][i].in_four_week_from_day for i in 1:length(out["Default3"])]) + unique([out["Default3"][i].inputs_agreement for i in 1:length(out["Default3"])]) \ No newline at end of file From 2d343f666a99fc04bebe225fac861fb03a448bbe Mon Sep 17 00:00:00 2001 From: Samuel-AMAP Date: Tue, 28 Oct 2025 15:13:02 +0100 Subject: [PATCH 4/8] Addendum : prototype tests moved to a new file, remove them from the old file. --- test/test-simulation.jl | 151 +--------------------------------------- 1 file changed, 1 insertion(+), 150 deletions(-) diff --git a/test/test-simulation.jl b/test/test-simulation.jl index c4439308..4ff83a64 100644 --- a/test/test-simulation.jl +++ b/test/test-simulation.jl @@ -292,153 +292,4 @@ end; end end end -end - - - - -using PlantSimEngine -# Include the example dummy processes: -using PlantSimEngine.Examples -using Test, Aqua -using Tables, DataFrames, CSV -using MultiScaleTreeGraph -using PlantMeteo, Statistics -using Documenter # for doctests - -using PlantMeteo.Dates -include("helper-functions.jl") - - - -# These models might be worth exposing in the future ? -PlantSimEngine.@process "basic_current_timestep" verbose = false - -struct HelperCurrentTimestepModel <: AbstractBasic_Current_TimestepModel -end - -PlantSimEngine.inputs_(::HelperCurrentTimestepModel) = (next_timestep=1,) -PlantSimEngine.outputs_(m::HelperCurrentTimestepModel) = (current_timestep=1,) - -function PlantSimEngine.run!(m::HelperCurrentTimestepModel, models, status, meteo, constants=nothing, extra=nothing) - status.current_timestep = status.next_timestep - end - - PlantSimEngine.ObjectDependencyTrait(::Type{<:HelperCurrentTimestepModel}) = PlantSimEngine.IsObjectDependent() - PlantSimEngine.TimeStepDependencyTrait(::Type{<:HelperCurrentTimestepModel}) = PlantSimEngine.IsTimeStepDependent() - -PlantSimEngine.timestep_range_(m::HelperCurrentTimestepModel) = Day(1) - - - PlantSimEngine.@process "basic_next_timestep" verbose = false - struct HelperNextTimestepModel <: AbstractBasic_Next_TimestepModel - end - - PlantSimEngine.inputs_(::HelperNextTimestepModel) = (current_timestep=1,) - PlantSimEngine.outputs_(m::HelperNextTimestepModel) = (next_timestep=1,) - - function PlantSimEngine.run!(m::HelperNextTimestepModel, models, status, meteo, constants=nothing, extra=nothing) - status.next_timestep = status.current_timestep + 1 - end - -PlantSimEngine.timestep_range_(m::HelperNextTimestepModel) = Day(1) - - - - - -PlantSimEngine.@process "ToyDay" verbose = false - -struct MyToyDayModel <: AbstractToydayModel end - -PlantSimEngine.inputs_(m::MyToyDayModel) = (a=1,) -PlantSimEngine.outputs_(m::MyToyDayModel) = (daily_temperature=-Inf,) - -function PlantSimEngine.run!(m::MyToyDayModel, models, status, meteo, constants=nothing, extra=nothing) - status.daily_temperature = meteo.T -end - -PlantSimEngine.@process "ToyWeek" verbose = false - -struct MyToyWeekModel <: AbstractToyweekModel - temperature_threshold::Float64 -end - -MyToyWeekModel() = MyToyWeekModel(30.0) -function PlantSimEngine.inputs_(::MyToyWeekModel) - (weekly_max_temperature=-Inf,) -end -PlantSimEngine.outputs_(m::MyToyWeekModel) = (hot = false,) - -function PlantSimEngine.run!(m::MyToyWeekModel, models, status, meteo, constants=nothing, extra=nothing) - status.hot = status.weekly_max_temperature > m.temperature_threshold -end - -PlantSimEngine.timestep_range_(m::MyToyWeekModel) = Week(1) - - - -PlantSimEngine.@process "DWConnector" verbose = false - -struct MyDwconnectorModel <: AbstractDwconnectorModel - T_daily::Array{Float64} -end - -MyDwconnectorModel() = MyDwconnectorModel(Array{Float64}(undef, 7)) - -function PlantSimEngine.inputs_(::MyDwconnectorModel) - (daily_temperature=-Inf, current_timestep=1,) -end -PlantSimEngine.outputs_(m::MyDwconnectorModel) = (weekly_max_temperature = 0.0,) - -function PlantSimEngine.run!(m::MyDwconnectorModel, models, status, meteo, constants=nothing, extra=nothing) - m.T_daily[1 + (status.current_timestep % 7)] = status.daily_temperature - - if(status.current_timestep % 7 == 1) - status.weekly_max_temperature = sum(m.T_daily)/7.0 - else - status.weekly_max_temperature = 0 - end -end - - PlantSimEngine.timestep_range_(m::MyDwconnectorModel) = Day(1) - - - - - -meteo_day = read_weather(joinpath(pkgdir(PlantSimEngine), "examples/meteo_day.csv"), duration=Day) - -m = Dict("Default" => ( - MyToyDayModel(), - MyToyWeekModel(), - MyDwconnectorModel(), - HelperNextTimestepModel(), - MultiScaleModel( - model=HelperCurrentTimestepModel(), - mapped_variables=[PreviousTimeStep(:next_timestep),], - ), - Status(a=1,))) - -to_initialize(m) - -models_timestep = Dict(MyToyDayModel=>1, MyDwconnectorModel => 1, MyToyWeekModel =>7, HelperNextTimestepModel => 1, HelperCurrentTimestepModel => 1) - -mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) - -out = run!(mtg, m, meteo_day, default_timestep=1, model_timesteps=models_timestep) - -@testset "Test varying timestep" begin - - - @test - @test - -end - - - # NOTE : replace_mapping_status_vectors_with_generated_models is assumed to have already run if used - # otherwise there might be vector length conflicts with timesteps - sim = @enter PlantSimEngine.GraphSimulation(mtg, m, nsteps=nothing, check=true, outputs=nothing, default_timestep=1, model_timesteps=models_timestep) - -using PlantSimEngine +end \ No newline at end of file From da65ff7399e331e6834d80b4a5765d394129f765 Mon Sep 17 00:00:00 2001 From: Samuel-AMAP Date: Fri, 31 Oct 2025 14:44:48 +0100 Subject: [PATCH 5/8] Fix little issues + attempt at addressing main blocker : disconnect timestep-mapped variables from their source, to avoid overwriting source (inputs to non-default-timestep models are changed by the accumulation function so can't be a simple Ref to source) --- src/dependencies/dependency_graph.jl | 78 -------------- src/dependencies/hard_dependencies.jl | 88 +++++++++++++++- src/dependencies/soft_dependencies.jl | 2 +- src/mtg/initialisation.jl | 50 ++++++--- src/mtg/mapping/compute_mapping.jl | 20 +++- src/mtg/save_results.jl | 10 +- src/run.jl | 7 +- src/timestep/timestep_mapping.jl | 2 +- test/test_multitimestep.jl | 142 ++++++++++++++++++-------- 9 files changed, 249 insertions(+), 150 deletions(-) diff --git a/src/dependencies/dependency_graph.jl b/src/dependencies/dependency_graph.jl index 33b5fc62..59ec29f5 100644 --- a/src/dependencies/dependency_graph.jl +++ b/src/dependencies/dependency_graph.jl @@ -89,82 +89,4 @@ function Base.show(io::IO, ::MIME"text/plain", t::DependencyGraph) else draw_dependency_graph(io, t) end -end - -""" - variables_multiscale(node, organ, mapping, st=NamedTuple()) - -Get the variables of a HardDependencyNode, taking into account the multiscale mapping, *i.e.* -defining variables as `MappedVar` if they are mapped to another scale. The default values are -taken from the model if not given by the user (`st`), and are marked as `UninitializedVar` if -they are inputs of the node. - -Return a NamedTuple with the variables and their default values. - -# Arguments - -- `node::HardDependencyNode`: the node to get the variables from. -- `organ::String`: the organ type, *e.g.* "Leaf". -- `vars_mapping::Dict{String,T}`: the mapping of the models (see details below). -- `st::NamedTuple`: an optional named tuple with default values for the variables. - -# Details - -The `vars_mapping` is a dictionary with the organ type as key and a dictionary as value. It is -computed from the user mapping like so: -""" -function variables_multiscale(node, organ, vars_mapping, st=NamedTuple()) - node_vars = variables(node) # e.g. (inputs = (:var1=-Inf, :var2=-Inf), outputs = (:var3=-Inf,)) - ins = node_vars.inputs - ins_variables = keys(ins) - outs_variables = keys(node_vars.outputs) - defaults = merge(node_vars...) - map((inputs=ins_variables, outputs=outs_variables)) do vars # Map over vars from :inputs and vars from :outputs - vars_ = Vector{Pair{Symbol,Any}}() - for var in vars # e.g. var = :carbon_biomass - if var in keys(st) - #If the user has given a status, we use it as default value. - default = st[var] - elseif var in ins_variables - # Otherwise, we use the default value given by the model: - # If the variable is an input, we mark it as uninitialized: - default = UninitializedVar(var, defaults[var]) - else - # If the variable is an output, we use the default value given by the model: - default = defaults[var] - end - - if haskey(vars_mapping[organ], var) - organ_mapped, organ_mapped_var = _node_mapping(vars_mapping[organ][var]) - push!(vars_, var => MappedVar(organ_mapped, var, organ_mapped_var, default)) - #* We still check if the variable also exists wrapped in PreviousTimeStep, because one model could use the current - #* values, and another one the previous values. - if haskey(vars_mapping[organ], PreviousTimeStep(var, node.process)) - organ_mapped, organ_mapped_var = _node_mapping(vars_mapping[organ][PreviousTimeStep(var, node.process)]) - push!(vars_, var => MappedVar(organ_mapped, PreviousTimeStep(var, node.process), organ_mapped_var, default)) - end - elseif haskey(vars_mapping[organ], PreviousTimeStep(var, node.process)) - # If not found in the current time step, we check if the variable is mapped to the previous time step: - organ_mapped, organ_mapped_var = _node_mapping(vars_mapping[organ][PreviousTimeStep(var, node.process)]) - push!(vars_, var => MappedVar(organ_mapped, PreviousTimeStep(var, node.process), organ_mapped_var, default)) - else - # Else we take the default value: - push!(vars_, var => default) - end - end - return (; vars_...,) - end -end - -function _node_mapping(var_mapping::Pair{String,Symbol}) - # One organ is mapped to the variable: - return SingleNodeMapping(first(var_mapping)), last(var_mapping) -end - -function _node_mapping(var_mapping) - # Several organs are mapped to the variable: - organ_mapped = MultiNodeMapping([first(i) for i in var_mapping]) - organ_mapped_var = [last(i) for i in var_mapping] - - return organ_mapped, organ_mapped_var end \ No newline at end of file diff --git a/src/dependencies/hard_dependencies.jl b/src/dependencies/hard_dependencies.jl index eda0573f..04fb16f4 100644 --- a/src/dependencies/hard_dependencies.jl +++ b/src/dependencies/hard_dependencies.jl @@ -110,6 +110,92 @@ function initialise_all_as_hard_dependency_node(models, scale) return dep_graph end +# Samuel : this requires the orchestrator, which requires the dependency graph +# Leaving it in dependency_graph.jl causes forward declaration issues, moving it here as a quick protoyping hack, it might not be the ideal spot +""" + variables_multiscale(node, organ, mapping, st=NamedTuple()) + +Get the variables of a HardDependencyNode, taking into account the multiscale mapping, *i.e.* +defining variables as `MappedVar` if they are mapped to another scale. The default values are +taken from the model if not given by the user (`st`), and are marked as `UninitializedVar` if +they are inputs of the node. + +Return a NamedTuple with the variables and their default values. + +# Arguments + +- `node::HardDependencyNode`: the node to get the variables from. +- `organ::String`: the organ type, *e.g.* "Leaf". +- `vars_mapping::Dict{String,T}`: the mapping of the models (see details below). +- `st::NamedTuple`: an optional named tuple with default values for the variables. + +# Details + +The `vars_mapping` is a dictionary with the organ type as key and a dictionary as value. It is +computed from the user mapping like so: +""" +function variables_multiscale(node, organ, vars_mapping, st=NamedTuple(), orchestrator::Orchestrator2=Orchestrator2()) + node_vars = variables(node) # e.g. (inputs = (:var1=-Inf, :var2=-Inf), outputs = (:var3=-Inf,)) + ins = node_vars.inputs + ins_variables = keys(ins) + outs_variables = keys(node_vars.outputs) + defaults = merge(node_vars...) + map((inputs=ins_variables, outputs=outs_variables)) do vars # Map over vars from :inputs and vars from :outputs + vars_ = Vector{Pair{Symbol,Any}}() + for var in vars # e.g. var = :carbon_biomass + if var in keys(st) + #If the user has given a status, we use it as default value. + default = st[var] + elseif var in ins_variables + # Otherwise, we use the default value given by the model: + # If the variable is an input, we mark it as uninitialized: + default = UninitializedVar(var, defaults[var]) + else + # If the variable is an output, we use the default value given by the model: + default = defaults[var] + end + + # TODO no idea how this meshes with refvector situations or previoustimestep + if is_timestep_mapped(organ => var, orchestrator, search_inputs_only=true) + + push!(vars_, var => default) + else + + if haskey(vars_mapping[organ], var) + organ_mapped, organ_mapped_var = _node_mapping(vars_mapping[organ][var]) + push!(vars_, var => MappedVar(organ_mapped, var, organ_mapped_var, default)) + #* We still check if the variable also exists wrapped in PreviousTimeStep, because one model could use the current + #* values, and another one the previous values. + if haskey(vars_mapping[organ], PreviousTimeStep(var, node.process)) + organ_mapped, organ_mapped_var = _node_mapping(vars_mapping[organ][PreviousTimeStep(var, node.process)]) + push!(vars_, var => MappedVar(organ_mapped, PreviousTimeStep(var, node.process), organ_mapped_var, default)) + end + elseif haskey(vars_mapping[organ], PreviousTimeStep(var, node.process)) + # If not found in the current time step, we check if the variable is mapped to the previous time step: + organ_mapped, organ_mapped_var = _node_mapping(vars_mapping[organ][PreviousTimeStep(var, node.process)]) + push!(vars_, var => MappedVar(organ_mapped, PreviousTimeStep(var, node.process), organ_mapped_var, default)) + else + # Else we take the default value: + push!(vars_, var => default) + end + end + end + return (; vars_...,) + end +end + +function _node_mapping(var_mapping::Pair{String,Symbol}) + # One organ is mapped to the variable: + return SingleNodeMapping(first(var_mapping)), last(var_mapping) +end + +function _node_mapping(var_mapping) + # Several organs are mapped to the variable: + organ_mapped = MultiNodeMapping([first(i) for i in var_mapping]) + organ_mapped_var = [last(i) for i in var_mapping] + + return organ_mapped, organ_mapped_var +end # When we use a mapping (multiscale), we return the set of soft-dependencies (we put the hard-dependencies as their children): function hard_dependencies(mapping::Dict{String,T}; verbose::Bool=true, orchestrator::Orchestrator2=Orchestrator2()) where {T} @@ -148,7 +234,7 @@ function hard_dependencies(mapping::Dict{String,T}; verbose::Bool=true, orchestr status_scale = Dict{Symbol,Vector{Pair{Symbol,NamedTuple}}}() for (procname, node) in hard_deps[organ].roots # procname = :leaf_surface ; node = hard_deps.roots[procname] var = Pair{Symbol,NamedTuple}[] - traverse_dependency_graph!(node, x -> variables_multiscale(x, organ, full_vars_mapping, st_scale_user), var) + traverse_dependency_graph!(node, x -> variables_multiscale(x, organ, full_vars_mapping, st_scale_user, orchestrator), var) push!(status_scale, procname => var) end diff --git a/src/dependencies/soft_dependencies.jl b/src/dependencies/soft_dependencies.jl index 241a63bf..135e32b6 100644 --- a/src/dependencies/soft_dependencies.jl +++ b/src/dependencies/soft_dependencies.jl @@ -357,7 +357,7 @@ function add_timestep_data_to_node(soft_dependency_node, orchestrator::Orchestra # now we can create the mapping for mtsm in orchestrator.non_default_timestep_mapping if mtsm.scale == soft_dependency_node.scale && (mtsm.model) == typeof(model_(soft_dependency_node.value)) - for (var_to, var_from) in mtsm.var_to_var + for (var_from, var_to) in mtsm.var_to_var if !isnothing(soft_dependency_node.parent) parent = nothing variable_mapping = nothing diff --git a/src/mtg/initialisation.jl b/src/mtg/initialisation.jl index eac56d0c..a4004ebe 100644 --- a/src/mtg/initialisation.jl +++ b/src/mtg/initialisation.jl @@ -35,18 +35,23 @@ function init_statuses(mtg, mapping, dependency_graph=first(hard_dependencies(ma # Note 3: we do it before `convert_reference_values!` because we need the variables to be MappedVar{MultiNodeMapping} to get the reverse mapping. # Convert the MappedVar{SelfNodeMapping} or MappedVar{SingleNodeMapping} to RefValues, and MappedVar{MultiNodeMapping} to RefVectors: - convert_reference_values!(mapped_vars) + convert_reference_values!(mapped_vars, orchestrator) # Get the variables that are not initialised or computed by other models in the output: vars_need_init = Dict(org => filter(x -> isa(last(x), UninitializedVar), vars) |> keys |> collect for (org, vars) in mapped_vars) |> filter(x -> length(last(x)) > 0) # Filter out variables that are timestep-mapped by the user, - # as those variables are initialized by another model, but are currently flagged as needing initialization + # Since we disconnect outputs from the source variable (as values are changed by the accumulation function, + # meaning they differ and we can't just Ref point to the source variable) + # they will be currently flagged as needing initialization # At this stage, data present in the orchestrator is expected to be valid, so we can take it into account # A model with a different timestep can still have unitialized vars found in a node, the meteo, or to be initialized by the user # in which case it'll be absent from the timestep mapping, but this needs testing - #filter_timestep_mapped_variables!(vars_need_init, orchestrator) + # TODO, *however*, this isn't the cleanest in its current state, + # there may be some user initialisation issues that are hidden by this approach + # Needs to be checked + filter_timestep_mapped_variables!(vars_need_init, orchestrator) # Note: these variables may be present in the MTG attributes, we check that below when traversing the MTG. @@ -54,7 +59,7 @@ function init_statuses(mtg, mapping, dependency_graph=first(hard_dependencies(ma # We traverse the MTG to initialise the statuses linked to the nodes: statuses = Dict(i => Status[] for i in collect(keys(mapped_vars))) MultiScaleTreeGraph.traverse!(mtg) do node # e.g.: node = MultiScaleTreeGraph.get_node(mtg, 5) - init_node_status!(node, statuses, mapped_vars, reverse_multiscale_mapping, vars_need_init, type_promotion, check=check) + init_node_status!(node, statuses, mapped_vars, reverse_multiscale_mapping, vars_need_init, type_promotion, check=check, orchestrator=orchestrator) end return (; statuses, mapped_vars, reverse_multiscale_mapping, vars_need_init) @@ -153,7 +158,7 @@ function init_node_status!(node, statuses, mapped_vars, reverse_multiscale_mappi end # Make the node status from the template: - st = status_from_template(st_template) + st = status_from_template(st_template, symbol(node), orchestrator) push!(statuses[symbol(node)], st) @@ -209,17 +214,23 @@ julia> b 2.0 ``` """ -function status_from_template(d::Dict{Symbol,T} where {T}) +function status_from_template(d::Dict{Symbol,T} where {T}, scale::String, orchestrator::Orchestrator2) # Sort vars to put the values that are of type PerStatusRef at the end (we need the pass on the other ones first): sorted_vars = Dict{Symbol,Any}(sort([pairs(d)...], by=v -> last(v) isa RefVariable ? 1 : 0)) # Note: PerStatusRef are used to reference variables in the same status for renaming. # We create the status with the right references for variables, and for PerStatusRef (we reference the reference variable): for (k, v) in sorted_vars - if isa(v, RefVariable) - sorted_vars[k] = sorted_vars[v.reference_variable] - else + if is_timestep_mapped((scale => k), orchestrator, search_inputs_only=true) + # avoid referring to the original variable sorted_vars[k] = ref_var(v) + else + + if isa(v, RefVariable) + sorted_vars[k] = sorted_vars[v.reference_variable] + else + sorted_vars[k] = ref_var(v) + end end end @@ -423,17 +434,22 @@ function filter_timestep_mapped_variables!(vars_need_init, orchestrator) end=# function filter_timestep_mapped_variables!(vars_need_init, orchestrator) - for (org, vars) in vars_need_init - for tmst in orchestrator.non_default_timestep_mapping + for tmst in orchestrator.non_default_timestep_mapping + for (org, vars) in vars_need_init if tmst.scale == org - - for (var_to, var_from) in tmst.var_to_var - filter!(o -> o == var_to.name || o == var_from.name, vars) - end - if isempty(vars) - delete!(vars_need_init, org) + for (var_from, var_to) in tmst.var_to_var + filter!(o -> o == var_to.name, vars) + end + end + for (var_from, var_to) in tmst.var_to_var + if var_from.scale == org + filter!(o -> o == var_from.name, vars) end end + + if isempty(vars) + delete!(vars_need_init, org) + end end end end \ No newline at end of file diff --git a/src/mtg/mapping/compute_mapping.jl b/src/mtg/mapping/compute_mapping.jl index 029b9c3d..872fa67e 100644 --- a/src/mtg/mapping/compute_mapping.jl +++ b/src/mtg/mapping/compute_mapping.jl @@ -279,6 +279,18 @@ function default_variables_from_mapping(mapped_vars, verbose=true) end +function is_timestep_mapped(key, orchestrator::Orchestrator2; search_inputs_only::Bool=false) + for mtsm in orchestrator.non_default_timestep_mapping + for (var_from, var_to) in mtsm.var_to_var + if ((first(key) == mtsm.scale) && last(key) == var_to.name) || + (!search_inputs_only && (first(key) == var_from.scale && last(key) == var_from.name)) + return true + end + end + end + return false +end + """ convert_reference_values!(mapped_vars::Dict{String,Dict{Symbol,Any}}) @@ -286,7 +298,7 @@ Convert the variables that are `MappedVar{SelfNodeMapping}` or `MappedVar{Single common value for the variable; and convert `MappedVar{MultiNodeMapping}` to RefVectors that reference the values for the variable in the source organs. """ -function convert_reference_values!(mapped_vars::Dict{String,Dict{Symbol,Any}}) +function convert_reference_values!(mapped_vars::Dict{String,Dict{Symbol,Any}}, orchestrator::Orchestrator2) # For the variables that will be RefValues, i.e. referencing a value that exists for different scales, we need to first # create a common reference to the value that we use wherever we need this value. These values are created in the dict_mapped_vars # Dict, and then referenced from there every time we point to it. @@ -303,7 +315,11 @@ function convert_reference_values!(mapped_vars::Dict{String,Dict{Symbol,Any}}) # First time we encounter this variable as a MappedVar, we create its value into the dict_mapped_vars Dict: if !haskey(dict_mapped_vars, key) - push!(dict_mapped_vars, key => Ref(mapped_default(vars[k]))) +# if is_timestep_mapped(key, orchestrator) +# push!(dict_mapped_vars, key => (mapped_default(vars[k]))) +# else + push!(dict_mapped_vars, key => Ref(mapped_default(vars[k]))) +# end end # Then we use the value for the particular variable to replace the MappedVar to a RefValue in the mapping: diff --git a/src/mtg/save_results.jl b/src/mtg/save_results.jl index 51baffcb..52a94ba8 100644 --- a/src/mtg/save_results.jl +++ b/src/mtg/save_results.jl @@ -109,8 +109,8 @@ julia> collect(keys(preallocated_vars["Leaf"])) :carbon_demand ``` """ - -function pre_allocate_outputs(statuses, statuses_template, reverse_multiscale_mapping, vars_need_init, outs, nsteps; type_promotion=nothing, check=true) +# TODO orchestrator prob shouldn't be a kwarg with a default +function pre_allocate_outputs(statuses, statuses_template, reverse_multiscale_mapping, vars_need_init, outs, nsteps; type_promotion=nothing, check=true, orchestrator=Orchestrator2()) outs_ = Dict{String,Vector{Symbol}}() # default behaviour : track everything @@ -211,18 +211,18 @@ function pre_allocate_outputs(statuses, statuses_template, reverse_multiscale_ma # I don't know if this function barrier is necessary preallocated_outputs = Dict{String,Vector}() - complete_preallocation_from_types!(preallocated_outputs, nsteps, outs_, node_type, statuses_template) + complete_preallocation_from_types!(preallocated_outputs, nsteps, outs_, node_type, statuses_template, orchestrator) return preallocated_outputs end -function complete_preallocation_from_types!(preallocated_outputs, nsteps, outs_, node_type, statuses_template) +function complete_preallocation_from_types!(preallocated_outputs, nsteps, outs_, node_type, statuses_template, orchestrator) types = Vector{DataType}() for organ in keys(outs_) outs_no_node = filter(x -> x != :node, outs_[organ]) #types = [typeof(status_from_template(statuses_template[organ])[var]) for var in outs[organ]] - values = [status_from_template(statuses_template[organ])[var] for var in outs_no_node] + values = [status_from_template(statuses_template[organ], organ, orchestrator)[var] for var in outs_no_node] #push!(types, node_type) diff --git a/src/run.jl b/src/run.jl index 8a62a200..3bac7595 100644 --- a/src/run.jl +++ b/src/run.jl @@ -454,7 +454,11 @@ function run_node_multiscale!( # Check if the model needs to run this timestep if (1 + (i-1) % node_ratio) != node_ratio - node.simulation_id[1] += 1 + + # TODO : This does prevent the node form being updated by two different parents but probably should be changed + if any([p.simulation_id[1] >= node.simulation_id[1] for p in node.parent]) + node.simulation_id[1] += 1 + end return nothing end @@ -468,6 +472,7 @@ function run_node_multiscale!( node_statuses = status(object)[node.scale] # Get the status of the nodes at the current scale models_at_scale = models[node.scale] + # TODO : is empty check should be pre simulation if isnothing(node.timestep_mapping_data) || isempty(node.timestep_mapping_data) # Samuel : this is the happy path, no further timestep mapping checks needed diff --git a/src/timestep/timestep_mapping.jl b/src/timestep/timestep_mapping.jl index 20ad038e..06d40d01 100644 --- a/src/timestep/timestep_mapping.jl +++ b/src/timestep/timestep_mapping.jl @@ -78,7 +78,7 @@ struct ModelTimestepMapping model scale::String timestep::Period - var_to_var::Dict{Var_to, Var_from} + var_to_var::Dict{Var_from, Var_to} end mutable struct Orchestrator2 diff --git a/test/test_multitimestep.jl b/test/test_multitimestep.jl index 9e25f1ce..b9be19da 100644 --- a/test/test_multitimestep.jl +++ b/test/test_multitimestep.jl @@ -297,14 +297,65 @@ mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 2)) to = PlantSimEngine.Var_to(:weekly_max_temperature) from = PlantSimEngine.Var_from(MyToyDayModel, "Default", :daily_temperature, maximum) -dict_to_from = Dict(to => from) +dict_to_from = Dict(from => to) mtsm = PlantSimEngine.ModelTimestepMapping(MyToyWeekModel, "Default2", Week(1), dict_to_from) orch2 = PlantSimEngine.Orchestrator2(Day(1), [mtsm,]) -out = @enter run!(mtg, m_multiscale, meteo_day, orchestrator=orch2) +out = run!(mtg, m_multiscale, meteo_day, orchestrator=orch2) +########################### +# Three timestep model that is single-scale, to circumvent refvector/refvalue overwriting +# and explore alternatives +# (eg filtering out timestep-mapped variables from vars_need_init and storing the values elsewhere) +########################### + + m_singlescale = Dict("Default" => ( + MyToyDay2Model(), + MyToyWeek2Model(), + MyToyFourWeek2Model(), + ),) + + + model_timesteps_defaultscale = Dict(MyToyWeek2Model =>Week(1), MyToyFourWeek2Model =>Week(4), ) + tsm_d = PlantSimEngine.TimestepMapper(:out_day, Day(1), sum, nothing) + sth_d = PlantSimEngine.SimulationTimestepHandler(model_timesteps_defaultscale, Dict(:in_week => tsm_d )) + + tsm_w = PlantSimEngine.TimestepMapper(:out_week, Week(1), sum, nothing) + sth_w = PlantSimEngine.SimulationTimestepHandler(model_timesteps_defaultscale, Dict(:in_four_week_from_day => tsm_d, :in_four_week_from_week => tsm_w )) + + to_w = PlantSimEngine.Var_to(:in_week) + from_d = PlantSimEngine.Var_from(MyToyDay2Model, "Default", :out_day, sum) + dict_to_from_w = Dict(to_w => from_d) + + to_w4_d = PlantSimEngine.Var_to(:in_four_week_from_day) + to_w4_w = PlantSimEngine.Var_to(:in_four_week_from_week) + from_w = PlantSimEngine.Var_from(MyToyWeek2Model, "Default", :out_week, sum) + + dict_to_from_w4 = Dict(from_d => to_w4_d, from_w => to_w4_w) + + mtsm_w = PlantSimEngine.ModelTimestepMapping(MyToyWeek2Model, "Default", Week(1), dict_to_from_w) + mtsm_w4 = PlantSimEngine.ModelTimestepMapping(MyToyFourWeek2Model, "Default", Week(4), dict_to_from_w4) + + orch2 = PlantSimEngine.Orchestrator2(Day(1), [mtsm_w, mtsm_w4]) + +mtg_single = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) +out = @run run!(mtg_single, m_singlescale, df, orchestrator=orch2) + + + + + + + + + + + + + + ########################### # Test with three timesteps, multiscale # Issues with data overwriting (refvector/refvalue) @@ -326,6 +377,28 @@ function PlantSimEngine.run!(m::MyToyDay2Model, models, status, meteo, constants status.out_day = meteo.data end +#=PlantSimEngine.@process "ToyDay3" verbose = false + +struct MyToyDay3Model <: AbstractToyday3Model end + +PlantSimEngine.inputs_(m::MyToyDay3Model) = (in_day=-Inf, in_day_summed_prev_timestep=-Inf) +PlantSimEngine.outputs_(m::MyToyDay3Model) = (out_day_summed=-Inf,) + +function PlantSimEngine.run!(m::MyToyDay3Model, models, status, meteo, constants=nothing, extra=nothing) + status.out_day_summed = status.in_day + status.in_day_summed_prev_timestep +end + +PlantSimEngine.@process "ToyDay4" verbose = false + +struct MyToyDay4Model <: AbstractToyday4Model end + +PlantSimEngine.inputs_(m::MyToyDay4Model) = (in_day_summed=-Inf,) +PlantSimEngine.outputs_(m::MyToyDay4Model) = (out_day_summed_2= -Inf,) + +function PlantSimEngine.run!(m::MyToyDay4Model, models, status, meteo, constants=nothing, extra=nothing) + status.out_day_summed_2 = status.in_day_summed +end=# + PlantSimEngine.@process "ToyWeek2" verbose = false struct MyToyWeek2Model <: AbstractToyweek2Model end @@ -367,13 +440,13 @@ df = DataFrame(:data => [1 for i in 1:365], ) to_w = PlantSimEngine.Var_to(:in_week) from_d = PlantSimEngine.Var_from(MyToyDay2Model, "Default", :out_day, sum) - dict_to_from_w = Dict(to_w => from_d) + dict_to_from_w = Dict(from_d => to_w) to_w4_d = PlantSimEngine.Var_to(:in_four_week_from_day) to_w4_w = PlantSimEngine.Var_to(:in_four_week_from_week) from_w = PlantSimEngine.Var_from(MyToyWeek2Model, "Default2", :out_week, sum) - dict_to_from_w4 = Dict(to_w4_d => from_d, to_w4_w => from_w) + dict_to_from_w4 = Dict(from_d => to_w4_d, from_w => to_w4_w) mtsm_w = PlantSimEngine.ModelTimestepMapping(MyToyWeek2Model, "Default2", Week(1), dict_to_from_w) mtsm_w4 = PlantSimEngine.ModelTimestepMapping(MyToyFourWeek2Model, "Default3", Week(4), dict_to_from_w4) @@ -395,55 +468,36 @@ m_multiscale = Dict("Default" => ( :in_four_week_from_day => "Default" => :out_day, :in_four_week_from_week => "Default2" => :out_week, ], + ),), + #="Default4"=> ( + MultiScaleModel( + model=MyToyDay3Model(), + mapped_variables=[ + PlantSimEngine.PreviousTimeStep(:in_day_summed_prev_timestep) => "Default5" => :out_day_summed_2, + :in_day => "Default" => :out_day, + ]), + Status(in_day_summed_prev_timestep=0,) ), - ),) + "Default5" => ( + MultiScaleModel(model=MyToyDay4Model(), + mapped_variables= [:in_day_summed => "Default4" => :out_day_summed], + ), + Status(in_day_summed=0,out_day_summed_2=0) + ),=# + ) # TODO test with multiple nodes mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 2)) mtg3 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default3", 1, 3)) -#out = @enter run!(mtg, m_multiscale, df, orchestrator=orch2) - - - -########################### -# Three timestep model that is single-scale, to circumvent refvector/refvalue overwriting -# and explore alternatives -# (eg filtering out timestep-mapped variables from vars_need_init and storing the values elsewhere) -########################### - - m_singlescale = Dict("Default" => ( - MyToyDay2Model(), - MyToyWeek2Model(), - MyToyFourWeek2Model(), - ),) - - - model_timesteps_defaultscale = Dict(MyToyWeek2Model =>Week(1), MyToyFourWeek2Model =>Week(4), ) - tsm_d = PlantSimEngine.TimestepMapper(:out_day, Day(1), sum, nothing) - sth_d = PlantSimEngine.SimulationTimestepHandler(model_timesteps_defaultscale, Dict(:in_week => tsm_d )) +#mtg4 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default4", 1, 4)) +#mtg5 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default5", 1, 5)) - tsm_w = PlantSimEngine.TimestepMapper(:out_week, Week(1), sum, nothing) - sth_w = PlantSimEngine.SimulationTimestepHandler(model_timesteps_defaultscale, Dict(:in_four_week_from_day => tsm_d, :in_four_week_from_week => tsm_w )) - - to_w = PlantSimEngine.Var_to(:in_week) - from_d = PlantSimEngine.Var_from(MyToyDay2Model, "Default", :out_day, sum) - dict_to_from_w = Dict(to_w => from_d) + #orch2 = PlantSimEngine.Orchestrator2() - to_w4_d = PlantSimEngine.Var_to(:in_four_week_from_day) - to_w4_w = PlantSimEngine.Var_to(:in_four_week_from_week) - from_w = PlantSimEngine.Var_from(MyToyWeek2Model, "Default", :out_week, sum) - - dict_to_from_w4 = Dict(to_w4_d => from_d, to_w4_w => from_w) - - mtsm_w = PlantSimEngine.ModelTimestepMapping(MyToyWeek2Model, "Default", Week(1), dict_to_from_w) - mtsm_w4 = PlantSimEngine.ModelTimestepMapping(MyToyFourWeek2Model, "Default", Week(4), dict_to_from_w4) - - orch2 = PlantSimEngine.Orchestrator2(Day(1), [mtsm_w, mtsm_w4]) - -mtg_single = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) -out = @run run!(mtg_single, m_singlescale, df, orchestrator=orch2) +out = @run run!(mtg, m_multiscale, df, orchestrator=orch2) +#out = run!(mtg, m_multiscale, df, orchestrator=orch2) From d964abdb820d537dbdbd13ec6a2c6ac9c047e75e Mon Sep 17 00:00:00 2001 From: Samuel-AMAP Date: Fri, 31 Oct 2025 15:43:36 +0100 Subject: [PATCH 6/8] Fix issues with simulation id and models not running properly. Many loose ends to investigate --- src/run.jl | 119 ++++++++++++++++++++++++++--------------------------- 1 file changed, 59 insertions(+), 60 deletions(-) diff --git a/src/run.jl b/src/run.jl index 3bac7595..42ee7355 100644 --- a/src/run.jl +++ b/src/run.jl @@ -453,77 +453,76 @@ function run_node_multiscale!( node_ratio = node.timestep / Day(1) # Check if the model needs to run this timestep - if (1 + (i-1) % node_ratio) != node_ratio + if (1 + (i - 1) % node_ratio) != node_ratio # TODO : This does prevent the node form being updated by two different parents but probably should be changed - if any([p.simulation_id[1] >= node.simulation_id[1] for p in node.parent]) + if any([p.simulation_id[1] > node.simulation_id[1] for p in node.parent]) node.simulation_id[1] += 1 end - return nothing - end - - # run!(status(object), dependency_node, meteo, constants, extra) - # Check if all the parents have been called before the child: - if !AbstractTrees.isroot(node) && any([p.simulation_id[1] <= node.simulation_id[1] for p in node.parent]) - # If not, this node should be called via another parent - return nothing - end + else - node_statuses = status(object)[node.scale] # Get the status of the nodes at the current scale - models_at_scale = models[node.scale] + # run!(status(object), dependency_node, meteo, constants, extra) + # Check if all the parents have been called before the child: + if !AbstractTrees.isroot(node) && any([p.simulation_id[1] <= node.simulation_id[1] for p in node.parent]) + # If not, this node should be called via another parent + return nothing + end - # TODO : is empty check should be pre simulation - if isnothing(node.timestep_mapping_data) || isempty(node.timestep_mapping_data) - # Samuel : this is the happy path, no further timestep mapping checks needed + node_statuses = status(object)[node.scale] # Get the status of the nodes at the current scale + models_at_scale = models[node.scale] - for st in node_statuses # for each node status at the current scale (potentially in parallel over nodes) - # Actual call to the model: - run!(node.value, models_at_scale, st, meteo, constants, extra) - end - node.simulation_id[1] += 1 # increment the simulation id, to remember that the model has been called already + # TODO : is empty check should be pre simulation + if isnothing(node.timestep_mapping_data) || isempty(node.timestep_mapping_data) + # Samuel : this is the happy path, no further timestep mapping checks needed - # Recursively visit the children (soft dependencies only, hard dependencies are handled by the model itself): - for child in node.children - #! check if we can run this safely in a @floop loop. I would say no, - #! because we are running a parallel computation above already, modifying the node.simulation_id, - #! which is not thread-safe yet. - run_node_multiscale!(object, child, i, models, meteo, constants, extra, check, executor) - end - else - # we have a child with a different timestep than the parent, - # which requires accumulation/reduce and running only some of the time - - for st in node_statuses - run!(node.value, models_at_scale, st, meteo, constants, extra) - # TODO do the accumulation - # then write into the child's status if need be ? - for tmst in node.timestep_mapping_data - ratio = Int(tmst.node_to.timestep / node.timestep) - # TODO assert etc. - # do the accumulation for each variable - accumulation_index = 1 + ((i-1)%ratio) - tmst.mapping_data[node_id(st.node)][accumulation_index] = st[tmst.variable_from] - - # if we have completed a *full* cycle, then presumably we need to write out the value to - # the mapped model - # A full cycle isn't just the ratio to the parent, it's the ratio to the finest-grained timestep - if accumulation_index == ratio - node_statuses_to = status(object)[tmst.node_to.scale] - - # TODO : INCORRECT in a scale with multiple mtg nodes - for st_to in node_statuses_to - # TODO might be able to catch mapping_function type incompatibility errors and make them clearer - st_to[tmst.variable_to] = tmst.mapping_function(tmst.mapping_data[node_id(st.node)]) + for st in node_statuses # for each node status at the current scale (potentially in parallel over nodes) + # Actual call to the model: + run!(node.value, models_at_scale, st, meteo, constants, extra) + end + node.simulation_id[1] += 1 # increment the simulation id, to remember that the model has been called already + + # Recursively visit the children (soft dependencies only, hard dependencies are handled by the model itself): + for child in node.children + #! check if we can run this safely in a @floop loop. I would say no, + #! because we are running a parallel computation above already, modifying the node.simulation_id, + #! which is not thread-safe yet. + run_node_multiscale!(object, child, i, models, meteo, constants, extra, check, executor) + end + else + # we have a child with a different timestep than the parent, + # which requires accumulation/reduce and running only some of the time + + for st in node_statuses + run!(node.value, models_at_scale, st, meteo, constants, extra) + # TODO do the accumulation + # then write into the child's status if need be ? + for tmst in node.timestep_mapping_data + ratio = Int(tmst.node_to.timestep / node.timestep) + # TODO assert etc. This is all assuming the ratio is an integer, whereas it can be, like 1/7 + # do the accumulation for each variable + index = Int(i*Day(1) / node.timestep) + accumulation_index = 1 + ((index - 1) % ratio) + tmst.mapping_data[node_id(st.node)][accumulation_index] = st[tmst.variable_from] + + # if we have completed a *full* cycle, then presumably we need to write out the value to + # the mapped model + # A full cycle isn't just the ratio to the parent, it's the ratio to the finest-grained timestep + if accumulation_index == ratio + node_statuses_to = status(object)[tmst.node_to.scale] + + # TODO : INCORRECT in a scale with multiple mtg nodes + for st_to in node_statuses_to + # TODO might be able to catch mapping_function type incompatibility errors and make them clearer + st_to[tmst.variable_to] = tmst.mapping_function(tmst.mapping_data[node_id(st.node)]) + end end end - end - end - - node.simulation_id[1] += 1 + end - for child in node.children - - run_node_multiscale!(object, child, i, models, meteo, constants, extra, check, executor) + node.simulation_id[1] += 1 end end + for child in node.children + run_node_multiscale!(object, child, i, models, meteo, constants, extra, check, executor) + end end \ No newline at end of file From 80e12de17a5a9a3eafae9c04f149df5db68369ca Mon Sep 17 00:00:00 2001 From: Samuel-AMAP Date: Fri, 31 Oct 2025 16:29:39 +0100 Subject: [PATCH 7/8] Reintroduce error message which I may have unintentionally removed earlier --- src/mtg/initialisation.jl | 1 + 1 file changed, 1 insertion(+) diff --git a/src/mtg/initialisation.jl b/src/mtg/initialisation.jl index a4004ebe..1f604840 100644 --- a/src/mtg/initialisation.jl +++ b/src/mtg/initialisation.jl @@ -126,6 +126,7 @@ function init_node_status!(node, statuses, mapped_vars, reverse_multiscale_mappi end continue end + error("Variable `$(var)` is not computed by any model, not initialised by the user in the status, and not found in the MTG at scale $(symbol(node)) (checked for MTG node $(node_id(node))).") end # Applying the type promotion to the node attribute if needed: if isnothing(type_promotion) From 2ee5f85ee7dec5d948f1c8e1b100461e1b07fb22 Mon Sep 17 00:00:00 2001 From: Samuel-AMAP Date: Mon, 3 Nov 2025 15:46:29 +0100 Subject: [PATCH 8/8] Fix bug in timestep-mapped variables filtering from var_need_init --- src/mtg/initialisation.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/mtg/initialisation.jl b/src/mtg/initialisation.jl index 1f604840..f738ee42 100644 --- a/src/mtg/initialisation.jl +++ b/src/mtg/initialisation.jl @@ -439,12 +439,12 @@ function filter_timestep_mapped_variables!(vars_need_init, orchestrator) for (org, vars) in vars_need_init if tmst.scale == org for (var_from, var_to) in tmst.var_to_var - filter!(o -> o == var_to.name, vars) + filter!(o -> o != var_to.name, vars) end end for (var_from, var_to) in tmst.var_to_var if var_from.scale == org - filter!(o -> o == var_from.name, vars) + filter!(o -> o != var_from.name, vars) end end