Merge pull request #169 from VirtualPlantLab/graphsim-initialisation-changes

Streamline graphsim initialisation

Author: Samuel-amap

src/dependencies/dependencies.jl

@@ -105,11 +105,15 @@ function dep(mapping::Dict{String,T}; verbose::Bool=true) where {T}
     # 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)
+    
+    mapped_vars = mapped_variables(mapping, soft_dep_graphs_roots, verbose=false)
+    reverse_multiscale_mapping = reverse_mapping(mapped_vars, all=false)
+    
     # Second step, compute the soft-dependency graph between SoftDependencyNodes computed in the first step. To do so, we search the 
     # inputs of each process into the outputs of the other processes, at the same scale, but also between scales. Then we keep only the
     # nodes that have no soft-dependencies, and we set them as root nodes of the soft-dependency graph. The other nodes are set as children
     # of the nodes that they depend on.
-    dep_graph = soft_dependencies_multiscale(soft_dep_graphs_roots, mapping, hard_dep_dict)
+    dep_graph = soft_dependencies_multiscale(soft_dep_graphs_roots, 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.

src/dependencies/soft_dependencies.jl

@@ -139,10 +139,8 @@ 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 (proc, i) in soft_dep_graph
@@ -158,7 +156,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".
 

src/mtg/initialisation.jl

@@ -21,7 +21,7 @@ a dictionary of variables that need to be initialised or computed by other model
 
 `(;statuses, status_templates, reverse_multiscale_mapping, vars_need_init, nodes_with_models)`
 """
-function init_statuses(mtg, mapping, dependency_graph=first(hard_dependencies(mapping; verbose=false)); type_promotion=nothing, verbose=false, check=true)
+function init_statuses(mtg, mapping, dependency_graph; type_promotion=nothing, verbose=false, check=true)
     # We compute the variables mapping for each scale:
     mapped_vars = mapped_variables(mapping, dependency_graph, verbose=verbose)
 
@@ -314,10 +314,32 @@ 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
 
+    soft_dep_graphs_roots, hard_dep_dict = hard_dependencies(mapping; verbose=false)
+
     # 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_graphs_roots; type_promotion=type_promotion, verbose=verbose, check=check)
+
+    # 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.
+    # 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, 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(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 
+    
     # 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)
         model_no_node = join(findall(x -> length(x) == 0, statuses), ", ")
@@ -329,5 +351,5 @@ 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)
+    return (; mtg, statuses, status_templates, reverse_multiscale_mapping, vars_need_init, dependency_graph=dep_graph, models, outputs, outputs_index)
 end

src/run.jl

@@ -399,7 +399,7 @@ function run!(
     executor=ThreadedEx()
 )
 
-    dep_graph = dep(object)
+    dep_graph = object.dependency_graph
     models = get_models(object)
     # st = status(object)
 

test/test-mtg-multiscale-cyclic-dep.jl

@@ -48,7 +48,11 @@ out_vars = Dict(
     @test_throws "Cyclic dependency detected in the graph. Cycle:" dep(mapping_cyclic)
 
     soft_dep_graphs_roots, hard_dep_dict = PlantSimEngine.hard_dependencies(mapping_cyclic)
-    dep_graph = PlantSimEngine.soft_dependencies_multiscale(soft_dep_graphs_roots, mapping_cyclic, hard_dep_dict)
+    
+    mapped_vars = PlantSimEngine.mapped_variables(mapping_cyclic, soft_dep_graphs_roots, verbose=false)
+    reverse_mapping_cyclic = PlantSimEngine.reverse_mapping(mapped_vars, all=false)
+    
+    dep_graph = PlantSimEngine.soft_dependencies_multiscale(soft_dep_graphs_roots, reverse_mapping_cyclic, hard_dep_dict)
     iscyclic, cycle_vec = PlantSimEngine.is_graph_cyclic(dep_graph; warn=false)
 
     @test iscyclic
@@ -97,7 +101,10 @@ end
 
     soft_dep_graphs_roots, hard_dep_dict = PlantSimEngine.hard_dependencies(mapping_nocyclic)
     # soft_dep_graphs_roots.roots["Leaf"].inputs
-    dep_graph = PlantSimEngine.soft_dependencies_multiscale(soft_dep_graphs_roots, mapping_nocyclic, hard_dep_dict)
+    mapped_vars = PlantSimEngine.mapped_variables(mapping_nocyclic, soft_dep_graphs_roots, verbose=false)
+    reverse_mapping_nocyclic = PlantSimEngine.reverse_mapping(mapped_vars, all=false)
+
+    dep_graph = PlantSimEngine.soft_dependencies_multiscale(soft_dep_graphs_roots, reverse_mapping_nocyclic, hard_dep_dict)
     iscyclic, cycle_vec = PlantSimEngine.is_graph_cyclic(dep_graph; warn=false)
 
     @test !iscyclic

test/test-mtg-multiscale.jl

@@ -127,8 +127,10 @@ end
 end
 
 @testset "Status initialisation" begin
-    @test_throws "Variable `carbon_biomass` is not computed by any model, not initialised by the user in the status, and not found in the MTG at scale Internode (checked for MTG node 4)." PlantSimEngine.init_statuses(mtg, mapping_1)
-    organs_statuses, others = PlantSimEngine.init_statuses(mtg_init, mapping_1)
+    # Samuel : internal function, suppressing REPL display errors
+    hard_dep_graph = first(PlantSimEngine.hard_dependencies(mapping_1; verbose=false));
+    @test_throws "Variable `carbon_biomass` is not computed by any model, not initialised by the user in the status, and not found in the MTG at scale Internode (checked for MTG node 4)." PlantSimEngine.init_statuses(mtg, mapping_1, hard_dep_graph)  
+    organs_statuses, others = PlantSimEngine.init_statuses(mtg_init, mapping_1, hard_dep_graph)
 
     @test collect(keys(organs_statuses)) == ["Soil", "Internode", "Plant", "Leaf"]
     # Check that the soil_water_content is linked between the soil and the leaves:
@@ -148,7 +150,7 @@ end
     @test organs_statuses["Leaf"][1][:carbon_demand] == -Inf
 
     # Testing with a different type:
-    organs_statuses, others = PlantSimEngine.init_statuses(mtg_init, mapping_1, type_promotion=Dict(Float64 => Float32, Vector{Float64} => Vector{Float32}))
+    organs_statuses, others = PlantSimEngine.init_statuses(mtg_init, mapping_1, hard_dep_graph, type_promotion=Dict(Float64 => Float32, Vector{Float64} => Vector{Float32}))
 
     @test isa(organs_statuses["Plant"][1][:carbon_assimilation], PlantSimEngine.RefVector{Float32})
     @test isa(organs_statuses["Plant"][1][:carbon_allocation], PlantSimEngine.RefVector{Float32})
@@ -168,7 +170,8 @@ end
     type_promotion = nothing
     nsteps = 2
     dependency_graph = dep(mapping_1)
-    organs_statuses, others = PlantSimEngine.init_statuses(mtg_init, mapping_1, first(PlantSimEngine.hard_dependencies(mapping_1; verbose=false)); type_promotion=type_promotion)
+    hard_dep_graph = first(PlantSimEngine.hard_dependencies(mapping_1; verbose=false))
+    organs_statuses, others = PlantSimEngine.init_statuses(mtg_init, mapping_1, hard_dep_graph; type_promotion=type_promotion)
 
     @test collect(keys(organs_statuses)) == ["Soil", "Internode", "Plant", "Leaf"]
     @test collect(keys(organs_statuses["Soil"][1])) == [:node, :soil_water_content]
@@ -194,7 +197,7 @@ end
     @test PlantSimEngine.to_initialize(mapping_1, mtg) == Dict("Internode" => [:carbon_biomass], "Leaf" => [:carbon_biomass])
     @test PlantSimEngine.to_initialize(mapping_1, mtg_init) == Dict{String,Symbol}()
 
-    statuses, status_templates, reverse_multiscale_mapping, vars_need_init = PlantSimEngine.init_statuses(mtg_init, mapping_1)
+    statuses, status_templates, reverse_multiscale_mapping, vars_need_init = PlantSimEngine.init_statuses(mtg_init, mapping_1, hard_dep_graph)
     @test collect(keys(statuses)) == ["Soil", "Internode", "Plant", "Leaf"]
 
     @test length(statuses["Internode"]) == length(statuses["Leaf"]) == 2