diff --git a/src/Abstract_model_structs.jl b/src/Abstract_model_structs.jl index 9cef7b08..68f91986 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) +""" +timestep_valid(tsr::TimestepRange) = tsr.lower_bound <= tsr.upper_bound + +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 == Second(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 ? + +# TODO properly test on models with default timestep \ No newline at end of file diff --git a/src/PlantSimEngine.jl b/src/PlantSimEngine.jl index ebf59dcb..a4ff5c83 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") @@ -56,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") @@ -102,7 +106,9 @@ include("examples_import.jl") export PreviousTimeStep export AbstractModel -export ModelList, MultiScaleModel +export ModelList, MultiScaleModel, TimestepMappedVariable +export MultiScaleMapping +export Orchestrator, TimestepRange, 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..6672add8 100644 --- a/src/dependencies/dependencies.jl +++ b/src/dependencies/dependencies.jl @@ -99,27 +99,3 @@ end 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} - # 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) - # 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) - # 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. - - # Checking that the graph is acyclic: - iscyclic, cycle_vec = is_graph_cyclic(dep_graph; warn=false) - # Note: we could do that in `soft_dependencies_multiscale` but we prefer to keep the function as simple as possible, and - # usable on its own. - - iscyclic && error("Cyclic dependency detected in the graph. Cycle: \n $(print_cycle(cycle_vec)) \n You can break the cycle using the `PreviousTimeStep` variable in the mapping.") - # Third step, we identify which - return dep_graph -end diff --git a/src/dependencies/dependency_graph.jl b/src/dependencies/dependency_graph.jl index 2be06b64..e31ea4aa 100644 --- a/src/dependencies/dependency_graph.jl +++ b/src/dependencies/dependency_graph.jl @@ -6,23 +6,35 @@ mutable struct HardDependencyNode{T} <: AbstractDependencyNode dependency::NamedTuple missing_dependency::Vector{Int} scale::String - inputs - outputs + #inputs + #outputs parent::Union{Nothing,<:AbstractDependencyNode} children::Vector{HardDependencyNode} end +mutable struct TimestepMapping + variable_from::Symbol + variable_to::Symbol + timestep_to::Period + mapping_function::Function + mapping_data_template + mapping_data::Dict{Int, Any} # TODO fix type stability : Int is the node id, Any is a vector of n elements of the variable's type, n being the # of required timesteps +end + +# can hard dependency nodes also handle timestep mapped variables... ? mutable struct SoftDependencyNode{T} <: AbstractDependencyNode value::T process::Symbol scale::String - inputs - outputs + #inputs + #outputs hard_dependency::Vector{HardDependencyNode} parent::Union{Nothing,Vector{SoftDependencyNode}} 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: @@ -77,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 0d1fa5b3..ae596d1e 100644 --- a/src/dependencies/hard_dependencies.jl +++ b/src/dependencies/hard_dependencies.jl @@ -100,8 +100,8 @@ function initialise_all_as_hard_dependency_node(models, scale) NamedTuple(), Int[], scale, - inputs_(i), - outputs_(i), + #inputs_(i), + #outputs_(i), nothing, HardDependencyNode[] ) for (p, i) in pairs(models) @@ -110,9 +110,100 @@ 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::Orchestrator=Orchestrator()) + 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 +# 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 + +function extract_timestep_mapped_outputs(m::MultiScaleModel, organ::String, outputs_process, timestep_mapped_outputs_process) + if length(m.timestep_mapped_variables) > 0 + if !haskey(timestep_mapped_outputs_process, organ) + timestep_mapped_outputs_process[organ] = Dict{Symbol,Vector}() + end + key = process(m.model) + extra_outputs = timestep_mapped_outputs_(m) + #ind = findfirst(x -> first(x) == key, outputs_process[organ][key]) + timestep_mapped_outputs_process[organ][key] = extra_outputs #TODO + end +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::Orchestrator=Orchestrator()) 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}() @@ -136,6 +227,7 @@ function hard_dependencies(mapping::Dict{String,T}; verbose::Bool=true) where {T #* dependency may not appear in its own scale, but this is treated in the soft-dependency computation inputs_process = Dict{String,Dict{Symbol,Vector}}() outputs_process = Dict{String,Dict{Symbol,Vector}}() + timestep_mapped_outputs_process = Dict{String,Dict{Symbol,Vector}}() for (organ, model) in mapping # Get the status given by the user, that is used to set the default values of the variables in the mapping: st_scale_user = get_status(model) @@ -148,12 +240,38 @@ function hard_dependencies(mapping::Dict{String,T}; verbose::Bool=true) where {T 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 inputs_process[organ] = Dict(key => [j.first => j.second.inputs for j in val] for (key, val) in status_scale) outputs_process[organ] = Dict(key => [j.first => j.second.outputs for j in val] for (key, val) in status_scale) + + # Samuel : This if else loop is a bit awkward + # None of the other code works this way, it uses the dependency grpah + # but the hard_dep graph loses the multiscale model information... + if isa(model, AbstractModel) + elseif isa(model, MultiScaleModel) + extract_timestep_mapped_outputs(model, organ, outputs_process, timestep_mapped_outputs_process) + else + for m in model + if isa(m, MultiScaleModel) + extract_timestep_mapped_outputs(m, organ, outputs_process, timestep_mapped_outputs_process) + end + end + end + + #=for m in model + if isa(m, MultiScaleModel) + if length(m.timestep_mapped_variables) > 0 + key = process(m.model) + extra_outputs = timestep_mapped_outputs_(m) + ind = findfirst(x -> first(x) == key, outputs_process[organ][key]) + outputs_process[organ][key][ind] = first(outputs_process[organ][key][ind]) => (; last(outputs_process[organ][key][ind])..., extra_outputs...) + + end + end + end=# end # If some models needed as hard-dependency are not found in their own scale, check the other scales: @@ -194,8 +312,8 @@ function hard_dependencies(mapping::Dict{String,T}; verbose::Bool=true) where {T dep_node_model.dependency, dep_node_model.missing_dependency, dep_node_model.scale, - dep_node_model.inputs, - dep_node_model.outputs, + #dep_node_model.inputs, + #dep_node_model.outputs, parent_node, dep_node_model.children ) @@ -226,23 +344,72 @@ 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( soft_dep_vars.value, process_, # process name organ, # scale - inputs_process[organ][process_], # These are the inputs, potentially multiscale - outputs_process[organ][process_], # Same for outputs + #inputs_process[organ][process_], # These are the inputs, potentially multiscale + #outputs_process[organ][process_], # Same for outputs AbstractTrees.children(soft_dep_vars), # hard dependencies 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_] ) - + for (process_, soft_dep_vars) in hard_deps[organ].roots + # TODO this is not good enough for some model ranges, and doesn't check for inconsistencies errors for models that have a modeltimestepmapping + if timestep_range_(soft_dep_vars.value).lower_bound == timestep_range_(soft_dep_vars.value).upper_bound + timestep = timestep_range_(soft_dep_vars.value).lower_bound + + # if the model has infinite range, set it to the simulation timestep + if timestep == Second(0) + timestep = orchestrator.default_timestep + end + soft_dep_graph[process_].timestep = timestep + end + end # Update the parent node of the hard dependency nodes to be the new SoftDependencyNode instead of the old # HardDependencyNode. for (p, node) in soft_dep_graph @@ -251,7 +418,7 @@ function hard_dependencies(mapping::Dict{String,T}; verbose::Bool=true) where {T end end - soft_dep_graphs[organ] = (soft_dep_graph=soft_dep_graph, inputs=inputs_process[organ], outputs=outputs_process[organ]) + soft_dep_graphs[organ] = (soft_dep_graph=soft_dep_graph, inputs=inputs_process[organ], outputs=outputs_process[organ], timestep_mapped_outputs=haskey(timestep_mapped_outputs_process,organ) ? timestep_mapped_outputs_process[organ] : Dict{Symbol,Vector}()) not_found = merge(not_found, hard_deps[organ].not_found) end diff --git a/src/dependencies/soft_dependencies.jl b/src/dependencies/soft_dependencies.jl index b772a207..12b98c92 100644 --- a/src/dependencies/soft_dependencies.jl +++ b/src/dependencies/soft_dependencies.jl @@ -64,13 +64,15 @@ function soft_dependencies(d::DependencyGraph{Dict{Symbol,HardDependencyNode}}, soft_dep_vars.value, process_, # process name "", - inputs_(soft_dep_vars.value), - outputs_(soft_dep_vars.value), + #inputs_(soft_dep_vars.value), + #outputs_(soft_dep_vars.value), AbstractTrees.children(soft_dep_vars), # hard dependencies nothing, nothing, SoftDependencyNode[], - fill(0, nsteps) + fill(0, nsteps), + Day(1), # TODO should be orchestrator.default_timestep, + nothing ) for (process_, soft_dep_vars) in d.roots ) @@ -139,16 +141,14 @@ function soft_dependencies(d::DependencyGraph{Dict{Symbol,HardDependencyNode}}, 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} - mapped_vars = mapped_variables(mapping, soft_dep_graphs_roots, verbose=false) - rev_mapping = reverse_mapping(mapped_vars, all=false) - +function soft_dependencies_multiscale(soft_dep_graphs_roots::DependencyGraph{Dict{String,Any}}, reverse_multiscale_mapping, hard_dep_dict::Dict{Pair{Symbol,String},HardDependencyNode}) independant_process_root = Dict{Pair{String,Symbol},SoftDependencyNode}() - for (organ, (soft_dep_graph, ins, outs)) in soft_dep_graphs_roots.roots # e.g. organ = "Plant"; soft_dep_graph, ins, outs = soft_dep_graphs_roots.roots[organ] + for (organ, (soft_dep_graph, ins, outs, timestep_mapped_outs)) in soft_dep_graphs_roots.roots # e.g. organ = "Plant"; soft_dep_graph, ins, outs = soft_dep_graphs_roots.roots[organ] + for (proc, i) in soft_dep_graph # proc = :leaf_surface; i = soft_dep_graph[proc] # Search if the process has soft dependencies: - soft_deps = search_inputs_in_output(proc, ins, outs) + soft_deps = search_inputs_in_output(proc, ins, outs, timestep_mapped_outs) # Remove the hard dependencies from the soft dependencies: soft_deps_not_hard = drop_process(soft_deps, [hd.process for hd in i.hard_dependency]) @@ -158,7 +158,7 @@ function soft_dependencies_multiscale(soft_dep_graphs_roots::DependencyGraph{Dic # NB: if a node is already a hard dependency of the node, it cannot be a soft dependency # Check if the process has soft dependencies at other scales: - soft_deps_multiscale = search_inputs_in_multiscale_output(proc, organ, ins, soft_dep_graphs_roots.roots, rev_mapping, hard_dependencies_from_other_scale) + soft_deps_multiscale = search_inputs_in_multiscale_output(proc, organ, ins, soft_dep_graphs_roots.roots, reverse_multiscale_mapping, hard_dependencies_from_other_scale) # Example output: "Soil" => Dict(:soil_water=>[:soil_water_content]), which means that the variable :soil_water_content # is computed by the process :soil_water at the scale "Soil". @@ -321,14 +321,47 @@ function soft_dependencies_multiscale(soft_dep_graphs_roots::DependencyGraph{Dic end end end + + if haskey(timestep_mapped_outs, proc) + create_timestep_mapping(i, timestep_mapped_outs, proc) + end + 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) + + return dep_graph end +function create_timestep_mapping(soft_dependency_node, timestep_mapped_outs, proc) -""" + for (proc_mapped, mapping_entries) in timestep_mapped_outs + + if proc_mapped == proc + for (var_to, (timestep_mapped_var_data, default_value)) in mapping_entries + timestep_ratio = timestep_mapped_var_data.timestep_to / soft_dependency_node.timestep + if timestep_ratio > 1 # todo assert it's an int or a rational ? + @assert timestep_ratio == trunc(timestep_ratio) "Error : non-integer timestep ratio" + + tmvd = timestep_mapped_var_data + var_type = typeof(default_value) + + mapping_data_template = Vector{var_type}(undef, convert(Int64, timestep_ratio)) + tsm = TimestepMapping(tmvd.name_from, tmvd.name_to, timestep_mapped_var_data.timestep_to, tmvd.aggregation_function, mapping_data_template, Dict{MultiScaleTreeGraph.NodeMTG,Any}()) + + if isnothing(soft_dependency_node.timestep_mapping_data) + soft_dependency_node.timestep_mapping_data = Vector{TimestepMapping}() + end + push!(soft_dependency_node.timestep_mapping_data, tsm) + end + end + end + end +end + + + """ drop_process(proc_vars, process) Return a new `NamedTuple` with the process `process` removed from the `NamedTuple` `proc_vars`. @@ -398,7 +431,17 @@ search_inputs_in_output(:process3, in_, out_) (process4 = (:var1, :var2),) ``` """ -function search_inputs_in_output(process, inputs, outputs) + + +function extract_mapped_outputs(timestep_mapped_outputs) + extracted = Pair[] + for pair in timestep_mapped_outputs + push!(extracted, Pair(first(last(pair)).name_to, last(last(pair)))) + end + return extracted +end + +function search_inputs_in_output(process, inputs, outputs, timestep_mapped_outputs=Dict{Symbol,NamedTuple}()) # proc, ins, outs # get the inputs of the node: vars_input = flatten_vars(inputs[process]) @@ -407,7 +450,12 @@ function search_inputs_in_output(process, inputs, outputs) for (proc_output, pairs_vars_output) in outputs # e.g. proc_output = :carbon_biomass; pairs_vars_output = outs[proc_output] if process != proc_output vars_output = flatten_vars(pairs_vars_output) - inputs_in_outputs = vars_in_variables(vars_input, vars_output) + vars_all_outputs = vars_output + if haskey(timestep_mapped_outputs, proc_output) + vars_all_outputs = (; vars_output..., extract_mapped_outputs(timestep_mapped_outputs[proc_output])...) + end + + inputs_in_outputs = vars_in_variables(vars_input, vars_all_outputs) if any(inputs_in_outputs) ins_in_outs = [vars_input...][inputs_in_outputs] @@ -511,13 +559,20 @@ function search_inputs_in_multiscale_output(process, organ, inputs, soft_dep_gra return inputs_as_output_of_other_scale end - function add_input_as_output!(inputs_as_output_of_other_scale, soft_dep_graphs, organ_source, variable, value) + + timestep_mapped_outputs = soft_dep_graphs[organ_source][:timestep_mapped_outputs] for (proc_output, pairs_vars_output) in soft_dep_graphs[organ_source][:outputs] # e.g. proc_output = :maintenance_respiration; pairs_vars_output = soft_dep_graphs_roots.roots[organ_source][:outputs][proc_output] vars_output = flatten_vars(pairs_vars_output) + + vars_all_outputs = vars_output + + if haskey(timestep_mapped_outputs, proc_output) + vars_all_outputs = (; vars_output..., extract_mapped_outputs(timestep_mapped_outputs[proc_output])...) + end # If the variable is found in the outputs of the process at the other scale: - if variable in keys(vars_output) + if variable in keys(vars_all_outputs) # The variable is found at another scale: if haskey(inputs_as_output_of_other_scale, organ_source) if haskey(inputs_as_output_of_other_scale[organ_source], proc_output) diff --git a/src/mtg/GraphSimulation.jl b/src/mtg/GraphSimulation.jl index 93376240..f0fb652b 100644 --- a/src/mtg/GraphSimulation.jl +++ b/src/mtg/GraphSimulation.jl @@ -14,6 +14,7 @@ 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} @@ -26,10 +27,12 @@ struct GraphSimulation{T,S,U,O,V} models::Dict{String,U} outputs::Dict{String,O} outputs_index::Dict{String, Int} + orchestrator::Orchestrator + 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, orchestrator=Orchestrator()) + 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 @@ -144,4 +147,20 @@ end function convert_outputs(out::TimeStepTable{T} where T, sink) @assert Tables.istable(sink) "The sink argument must be compatible with the Tables.jl interface (`Tables.istable(sink)` must return `true`, *e.g.* `DataFrame`)" return sink(out) +end + + +struct MultiScaleMapping{T} + default_timestep::Date + mapping::Dict{String,T} +end + + +function MultiScaleMapping(mapping, mtg; nsteps=1, outputs=nothing, type_promotion=nothing, check=true, verbose=false, orchestrator=Orchestrator()) + GraphSimulation(init_simulation(mtg, mapping; nsteps=nsteps, outputs=outputs, type_promotion=type_promotion, check=check, verbose=verbose, orchestrator=orchestrator)...) +end + + +function dep(m::MultiScaleMapping) + return dep(m.graph) end \ No newline at end of file diff --git a/src/mtg/MultiScaleModel.jl b/src/mtg/MultiScaleModel.jl index 3243afc6..439ff53b 100644 --- a/src/mtg/MultiScaleModel.jl +++ b/src/mtg/MultiScaleModel.jl @@ -102,11 +102,26 @@ julia> PlantSimEngine.model_(multiscale_model) ToyCAllocationModel() ``` """ -struct MultiScaleModel{T<:AbstractModel,V<:AbstractVector{Pair{A,Union{Pair{S,Symbol},Vector{Pair{S,Symbol}}}}} where {A<:Union{Symbol,PreviousTimeStep},S<:AbstractString}} + +# timestep mapped variables : original variable name, aggregated variable name, timestep at which it operates, aggregation function +struct TimestepMappedVariable + name_from::Symbol + name_to::Symbol + timestep_to::Period + aggregation_function +end + +#=function TimestepMappedVariable(name::Symbol, name_to::Symbol, ts_to::Date, aggreg_fn) + TimestepMappedVariable(name, name_to, ts_to, aggreg_fn) +end=# + + +struct MultiScaleModel{T<:AbstractModel,V<:AbstractVector{Pair{A,Union{Pair{S,Symbol},Vector{Pair{S,Symbol}}}}} where {A<:Union{Symbol,PreviousTimeStep},S<:AbstractString}, W<:AbstractVector{TimestepMappedVariable}} model::T mapped_variables::V + timestep_mapped_variables::W - function MultiScaleModel{T}(model::T, mapped_variables) where {T<:AbstractModel} + function MultiScaleModel{T}(model::T, mapped_variables, timestep_mapped_variables) where {T<:AbstractModel} # Check that the variables in the mapping are variables of the model: model_variables = keys(variables(model)) for i in mapped_variables @@ -117,10 +132,28 @@ struct MultiScaleModel{T<:AbstractModel,V<:AbstractVector{Pair{A,Union{Pair{S,Sy var = isa(var, PreviousTimeStep) ? var.variable : var if !(var in model_variables) + # TODO check duplicates within the timestep_mapped_variables error("Mapping for model $model defines variable $var, but it is not a variable of the model.") end end + for i in timestep_mapped_variables + # TODO handle cases where name_from is a PreviousTimestep in the mapped variables + # var = isa(i, PreviousTimeStep) ? i.variable : first(i) + var = i.name_from + + # TODO ensure no conflicts with other mappings and refs on that variable... ? + if !(var in model_variables) + error("Timestep mapping for model $model requires variable $var, but it is not a variable of the model.") + end + + # Avoid name conflicts # TODO make sure no model at that scale causes name conflicts by having the same name (amongst its outputs), not just the model owning that variable + var_out = i.name_to + if var_out in model_variables + error("Timestep mapping for model $model defines an output variable $var, but that name is already used as a variable in the model.") + end + end + # If the name of the variable mapped from the other scale is not given, we add it as the same of the variable name in the model. Cases: # 1. `[:variable_name => "Plant"]` # We take one value from the Plant node # 2. `[:variable_name => ["Leaf"]]` # We take a vector of values from the Leaf nodes @@ -137,8 +170,12 @@ struct MultiScaleModel{T<:AbstractModel,V<:AbstractVector{Pair{A,Union{Pair{S,Sy push!(unfolded_mapping, _get_var(isa(i, PreviousTimeStep) ? i : Pair(i.first, i.second), process_)) # Note: We are using Pair(i.first, i.second) to make sure the Pair is specialized enough, because sometimes the vector in the mapping made the Pair not specialized enough e.g. [:v1 => "S" => :v2,:v3 => "S"] makes the pairs `Pair{Symbol, Any}`. end - - new{T,typeof(unfolded_mapping)}(model, unfolded_mapping) + + for i in timestep_mapped_variables + #TODO + end + + new{T,typeof(unfolded_mapping), typeof(timestep_mapped_variables)}(model, unfolded_mapping, timestep_mapped_variables) end end @@ -185,16 +222,25 @@ end -function MultiScaleModel(model::T, mapped_variables) where {T<:AbstractModel} - MultiScaleModel{T}(model, mapped_variables) +function MultiScaleModel(model::T, mapped_variables, timestep_mapped_variables) where {T<:AbstractModel} + MultiScaleModel{T}(model, mapped_variables, timestep_mapped_variables) end -MultiScaleModel(; model, mapped_variables) = MultiScaleModel(model, mapped_variables) + +MultiScaleModel(; model, mapped_variables, timestep_mapped_variables=TimestepMappedVariable[]) = MultiScaleModel(model, mapped_variables, timestep_mapped_variables) mapped_variables_(m::MultiScaleModel) = m.mapped_variables +#timestep_mapped_variables_(m::MultiScaleModel) = m.timestep_mapped_variables model_(m::MultiScaleModel) = m.model inputs_(m::MultiScaleModel) = inputs_(m.model) outputs_(m::MultiScaleModel) = outputs_(m.model) +function timestep_mapped_outputs_(m::MultiScaleModel) + # TODO outputs_(m.model)[i.name_from] is the default value of the source variable + # this is not going to be correct beyond simple examples + [i.name_to => Pair(i, outputs_(m.model)[i.name_from]) for i in m.timestep_mapped_variables] +end get_models(m::MultiScaleModel) = [model_(m)] # Get the models of a MultiScaleModel: # Note: it is returning a vector of models, because in this case the user provided a single MultiScaleModel instead of a vector of. get_status(m::MultiScaleModel) = nothing -get_mapped_variables(m::MultiScaleModel{T,S}) where {T,S} = mapped_variables_(m) \ No newline at end of file +function get_mapped_variables(m::MultiScaleModel{T,S}) where {T,S} + mapped_variables_(m) +end \ No newline at end of file diff --git a/src/mtg/add_organ.jl b/src/mtg/add_organ.jl index 1c2925ef..1147056a 100644 --- a/src/mtg/add_organ.jl +++ b/src/mtg/add_organ.jl @@ -33,5 +33,9 @@ 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 RefVector udpating with node insertion + # 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 57aef81c..51de0dce 100644 --- a/src/mtg/initialisation.jl +++ b/src/mtg/initialisation.jl @@ -21,9 +21,9 @@ 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=Orchestrator()))=#; type_promotion=nothing, verbose=false, check=true, orchestrator=Orchestrator()) # We compute the variables mapping for each scale: - mapped_vars = mapped_variables(mapping, dependency_graph, verbose=verbose) + mapped_vars = mapped_variables(mapping, dependency_graph, verbose=verbose,orchestrator=orchestrator) # Update the types of the variables as desired by the user: convert_vars!(mapped_vars, type_promotion) @@ -35,18 +35,31 @@ 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 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, + # 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 + # 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. # 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) @@ -91,7 +104,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=Orchestrator()) # Check if the node has a model defined for its symbol, if not, no need to compute symbol(node) ∉ collect(keys(mapped_vars)) && return @@ -146,7 +159,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) @@ -169,7 +182,7 @@ function init_node_status!(node, statuses, mapped_vars, reverse_multiscale_mappi end """ - status_from_template(d::Dict{Symbol,Any}) + status_from_template(d::Dict{Symbol,Any}, scale::String) Create a status from a template dictionary of variables and values. If the values are already RefValues or RefVectors, they are used as is, else they are converted to Refs. @@ -189,7 +202,7 @@ julia> using PlantSimEngine ``` ```jldoctest mylabel -julia> a, b = PlantSimEngine.status_from_template(Dict(:a => 1.0, :b => 2.0)); +julia> a, b = PlantSimEngine.status_from_template(Dict(:a => 1.0, :b => 2.0), "Dummy_Scale"); ``` ```jldoctest mylabel @@ -202,7 +215,7 @@ 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::Orchestrator=Orchestrator()) # 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. @@ -305,7 +318,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 @@ -314,9 +327,18 @@ 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) + + soft_dep_graph_roots, hard_dep_dict = hard_dependencies(mapping; verbose=false, orchestrator=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, soft_dep_graph_roots; type_promotion=type_promotion, verbose=verbose, check=check, orchestrator=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, soft_dep_graph_roots; 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 +351,125 @@ 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) + #dependency_graph = dep(mapping, soft_dep_graph_roots=soft_dep_graph_roots, hard_dep_dict=hard_dep_dict, orchestrator=orchestrator) + + # 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. + dependency_graph = soft_dependencies_multiscale(soft_dep_graph_roots, reverse_multiscale_mapping, hard_dep_dict) + # 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. + + # Checking that the graph is acyclic: + iscyclic, cycle_vec = is_graph_cyclic(dependency_graph; warn=false) + # Note: we could do that in `soft_dependencies_multiscale` but we prefer to keep the function as simple as possible, and + # usable on its own. + + iscyclic && error("Cyclic dependency detected in the graph. Cycle: \n $(print_cycle(cycle_vec)) \n You can break the cycle using the `PreviousTimeStep` variable in the mapping.") + # Third step, we identify which + + # 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 tmst in orchestrator.non_default_timestep_mapping + 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) + 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 21d64008..94130a3b 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; verbose=false, orchestrator=Orchestrator()) # Initialise a dict that defines the multiscale variables for each organ type: mapped_vars = mapped_variables_no_outputs_from_other_scale(mapping, dependency_graph) @@ -22,7 +22,7 @@ function mapped_variables(mapping, dependency_graph=first(hard_dependencies(mapp # Find variables that are inputs to other scales as a `SingleNodeMapping` and declare them as MappedVar from themselves in the source scale. # This helps us declare it as a reference when we create the template status objects. - transform_single_node_mapped_variables_as_self_node_output!(mapped_vars) + transform_single_node_mapped_variables_as_self_node_output!(mapped_vars, orchestrator) # We now merge inputs and outputs into a single dictionary: mapped_vars_per_organ = merge(merge, mapped_vars[:inputs], mapped_vars[:outputs]) @@ -54,10 +54,24 @@ 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))) - 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) + +function combine_outputs(organ, outs, timestep_mapped_outputs_process) + if length(timestep_mapped_outputs_process) > 0 + timestep_mapped_organ = NamedTuple() + for (proc, data) in timestep_mapped_outputs_process + timestep_mapped_organ = (; timestep_mapped_organ..., extract_mapped_outputs(data)...) + end + #timestep_mapped_organ = collect(values(timestep_mapped_outputs_process))[1] + (; flatten_vars(vcat(values(outs)...))..., timestep_mapped_organ...) + else + flatten_vars(vcat(values(outs)...)) + end +end + +function mapped_variables_no_outputs_from_other_scale(mapping, dependency_graph) + nodes_insouts = Dict(organ => (inputs=ins, outputs=outs, timestep_outs=timestep_mapped_outs) for (organ, (soft_dep_graph, ins, outs, timestep_mapped_outs)) in dependency_graph.roots) + ins = Dict{String,NamedTuple}(organ => flatten_vars(vcat(values(ins)...)) for (organ, (ins, outs, timestep_mapped_outs)) in nodes_insouts) + outs = Dict{String,NamedTuple}(organ => combine_outputs(organ, outs, timestep_mapped_outs) for (organ, (ins, outs, timestep_mapped_outs)) in nodes_insouts) return Dict(:inputs => ins, :outputs => outs) end @@ -159,18 +173,21 @@ This helps us declare it as a reference when we create the template status objec These node are found in the mapping as `[:variable_name => "Plant"]` (notice that "Plant" is a scalar value). """ -function transform_single_node_mapped_variables_as_self_node_output!(mapped_vars) +function transform_single_node_mapped_variables_as_self_node_output!(mapped_vars, orchestrator=Orchestrator()) for (organ, vars) in mapped_vars[:inputs] # e.g. organ = "Leaf"; vars = mapped_vars[:inputs][organ] for (var, mapped_var) in pairs(vars) # e.g. var = :carbon_biomass; mapped_var = vars[var] if isa(mapped_var, MappedVar{SingleNodeMapping}) source_organ = mapped_organ(mapped_var) source_organ == "" && continue # We skip the variables that are mapped to themselves (e.g. [PreviousTimeStep(:variable_name)], or just renaming a variable) + # TODO dirty prototyping to see how timestep mapped variables work at the same scale + # not good to special-case them + #source_organ == organ && continue @assert source_organ != organ "Variable `$var` is mapped to its own scale in organ $organ. This is not allowed." @assert haskey(mapped_vars[:outputs], source_organ) "Scale $source_organ not found in the mapping, but mapped to the $organ scale." + @assert haskey(mapped_vars[:outputs][source_organ], source_variable(mapped_var)) "The variable `$(source_variable(mapped_var))` is mapped from scale `$source_organ` to " * "scale `$organ`, but is not computed by any model at `$source_organ` scale." - # If the source variable was already defined as a `MappedVar{SelfNodeMapping}` by another scale, we skip it: isa(mapped_vars[:outputs][source_organ][source_variable(mapped_var)], MappedVar{SelfNodeMapping}) && continue # Note: this happens when a variable is mapped to several scales, e.g. soil_water_content computed at soil scale can be @@ -186,6 +203,7 @@ function transform_single_node_mapped_variables_as_self_node_output!(mapped_vars mapped_vars[:outputs][source_organ][source_variable(mapped_var)], ) ) + # end mapped_vars[:outputs][source_organ] = merge(mapped_vars[:outputs][source_organ], self_mapped_var) # Note: merge overwrites the LHS values with the RHS values if they have the same key. end @@ -286,7 +304,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::Orchestrator) # 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 +321,7 @@ 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]))) + push!(dict_mapped_vars, key => Ref(mapped_default(vars[k]))) 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/mapping/reverse_mapping.jl b/src/mtg/mapping/reverse_mapping.jl index 059b4294..7b07135e 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=Orchestrator())), verbose=false) reverse_mapping(mapped_vars, all=all) end diff --git a/src/mtg/save_results.jl b/src/mtg/save_results.jl index 51baffcb..1845bcf1 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=Orchestrator()) 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/processes/model_initialisation.jl b/src/processes/model_initialisation.jl index 95a8cc76..d944b16c 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=Orchestrator())), 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 d3312473..127fa323 100644 --- a/src/run.jl +++ b/src/run.jl @@ -368,14 +368,15 @@ function run!( nsteps=nothing, tracked_outputs=nothing, check=true, - executor=ThreadedEx() + executor=ThreadedEx(), + orchestrator::Orchestrator=Orchestrator(), ) 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, orchestrator=orchestrator) run!( sim, meteo_adjusted, @@ -448,28 +449,85 @@ function run_node_multiscale!( executor ) where {T<:GraphSimulation} # T is the status of each node by organ type - # 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 + # Hack until default timestep is inferred from the meteo + node_ratio = node.timestep / Day(1) - node_statuses = status(object)[node.scale] # Get the status of the nodes at the current scale - models_at_scale = models[node.scale] + # Check if the model needs to run this timestep - 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 + skip = (1 + (i - 1) % node_ratio) != node_ratio - node.simulation_id[1] += 1 # increment the simulation id, to remember that the model has been called already + if skip - # 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. + # TODO : This prevents a non-default timestep node form being updated by two different parents + # but probably should be changed, it is bug-prone + if isnothing(node.parent) || any([p.simulation_id[1] > node.simulation_id[1] for p in node.parent]) + node.simulation_id[1] += 1 + end + else + + # 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 + + 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 actually be checkde pre-simulation, it is incorrect behaviour + if isnothing(node.timestep_mapping_data) || isempty(node.timestep_mapping_data) + # Samuel : this is the happy path, no further timestep mapping checks needed + + 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 remove this bit + # 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) + + # Do the accumulation then write into the child's status if need be + for tmst in node.timestep_mapping_data + ratio = Int(tmst.timestep_to / 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] + st[tmst.variable_to] = tmst.mapping_function(tmst.mapping_data[node_id(st.node)]) + # 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 + 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 diff --git a/src/timestep/timestep_mapping.jl b/src/timestep/timestep_mapping.jl new file mode 100644 index 00000000..d9ac170d --- /dev/null +++ b/src/timestep/timestep_mapping.jl @@ -0,0 +1,80 @@ + +# 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 + +# Many shortcuts will be taken, I'll try and comment what's missing/implicit/etc. + +# 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 + +# 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 +end + +struct Var_to + name::Symbol +end + +struct ModelTimestepMapping + model + scale::String + timestep::Period +end + +mutable struct Orchestrator + default_timestep::Period + non_default_timestep_mapping::Vector{ModelTimestepMapping} + + function Orchestrator(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 + +Orchestrator() = Orchestrator(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..e4563527 --- /dev/null +++ b/test/test-multitimestep.jl @@ -0,0 +1,1231 @@ + +########################### +# Test with three timesteps, multiscale +########################### + +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 "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 + +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), ) + + to_w = PlantSimEngine.Var_to(:in_week) + from_d = PlantSimEngine.Var_from(MyToyDay2Model, "Default", :out_day, sum) + 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(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) + + orch2 = PlantSimEngine.Orchestrator(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, + ], + ),), + #="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)) +#mtg4 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default4", 1, 4)) +#mtg5 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default5", 1, 5)) + + #orch2 = PlantSimEngine.Orchestrator() + +#out = @run run!(mtg, m_multiscale, df, orchestrator=orch2) +out = run!(mtg, m_multiscale, df, orchestrator=orch2) + + + +using Test + @test unique([out["Default3"][i].in_four_week_from_day for i in 1:length(out["Default3"])]) == [-Inf, 28.0] + @test unique([out["Default3"][i].in_four_week_from_week for i in 1:length(out["Default3"])]) == [-Inf, 28.0] + + # Note : until the models actually run, inputs_agreement defaults to false, so it's only expected to be true + # from day 28 onwards + @test unique([out["Default3"][i].inputs_agreement for i in 28:length(out["Default3"])]) == [1] + + ########################### +# Three timestep model that is single-scale, to circumvent refvector/refvalue overwriting +# (eg filtering out timestep-mapped variables from vars_need_init and storing the values elsewhere) +# and check mapping at the same scale +########################### + + m_singlescale = Dict("Default" => ( + MyToyDay2Model(), + MyToyWeek2Model(), + MyToyFourWeek2Model(), + ),) + + + model_timesteps_defaultscale = Dict(MyToyWeek2Model =>Week(1), MyToyFourWeek2Model =>Week(4), ) + to_w = PlantSimEngine.Var_to(:in_week) + from_d = PlantSimEngine.Var_from(MyToyDay2Model, "Default", :out_day, sum) + 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, "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.Orchestrator(Day(1), [mtsm_w, mtsm_w4]) + + mtg_single = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) + out = run!(mtg_single, m_singlescale, df, orchestrator=orch2) + + + +########################### +# Test with three timesteps, multiscale + previoustimestep +########################### + +# note the daily models don't specify a timestep range +# if you copy-paste this elsewhere, bear in mind that you might need to specify it + +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 "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 + +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 "ToyPreviousWeek2" verbose = false + +struct MyToyPreviousWeek2Model <: AbstractToypreviousweek2Model end + +# TODO initialisation issue +PlantSimEngine.inputs_(::MyToyPreviousWeek2Model) = (in_last_week=-Inf,) +PlantSimEngine.outputs_(m::MyToyPreviousWeek2Model) = (out_last_week=-Inf,) + +function PlantSimEngine.run!(m::MyToyPreviousWeek2Model, models, status, meteo, constants=nothing, extra=nothing) + status.out_last_week += status.in_last_week/7.0 +end + +PlantSimEngine.timestep_range_(m::MyToyPreviousWeek2Model) = 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), ) + + to_w = PlantSimEngine.Var_to(:in_week) + from_d = PlantSimEngine.Var_from(MyToyDay2Model, "Default", :out_day, sum) + 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(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) + + orch2 = PlantSimEngine.Orchestrator(Day(1), [mtsm_w, mtsm_w4]) + + +m_multiscale = Dict( +"Default6" => + ( + MultiScaleModel(model=MyToyPreviousWeek2Model(), + mapped_variables=[PlantSimEngine.PreviousTimeStep(:in_last_week) => "Default2" => :out_week], + ), + Status(in_last_week=0.0, out_last_week=0.0) + ), +"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, + ], + ),), + "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)) +#mtg4 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default4", 1, 4)) +#mtg5 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default5", 1, 5)) +mtg6 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default6", 1, 6)) + +out = @run run!(mtg, m_multiscale, df, orchestrator=orch2) +out = run!(mtg, m_multiscale, df, orchestrator=orch2) + +unique!([out["Default6"][i].out_last_week for i in 1:length(out["Default6"])]) + +# TODO : out_last_week is an output at the weekly scale, it isn't mapped +# Doesn't mesh well with current implementation : +# non_default_timestep_mapping requires input + output, so the model and variable aren't declared +# I can infer the timestep from the timestep_range in this situation, but not for a model with a wider range +# This means the model needs to be declared somewhere + + +########################### +# Test with one timestep, multiscale + previoustimestep +########################### + +# note the daily models don't specify a timestep range +# if you copy-paste this elsewhere, bear in mind that you might need to specify it + +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(Day(1)) + +PlantSimEngine.@process "ToyPreviousWeek2" verbose = false + +struct MyToyPreviousWeek2Model <: AbstractToypreviousweek2Model end + +# TODO initialisation issue +PlantSimEngine.inputs_(::MyToyPreviousWeek2Model) = (in_last_week=-Inf,) +PlantSimEngine.outputs_(m::MyToyPreviousWeek2Model) = (out_last_week=-Inf,) + +function PlantSimEngine.run!(m::MyToyPreviousWeek2Model, models, status, meteo, constants=nothing, extra=nothing) + status.out_last_week += status.in_last_week +end + +PlantSimEngine.timestep_range_(m::MyToyPreviousWeek2Model) = TimestepRange(Day(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(Day(1)) + + + +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), ) + + to_w = PlantSimEngine.Var_to(:in_week) + from_d = PlantSimEngine.Var_from(MyToyDay2Model, "Default", :out_day, sum) + 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(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) + + orch2 = PlantSimEngine.Orchestrator(Day(1), [mtsm_w, mtsm_w4])=# +orch2 = PlantSimEngine.Orchestrator() + +m_multiscale = Dict( +"Default6" => + ( + MultiScaleModel(model=MyToyPreviousWeek2Model(), + mapped_variables=[PlantSimEngine.PreviousTimeStep(:in_last_week) => "Default2" => :out_week], + ), + Status(in_last_week=0.0, out_last_week=0.0) + ), +"Default" => ( + MyToyDay2Model(), + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeek2Model(), + mapped_variables=[:in_week => "Default" => :out_day], + ), + ), + + ) + + +# TODO test with multiple nodes +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) +mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 2)) +mtg6 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default6", 1, 6)) + +out = @run run!(mtg, m_multiscale, df, orchestrator=orch2) +out = run!(mtg, m_multiscale, df, orchestrator=orch2) + +unique!([out["Default6"][i].out_last_week for i in 1:length(out["Default6"])]) + +########################### +# Previous timestep debugging, not useful for testing timestep mapping atm +########################### +#= +using MultiScaleTreeGraph +using PlantSimEngine +using PlantMeteo +using PlantMeteo.Dates + +PlantSimEngine.@process "current_timestep" verbose = false + +struct HelperCurrentTimestepModel <: AbstractCurrent_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.@process "next_timestep" verbose = false + struct HelperNextTimestepModel <: AbstractNext_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 + + df = DataFrame(:data => [1 for i in 1:365], ) + + m_ms = Dict( + "B" => ( + MultiScaleModel( + model=HelperCurrentTimestepModel(), + mapped_variables=[PreviousTimeStep(:next_timestep),], + ), + Status(next_timestep=2) + + ), + "A" => ( + HelperNextTimestepModel(), + Status(current_timestep=1) + ), + ) + +m_ss = Dict( + "A" => ( + HelperNextTimestepModel(), + MultiScaleModel( + model=HelperCurrentTimestepModel(), + mapped_variables=[PreviousTimeStep(:next_timestep),], + ), + Status(current_timestep=1, next_timestep=1)), +) + + mtg_ = Node(MultiScaleTreeGraph.NodeMTG("/", "A", 1, 1)) +#mtg2 = Node(mtg_, MultiScaleTreeGraph.NodeMTG("/", "B", 1, 2)) +out = @run run!(mtg_, m_ss, df) +=# + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +################################ +## API change : integrate timestep mapping into multiscalemodels +################################ + + +########################### +# Simple test with an orchestrator +########################### +using MultiScaleTreeGraph +using PlantSimEngine +using PlantMeteo +using PlantMeteo.Dates +using Test + +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(15.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) + +m_multiscale = Dict("Default" => ( + MultiScaleModel(model=MyToyDayModel(), + mapped_variables=[], + timestep_mapped_variables=[TimestepMappedVariable(:daily_temperature, :weekly_max_temperature, Week(1), maximum),] + ), + Status(a=1,) + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeekModel(), + #mapped_variables=[:weekly_max_temperature => ["Default" => :daily_temperature]], # TODO test this + mapped_variables=[:weekly_max_temperature => "Default" => :weekly_max_temperature], + timestep_mapped_variables=PlantSimEngine.TimestepMappedVariable[], #TODO avoid this + ), + ),) + + +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) +mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 2)) + +mtsm = PlantSimEngine.ModelTimestepMapping(MyToyWeekModel, "Default2", Week(1)) + +orch2 = PlantSimEngine.Orchestrator(Day(1), [mtsm,]) + +#out = @run run!(mtg, m_multiscale, meteo_day, orchestrator=orch2) +out = run!(mtg, m_multiscale, meteo_day, orchestrator=orch2) + +temps = [out["Default"][i].daily_temperature for i in 1:365] +temp_m = maximum(temps) + +# At least one week should have max temp > 28 +@test temp_m > 28 && unique!([out["Default2"][i].hot for i in 1:365]) == [0,1] + + +########################### +# Test with three timesteps, multiscale +########################### + +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], ) + + mtsm_w = PlantSimEngine.ModelTimestepMapping(MyToyWeek2Model, "Default2", Week(1)) + mtsm_w4 = PlantSimEngine.ModelTimestepMapping(MyToyFourWeek2Model, "Default3", Week(4)) + + orch2 = PlantSimEngine.Orchestrator(Day(1), [mtsm_w, mtsm_w4]) + + +m_multiscale = Dict("Default" => ( + MultiScaleModel(model=MyToyDay2Model(), + mapped_variables=[], + timestep_mapped_variables=[TimestepMappedVariable(:out_day, :out_week_from_day, Week(1), sum), + TimestepMappedVariable(:out_day, :out_four_week_from_day, Week(4), sum),] + ), + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeek2Model(), + mapped_variables=[:in_week => "Default" => :out_week_from_day], + timestep_mapped_variables=[TimestepMappedVariable(:out_week, :out_four_week_from_week, Week(4), sum),] + ), + ), + "Default3" => ( + MultiScaleModel(model=MyToyFourWeek2Model(), + mapped_variables=[ + :in_four_week_from_day => "Default" => :out_four_week_from_day, + :in_four_week_from_week => "Default2" => :out_four_week_from_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 = run!(mtg, m_multiscale, df, orchestrator=orch2) + + + +using Test + @test unique([out["Default3"][i].in_four_week_from_day for i in 1:length(out["Default3"])]) == [-Inf, 28.0] + @test unique([out["Default3"][i].in_four_week_from_week for i in 1:length(out["Default3"])]) == [-Inf, 28.0] + + # Note : until the models actually run, inputs_agreement defaults to false, so it's only expected to be true + # from day 28 onwards + @test unique([out["Default3"][i].inputs_agreement for i in 28:length(out["Default3"])]) == [1] + +########################### +# Three timestep model that is single-scale, to circumvent refvector/refvalue overwriting +# (eg filtering out timestep-mapped variables from vars_need_init and storing the values elsewhere) +# and check mapping at the same scale +########################### + +# This example has variable renaming at the same scale + + m_singlescale_mapped = Dict("Default" => ( + MultiScaleModel(model=MyToyDay2Model(), + mapped_variables=[], + timestep_mapped_variables=[TimestepMappedVariable(:out_day, :out_week_from_day, Week(1), sum), + TimestepMappedVariable(:out_day, :out_four_week_from_day, Week(4), sum),] + ), + MultiScaleModel(model=MyToyWeek2Model(), + mapped_variables=[:in_week => "Default" => :out_week_from_day], + timestep_mapped_variables=[TimestepMappedVariable(:out_week, :out_four_week_from_week, Week(4), sum),] + ), + MultiScaleModel(model=MyToyFourWeek2Model(), + mapped_variables=[ + :in_four_week_from_day => "Default" => :out_four_week_from_day, + :in_four_week_from_week => "Default" => :out_four_week_from_week, + ], + ),)) + + # This one reuses the variable names directly, so requires only timestep mapping + m_singlescale = Dict("Default" => ( + MultiScaleModel(model=MyToyDay2Model(), + mapped_variables=[], + timestep_mapped_variables=[TimestepMappedVariable(:out_day, :in_week, Week(1), sum), + TimestepMappedVariable(:out_day, :in_four_week_from_day, Week(4), sum),] + ), + MultiScaleModel(model=MyToyWeek2Model(), + mapped_variables=[], + timestep_mapped_variables=[TimestepMappedVariable(:out_week, :in_four_week_from_week, Week(4), sum),] + ), + MyToyFourWeek2Model(), + )) + + mtsm_w = PlantSimEngine.ModelTimestepMapping(MyToyWeek2Model, "Default", Week(1)) + mtsm_w4 = PlantSimEngine.ModelTimestepMapping(MyToyFourWeek2Model, "Default", Week(4)) + + orch2 = PlantSimEngine.Orchestrator(Day(1), [mtsm_w, mtsm_w4]) + +mtg_single = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) +out = run!(mtg_single, m_singlescale, df, orchestrator=orch2) +out = run!(mtg_single, m_singlescale_mapped, df, orchestrator=orch2) + +########################### +# Test with a D -> W -> D configuration, with multiple variables mapped between timesteps +########################### + +using MultiScaleTreeGraph +using PlantSimEngine +using PlantMeteo +using PlantMeteo.Dates + +PlantSimEngine.@process "ToyDayDWD" verbose = false + +struct MyToyDayDWDModel <: AbstractToydaydwdModel end + +PlantSimEngine.inputs_(m::MyToyDayDWDModel) = (a=1,) +PlantSimEngine.outputs_(m::MyToyDayDWDModel) = (daily_temperature=-Inf,) + +function PlantSimEngine.run!(m::MyToyDayDWDModel, models, status, meteo, constants=nothing, extra=nothing) + status.daily_temperature = meteo.data +end + +PlantSimEngine.@process "ToyWeekDWD" verbose = false + +struct MyToyWeekDWDModel <: AbstractToyweekdwdModel + temperature_threshold::Float64 +end + +MyToyWeekDWDModel() = MyToyWeekDWDModel(30.0) +function PlantSimEngine.inputs_(::MyToyWeekDWDModel) + (weekly_max_temperature=-Inf, weekly_sum_temperature=-Inf) +end +PlantSimEngine.outputs_(m::MyToyWeekDWDModel) = (hot = false, sum=-Inf) + +function PlantSimEngine.run!(m::MyToyWeekDWDModel, models, status, meteo, constants=nothing, extra=nothing) + status.hot = status.weekly_max_temperature > m.temperature_threshold + status.sum += status.weekly_sum_temperature +end + +PlantSimEngine.timestep_range_(m::MyToyWeekDWDModel) = TimestepRange(Week(1)) + +PlantSimEngine.@process "ToyDayDWDOut" verbose = false + +struct MyToyDayDWDOutModel <: AbstractToydaydwdoutModel end + +PlantSimEngine.inputs_(m::MyToyDayDWDOutModel) = (sum=-Inf,weekly_sum_temperature=-Inf,) +PlantSimEngine.outputs_(m::MyToyDayDWDOutModel) = (out=-Inf,) + +function PlantSimEngine.run!(m::MyToyDayDWDOutModel, models, status, meteo, constants=nothing, extra=nothing) + status.out = status.sum - status.weekly_sum_temperature +end + +df = DataFrame(:data => [1 for i in 1:365], ) + +# TODO check that DWDOUT properly uses the variables from Default2 and not Default +m_dwd = Dict("Default" => ( + MultiScaleModel( + model=MyToyDayDWDModel(), + mapped_variables=[], + timestep_mapped_variables=[TimestepMappedVariable(:daily_temperature, :weekly_max_temperature, Week(1), maximum), + TimestepMappedVariable(:daily_temperature, :weekly_sum_temperature, Week(1), sum), + ] ), + MultiScaleModel( + model=MyToyDayDWDOutModel(), + mapped_variables=[:sum => "Default2",]# :weekly_sum_temperature => "Default2"] + ), + #MyToyDayDWDOutModel(), + Status(a=1,out=0.0) + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeekDWDModel(), + #mapped_variables=[:weekly_max_temperature => ["Default" => :daily_temperature]], # TODO test this + mapped_variables=[:weekly_max_temperature => "Default", :weekly_sum_temperature => "Default"], + ), + Status(weekly_max_temperature=0.0, weekly_sum_temperature=0.0, sum=0.0) + ), +) + + +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) +mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 2)) + +mtsm_dwd = PlantSimEngine.ModelTimestepMapping(MyToyWeekDWDModel, "Default2", Week(1)) + +orch_dwd = PlantSimEngine.Orchestrator(Day(1), [mtsm_dwd,])#mtsm2]) + +out = @run run!(mtg, m_dwd, df, orchestrator=orch_dwd) +out = run!(mtg, m_dwd, df, orchestrator=orch_dwd) + + +################################## +# Two variables mapped +################################## + +using MultiScaleTreeGraph +using PlantSimEngine +using PlantMeteo +using PlantMeteo.Dates + +PlantSimEngine.@process "ToyDayDWD" verbose = false + +struct MyToyDayDWDModel <: AbstractToydaydwdModel end + +PlantSimEngine.inputs_(m::MyToyDayDWDModel) = (a=1,) +PlantSimEngine.outputs_(m::MyToyDayDWDModel) = (daily_temperature=-Inf,) + +function PlantSimEngine.run!(m::MyToyDayDWDModel, models, status, meteo, constants=nothing, extra=nothing) + status.daily_temperature = meteo.T +end + +PlantSimEngine.@process "ToyWeekDWD" verbose = false + +struct MyToyWeekDWDModel <: AbstractToyweekdwdModel + temperature_threshold::Float64 +end + +MyToyWeekDWDModel() = MyToyWeekDWDModel(30.0) +function PlantSimEngine.inputs_(::MyToyWeekDWDModel) + (weekly_max_temperature=-Inf, weekly_sum_temperature=-Inf) +end +PlantSimEngine.outputs_(m::MyToyWeekDWDModel) = (hot = false, sum=-Inf) + +function PlantSimEngine.run!(m::MyToyWeekDWDModel, models, status, meteo, constants=nothing, extra=nothing) + status.hot = status.weekly_max_temperature > m.temperature_threshold + status.sum += status.weekly_sum_temperature +end + +PlantSimEngine.timestep_range_(m::MyToyWeekDWDModel) = TimestepRange(Week(1)) + +#df = DataFrame(:data => [1 for i in 1:365], ) +meteo_day = read_weather(joinpath(pkgdir(PlantSimEngine), "examples/meteo_day.csv"), duration=Day) + +m_dwd = Dict("Default" => ( + MultiScaleModel(model=MyToyDayDWDModel(), + mapped_variables=[], + timestep_mapped_variables=[TimestepMappedVariable(:daily_temperature, :weekly_max_temperature, Week(1), maximum), + TimestepMappedVariable(:daily_temperature, :weekly_sum_temperature, Week(1), sum), + ] + ), + Status(a=1,) + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeekDWDModel(), + #mapped_variables=[:weekly_max_temperature => ["Default" => :daily_temperature]], # TODO test this + mapped_variables=[:weekly_max_temperature => "Default", :weekly_sum_temperature => "Default" ], + ), + Status(weekly_max_temperature=0.0, weekly_sum_temperature=0.0, sum =0.0) + ), +) + + +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default", 1, 1)) +mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 2)) + +mtsm_dwd = PlantSimEngine.ModelTimestepMapping(MyToyWeekDWDModel, "Default2", Week(1)) + +orch_dwd = PlantSimEngine.Orchestrator(Day(1), [mtsm_dwd,]) + +out = run!(mtg, m_dwd, meteo_day, orchestrator=orch_dwd) + +# TODO previous timestep, timestep-mapping to the same variable name + + +#TODO should timestep mapped vars also be part of a model's outputs ? + + + + +#TODO +########################## +# Two models, D -> W, but D has two MTG nodes +########################## + +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) + +m_multiscale = Dict("Default" => ( + MultiScaleModel( + model=MyToyDayModel(), + mapped_variables=[], + timestep_mapped_variables=[TimestepMappedVariable(:daily_temperature, :weekly_temperature, Week(1), maximum)], + ), + Status(a=1,) + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeekModel(), + mapped_variables=[:weekly_max_temperature => "Default" => :weekly_temperature], # TODO test this + #mapped_variables=[:weekly_max_temperature => "Default" => :daily_temperature], + ), + ),) + + +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 1)) +mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("+", "Default", 1, 1)) +mtg3 = Node(mtg, MultiScaleTreeGraph.NodeMTG("+", "Default", 1, 2)) + +mtsm = PlantSimEngine.ModelTimestepMapping(MyToyWeekModel, "Default2", Week(1)) + +orch2 = PlantSimEngine.Orchestrator(Day(1), [mtsm,]) + +out = run!(mtg, m_multiscale, meteo_day, orchestrator=orch2) + + +########################## +# Two models, D -> W, but D has two MTG nodes, and we map as a refvector +########################## + +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) = (refvector = false,) + +function PlantSimEngine.run!(m::MyToyWeekModel, models, status, meteo, constants=nothing, extra=nothing) + status.refvector = status.weekly_max_temperature[1] == status.weekly_max_temperature[2] +end + +PlantSimEngine.timestep_range_(m::MyToyWeekModel) = TimestepRange(Week(1)) + + +meteo_day = read_weather(joinpath(pkgdir(PlantSimEngine), "examples/meteo_day.csv"), duration=Day) + +m_multiscale = Dict("Default" => ( + MultiScaleModel( + model=MyToyDayModel(), + mapped_variables=[], + timestep_mapped_variables=[TimestepMappedVariable(:daily_temperature, :weekly_temperature, Week(1), maximum)], + ), + Status(a=1,) + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeekModel(), + mapped_variables=[:weekly_max_temperature => ["Default" => :weekly_temperature]], # TODO test this + #mapped_variables=[:weekly_max_temperature => "Default" => :daily_temperature], + ), + ),) + + +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 1)) +mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("+", "Default", 1, 1)) +mtg3 = Node(mtg, MultiScaleTreeGraph.NodeMTG("+", "Default", 1, 2)) + +mtsm = PlantSimEngine.ModelTimestepMapping(MyToyWeekModel, "Default2", Week(1)) + +orch2 = PlantSimEngine.Orchestrator(Day(1), [mtsm,]) + +# The RefVector will be in the outputs, so intermediate data is lost for such timestep-mapped variables, and it makes the outputs confusing +out = run!(mtg, m_multiscale, meteo_day, orchestrator=orch2) + +using Test +#@test out["Default2"][1] + + +########################## +# Two models, D -> W, but both D and W have two MTG nodes +########################## + +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 + node_id(status.node) +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) = (refvector = false,) + +function PlantSimEngine.run!(m::MyToyWeekModel, models, status, meteo, constants=nothing, extra=nothing) + status.refvector = status.weekly_max_temperature[1] + 1== status.weekly_max_temperature[2] +end + +PlantSimEngine.timestep_range_(m::MyToyWeekModel) = TimestepRange(Week(1)) + + +meteo_day = read_weather(joinpath(pkgdir(PlantSimEngine), "examples/meteo_day.csv"), duration=Day) + +m_multiscale = Dict("Default" => ( + MultiScaleModel( + model=MyToyDayModel(), + mapped_variables=[], + timestep_mapped_variables=[TimestepMappedVariable(:daily_temperature, :weekly_temperature, Week(1), maximum)], + ), + Status(a=1,) + ), + "Default2" => ( + MultiScaleModel(model=MyToyWeekModel(), + mapped_variables=[:weekly_max_temperature => ["Default" => :weekly_temperature]], # TODO test this + #mapped_variables=[:weekly_max_temperature => "Default" => :daily_temperature], + ), + ),) + + +mtg = Node(MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 1)) +mtg2 = Node(mtg, MultiScaleTreeGraph.NodeMTG("/", "Default2", 1, 1)) +mtg3 = Node(mtg2, MultiScaleTreeGraph.NodeMTG("+", "Default", 1, 1)) +mtg4 = Node(mtg2, MultiScaleTreeGraph.NodeMTG("+", "Default", 1, 2)) + +mtsm = PlantSimEngine.ModelTimestepMapping(MyToyWeekModel, "Default2", Week(1)) + +orch2 = PlantSimEngine.Orchestrator(Day(1), [mtsm,]) + +# The RefVector will be in the outputs, so intermediate data is lost for such timestep-mapped variables +out = run!(mtg, m_multiscale, meteo_day, orchestrator=orch2) + +using Test +#@test out["Default2"][1] \ No newline at end of file