Add MultiScale test with status vectors, and fix issues that popped up... Make the function that checks whether vectors are still in mapping statuses return organ + boolean instead of asserting

Author: Samuel-amap

src/mtg/initialisation.jl

@@ -309,7 +309,10 @@ function init_simulation(mtg, mapping; nsteps=1, outputs=nothing, type_promotion
 
     # Ensure the user called the model generation function to handle vectors passed into a status
     # before we keep going
-    check_statuses_contain_no_remaining_vectors(mapping)
+    (organ_with_vector, no_vectors_found) = (check_statuses_contain_no_remaining_vectors(mapping))
+    if !no_vectors_found
+        @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
 
     # Get the status of each node by node type, pre-initialised considering multi-scale variables:
     statuses, status_templates, reverse_multiscale_mapping, vars_need_init =

src/mtg/mapping/model_generation_from_status_vectors.jl

@@ -49,30 +49,54 @@ PlantSimEngine.TimeStepDependencyTrait(::Type{<:HelperNextTimestepModel}) = Plan
 # Note : User specifies at which level they want the basic timestep model to be inserted at, as well as the meteo length
 function replace_mapping_status_vectors_with_generated_models(mapping_with_vectors_in_status, timestep_model_organ_level, nsteps)
 
+    (organ, check) = check_statuses_contain_no_remaining_vectors(mapping_with_vectors_in_status)
+        if check
+        @warn "No vectors, or types deriving from AbstractVector found in statuses, returning mapping as is."
+        return mapping_with_vectors_in_status
+    end
+
+    # we are now certain a model will be generated, and that the timestep models need to be inserted
     mapping = Dict(organ => models for (organ, models) in mapping_with_vectors_in_status)
     for (organ,models) in mapping
-        for status in models
-            if isa(status, Status)
-                new_status, generated_models = generate_model_from_status_vector_variable(mapping, timestep_scale, status, organ, nsteps)
+        for status in models           
+            if isa(status, Status)                
+                # Generate models and remove corresponding vectors from status
+                new_status, generated_models = generate_model_from_status_vector_variable(mapping, timestep_model_organ_level, status, organ, nsteps)
 
-                if length(generated_models) > 0
-                    if organ == timestep_model_organ_level
-                        mapping[organ] = (
-                            models..., 
-                            generated_models..., 
-                            HelperNextTimestepModel(),
-                            MultiScaleModel(
-                            model=HelperCurrentTimestepModel(),
-                            mapping=[PreviousTimeStep(:next_timestep),],
-                            ),
-                            new_status,)
-                    else
-                        mapping[organ] = (
-                            models..., 
-                            generated_models...,                          
-                            new_status,)
-                    end
+                # Avoid inserting empty named tuples into the mapping
+                models_and_new_status = [model for model in models if !isa(model, Status)]
+                if length(new_status) != 0
+                    models_and_new_status = [models_and_new_status..., new_status]
                 end
+
+                # The timestep models might be inserted elsewhere in the mapping, handle various cases
+                if length(generated_models) > 0
+                    mapping[organ] = (
+                        generated_models...,                          
+                        models_and_new_status...,)
+                end      
+            end         
+        end
+        
+        # insert timestep models wherever they're required
+        if organ == timestep_model_organ_level
+            # mapping at a given level can be a tuple or a single model
+            if isa(mapping[organ], AbstractModel) || isa(mapping[organ], MultiScaleModel)
+                mapping[organ] = (
+                    HelperNextTimestepModel(),
+                    MultiScaleModel(
+                    model=HelperCurrentTimestepModel(),
+                    mapping=[PreviousTimeStep(:next_timestep),],
+                    ),
+                    mapping[organ], )
+            else
+                mapping[organ] = (
+                HelperNextTimestepModel(),
+                MultiScaleModel(
+                model=HelperCurrentTimestepModel(),
+                mapping=[PreviousTimeStep(:next_timestep),],
+                ),
+                mapping[organ]..., )
             end
         end
     end
@@ -98,7 +122,7 @@ function generate_model_from_status_vector_variable(mapping, timestep_scale, sta
     # See the test code in test-mapping.jl : cumsum(meteo_day.TT) returns such a data structure
 
     global generated_models_534f1c161f91bb346feba1a84a55e8251f5ad446 = ()
-    global new_status_534f1c161f91bb346feba1a84a55e8251f5ad446 = NamedTuple()
+    global new_status_534f1c161f91bb346feba1a84a55e8251f5ad446 = Status(NamedTuple())
 
     for symbol in keys(status)
         global value_534f1c161f91bb346feba1a84a55e8251f5ad446 = getproperty(status, symbol)
@@ -137,14 +161,13 @@ function generate_model_from_status_vector_variable(mapping, timestep_scale, sta
             run_decl = "function PlantSimEngine.run!(m::$model_name, models, status, meteo, constants=nothing, extra_args=nothing)\nstatus.$symbol = m.$var_vector[status.current_timestep]\nend\n"
             eval(Meta.parse(run_decl))
 
-            # add name to vector of models
+            model_add_decl = "generated_models_534f1c161f91bb346feba1a84a55e8251f5ad446 = (generated_models_534f1c161f91bb346feba1a84a55e8251f5ad446..., $model_name($var_vector=$value_534f1c161f91bb346feba1a84a55e8251f5ad446),)"
+
+            # if :current_timestep is not in the same scale
             if timestep_scale != organ
-                mapping_decl = "mapping[\"($organ)\"] = MultiScaleModel($process_name($var_vector), mapping=\"($timestep_scale)\" => (:current_timestep,))"
-                eval(Meta.parse(mapping_decl))
-            else
+                model_add_decl = "generated_models_534f1c161f91bb346feba1a84a55e8251f5ad446 = (generated_models_534f1c161f91bb346feba1a84a55e8251f5ad446..., MultiScaleModel(model=$model_name($value_534f1c161f91bb346feba1a84a55e8251f5ad446), mapping=[:current_timestep=>\"$timestep_scale\"],),)"
             end 
 
-        model_add_decl = "generated_models_534f1c161f91bb346feba1a84a55e8251f5ad446 = (generated_models_534f1c161f91bb346feba1a84a55e8251f5ad446..., $model_name($var_vector=$value_534f1c161f91bb346feba1a84a55e8251f5ad446),)"
         eval(Meta.parse(model_add_decl))
         else
             new_status_decl = "new_status_534f1c161f91bb346feba1a84a55e8251f5ad446 = Status(; NamedTuple(new_status_534f1c161f91bb346feba1a84a55e8251f5ad446)..., $symbol=$value_534f1c161f91bb346feba1a84a55e8251f5ad446)"
@@ -231,10 +254,12 @@ function check_statuses_contain_no_remaining_vectors(mapping)
             if isa(status, Status)
                 for symbol in keys(status)
                     value = getproperty(status, symbol)
-                    @assert !isa(value, AbstractVector) "Error : Mapping status at $organ 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."
+                    if isa(value, AbstractVector)
+                        return (organ, false)
+                    end
                 end
             end
         end
     end
-    return true
+    return ("", true)
 end

test/test-mapping.jl

@@ -77,7 +77,9 @@ end
             ),
         )
 
-     @test_throws "call the function generate_models_from_status_vectors" PlantSimEngine.check_statuses_contain_no_remaining_vectors(mapping_with_vector)
+        mtg = import_mtg_example()
+        @test !last(PlantSimEngine.check_statuses_contain_no_remaining_vectors(mapping_with_vector))
+        @test_throws "call the function generate_models_from_status_vectors" PlantSimEngine.GraphSimulation(mtg, mapping_with_vector)
 
      mapping_with_empty_status = Dict(
         "Scale" =>
@@ -87,7 +89,7 @@ end
             ),
         )
 
