diff --git a/dev_notes/todo.txt b/dev_notes/todo.txt index f6cb73e..e7de47d 100644 --- a/dev_notes/todo.txt +++ b/dev_notes/todo.txt @@ -2,6 +2,12 @@ *** High pri *** + Index set dependency propagation along connections is wrong in some edge cases. + + e.g if it transports a quantity that is indexed over an index set, this index set dependency is not propagated by the connection aggregation. + + + Make functionality to patch in data from a sub-model in a data set. E.g. you have a setup for Simply+NIVAFjord and you want to patch in from a different file that just has the Simply part of the setup. Should not be an "include" on file level, because that causes too much ambiguity about how this should be saved if overwritten. @@ -168,15 +174,6 @@ *** Intermediate-pri *** - Got an (internal) error in the following case (nivafjord pipe boundary model) - flux(bnd.water, layer.water[vert.specific]) ... - - bnd.water.tox can be distributed over tox_index, but that dependency was ruled out due to lack of differing parametrization - layer.water.tox was distributed over tox_index with the actual dependency. - Caused problem when the dissolved flux were to be added to the connection aggregate (lack of index set). - Probably also a problem for vert.top etc. - It is a bit tricky for the framework to predetermine this.. But there is maybe something like this for graph aggregates already that could be applied here too. - If an input series is provided indexed over "subcatchment", but should be indexed over "water body" and "subcatchment" is a union member of "water body", that should be made to work. (and all similar cases) Just re-map the index. Store name of expected model file in data sets (?) diff --git a/models/demonstration_models/nivafjord_lake_simplycnp_glacier.txt b/models/demonstration_models/nivafjord_lake_simplycnp_glacier.txt index 1062060..48b18e4 100644 --- a/models/demonstration_models/nivafjord_lake_simplycnp_glacier.txt +++ b/models/demonstration_models/nivafjord_lake_simplycnp_glacier.txt @@ -27,4 +27,24 @@ model("NIVAFjord-SimplyCNP-Glacier") { # Add a module for the lake to discharge out of the system. load("modules/nivafjord/lake_discharge.txt", module("NIVAFjord lake discharge", basin, layer, water, dz, vert, dims)) + + ti : index_set("Tracer index") + tr : quantity("Tracer", ti) + + module("Extra things", version(0, 0, 0)) { + + par_group("Tracer", tr) { + init_tr : par_real("Initial tracer", [m g, l-1], 1) + } + + var(soil.water.tr, [k g, k m-2], [m g, l-1], "Soil water tracer") @initial_conc { init_tr } + + var(gw.water.tr, [k g, k m-2], [m g, l-1], "Groundwater tracer") + + var(river.water.tr, [k g], [m g, l-1], "River water tracer") + + var(layer.water.tr, [k g], [m g, l-1], "Layer water tracer") + + + } } \ No newline at end of file diff --git a/models/example_data/NIVAFjord/glacier_lake_example.dat b/models/example_data/NIVAFjord/glacier_lake_example.dat index ca2b62d..4767f19 100644 --- a/models/example_data/NIVAFjord/glacier_lake_example.dat +++ b/models/example_data/NIVAFjord/glacier_lake_example.dat @@ -19,6 +19,8 @@ data_set { phyt_type : index_set("Phytoplankton type") [ "Phyto" ] + ti : index_set("Tracer index") [ "Cl" "Ca" ] + position_map(layer) [! 0 : 1 ] diff --git a/src/model_compilation.cpp b/src/model_compilation.cpp index 0c31cc8..4414240 100644 --- a/src/model_compilation.cpp +++ b/src/model_compilation.cpp @@ -18,7 +18,10 @@ std::string Model_Instruction::debug_string(Model_Application *app) const { std::stringstream ss; - if(type == Model_Instruction::Type::compute_state_var) + if(type == Model_Instruction::Type::invalid) + // NOTE: This sometimes happens when it prints instruction dependencies. + ss << "(invalid)"; + else if(type == Model_Instruction::Type::compute_state_var) ss << "\"" << app->vars[var_id]->name << "\""; else if(type == Model_Instruction::Type::subtract_discrete_flux_from_source) ss << "\"" << app->vars[source_id]->name << "\" -= \"" << app->vars[var_id]->name << "\""; @@ -336,6 +339,25 @@ insert_var_order_depencencies(Model_Application *app, Model_Instruction *instr, } } +Index_Set_Tuple +get_quantity_index_sets(Model_Application *app, Var_Id var_id) { + Index_Set_Tuple result; + + auto var = app->vars[var_id]; + auto &loc = var->loc1; + if(!is_located(loc)) + fatal_error(Mobius_Error::internal, "Misuse of get_quantity_index_sets"); + // Iterate from 1 to not include the compartment. + for(int idx = 1; idx < loc.n_components; ++idx) { + auto quantity = app->model->components[loc.components[idx]]; + if(quantity->decl_type != Decl_Type::quantity) + fatal_error(Mobius_Error::internal, "Misuse of get_quantity_index_sets"); + for(auto index_set : quantity->index_sets) + result.insert(index_set); + } + return result; +} + void resolve_basic_dependencies(Model_Application *app, std::vector &instructions) { @@ -351,6 +373,13 @@ resolve_basic_dependencies(Model_Application *app, std::vectorvars[instr.var_id]; if(var->specific_target.get()) register_dependencies(var->specific_target.get(), &code_depends); + if(var->is_mass_balance_quantity()) { + // By default always distribute over quantity distributions. + // Otherwise there are several failure edge cases where these are not propagated (through connections mostly). + // And there is little to gain by optimizing them away + auto quantity_index_sets = get_quantity_index_sets(app, instr.var_id); + insert_dependecies(app, &instr, quantity_index_sets); + } } for(auto &dep : code_depends) { @@ -820,22 +849,21 @@ process_graph_connection_aggregation(Model_Application *app, std::vectoredge_index_set); } - if(!agg_var->is_out) { // TODO: should (something like) this also be done for the source aggregate in directed_graph? - - // TODO: Make a better explanation of what is going on in this block and why it is needed (What is the failure case otherwise). + if(!agg_var->is_out) { // If the target compartment (not just what the connection indexes over) has an index set shared with the source compartment, we must index the target variable over that. + // Since the target could get a different value from the connection depending on its own index, we have to force it to be computed per each of these indexes even if it were not to have an index set dependency on this otherwise. + auto target_comp_id = app->vars[agg_var->agg_for]->loc1.first(); auto find_target = app->find_connection_component(agg_var->connection, target_comp_id); auto target_comp = model->components[target_comp_id]; - auto target_index_sets = find_target->index_sets; // vector copy; + auto target_index_sets = find_target->index_sets; // intentional copy; for(auto index_set : find_source->index_sets) { if(std::find(target_comp->index_sets.begin(), target_comp->index_sets.end(), index_set) != target_comp->index_sets.end()) target_index_sets.push_back(index_set); } - // Since the target could get a different value from the connection depending on its own index, we have to force it to be computed per each of these indexes even if it were not to have an index set dependency on this otherwise. for(auto index_set : target_index_sets) insert_dependency(app, &instructions[agg_var->agg_for.id], index_set); } diff --git a/test/models/another_simple.txt b/test/models/another_simple.txt index ed58111..0b5d2b0 100644 --- a/test/models/another_simple.txt +++ b/test/models/another_simple.txt @@ -2,21 +2,34 @@ model("Another simple") { + i : index_set("Q index") + a : compartment("A") - q : quantity("Q") + b : compartment("B") + r : quantity("R") + q : quantity("Q", i) + + c : connection("Connection") @directed_graph { a b } module("So simple", version(2, 0, 0)) { - var(a.q, [], "A Queue!") @override { - time.step=>> - } @no_store + par_group("Q pars", q) { + init_q : par_real("Init Q", [], 1) + } + + var(a.r, [k g], "A water") @initial { 10[k g] } + + var(a.r.q, [k g], "A queue") @initial_conc { init_q } + + var(b.r, [k g], "B water") + + var(b.r.q, [k g], "B queue") - r : property("P") - var(a.r, [], "Arrr!") { a.q*2 } + flux(a.r, c, [k g, day-1], "A r flux") { a.r*0.01[day-1] } } sol : solver("solver", inca_dascru, [2, hr], 0.01) - solve(sol, a.q) + solve(sol, a.r, b.r) } \ No newline at end of file diff --git a/test/models/another_simple_data.dat b/test/models/another_simple_data.dat index 6009d18..5edaa59 100644 --- a/test/models/another_simple_data.dat +++ b/test/models/another_simple_data.dat @@ -4,6 +4,18 @@ data_set { The doc string """ + i : index_set("Q index") [ "Q1" "Q2" ] + + connection("Connection") { + a : compartment("A") + b : compartment("B") + + directed_graph [ + a[] -> b[] + ] + } + + par_group("System") { par_datetime("Start date")