-     @test PlantSimEngine.check_statuses_contain_no_remaining_vectors(mapping_with_empty_status)
+     @test last(PlantSimEngine.check_statuses_contain_no_remaining_vectors(mapping_with_empty_status))
 end
 
 # simple conversion to a mapping, with a manually written model
@@ -200,9 +202,7 @@ PlantSimEngine.ObjectDependencyTrait(::Type{<:ToyTestDegreeDaysCumulModel}) = Pl
     # TODO when outputs filtering is implemented, can test it with this function
     mtg, mapping, outputs_mapping = PlantSimEngine.modellist_to_mapping(models, st2; nsteps=nsteps, outputs=nothing)
 
-    # This test is now redundant as it is called during GraphSimulation initialisation
-    #@test PlantSimEngine.check_statuses_contain_no_remaining_vectors(mapping)
-    graphsim2 = PlantSimEngine.GraphSimulation(mtg, mapping, nsteps=nsteps, check=true, outputs=outputs_mapping)
+   graphsim2 = PlantSimEngine.GraphSimulation(mtg, mapping, nsteps=nsteps, check=true, outputs=outputs_mapping)
 
     sim2 = run!(graphsim2,
         meteo_day,
@@ -216,4 +216,65 @@ PlantSimEngine.ObjectDependencyTrait(::Type{<:ToyTestDegreeDaysCumulModel}) = Pl
 
 end
 
-#[getproperty(a,i) for i in fieldnames(typeof(a))]
+#[getproperty(a,i) for i in fieldnames(typeof(a))]
+
+
+@testset "Vector in status in a multiscale context" begin
+    meteo_day = CSV.read(joinpath(pkgdir(PlantSimEngine), "examples/meteo_day.csv"), DataFrame, header=18)
+    TT_v = Vector(meteo_day.TT)
+    TT_cu_vec = Vector(cumsum(meteo_day.TT))
+    nsteps = length(meteo_day.TT)
+
+    mapping_with_vector = Dict(
+    
+    "Plant" => (
+        MultiScaleModel(
+            model=ToyCAllocationModel(),
+            mapping=[
+                # inputs
+                :carbon_assimilation => ["Leaf"],
+                :carbon_demand => ["Leaf", "Internode"],
+                # outputs
+                :carbon_allocation => ["Leaf", "Internode"]
+            ],
+        ),
+        MultiScaleModel(
+            model=ToyPlantRmModel(),
+            mapping=[:Rm_organs => ["Leaf" => :Rm, "Internode" => :Rm],],
+        ),
+    ),
+    "Internode" => (
+        ToyCDemandModel(optimal_biomass=10.0, development_duration=200.0),
+        ToyMaintenanceRespirationModel(1.5, 0.06, 25.0, 0.6, 0.004),
+        Status(TT=TT_v, carbon_biomass=1.0)
+    ),
+    "Leaf" => (
+        MultiScaleModel(
+            model=ToyAssimModel(),
+            mapping=[:soil_water_content => "Soil",],
+            # Notice we provide "Soil", not ["Soil"], so a single value is expected here
+        ),
+        ToyCDemandModel(optimal_biomass=10.0, development_duration=200.0),
+        ToyMaintenanceRespirationModel(2.1, 0.06, 25.0, 1.0, 0.025),
+        Status(aPPFD=1300.0, carbon_biomass=2.0, TT=10.0), # TODO try calling the generated TT output through a variable mapping
+    ),
+    "Soil" => (
+        ToySoilWaterModel(),
+    ),
+)
+
+out_multiscale = Dict("Plant" => (:Rm_organs,),)
+mtg = import_mtg_example();
+
+mapping_without_vectors = PlantSimEngine.replace_mapping_status_vectors_with_generated_models(mapping_with_vector, "Soil", nsteps)
+
+ graph_sim_multiscale = @test_nowarn PlantSimEngine.GraphSimulation(mtg, mapping_without_vectors, nsteps=nsteps, check=true, outputs=out_multiscale)
+
+    sim_multiscale = run!(graph_sim_multiscale,
+        meteo_day,
+        PlantMeteo.Constants(),
+        nothing;
+        check=true,
+        executor=SequentialEx()
+    )
+end