|
| 1 | +# Simple two-system interaction |
| 2 | + |
| 3 | +Here, a two-node network including a river and a dam is represented. The example is |
| 4 | +based on the Lower Campaspe catchment - a small semi-arid basin in North-Central Victoria, |
| 5 | +Australia. |
| 6 | + |
| 7 | + |
| 8 | + |
| 9 | +The dam is the primary water store for farmers in the area. In this simplified example, we |
| 10 | +investigate the effect of changing water demands and related policies on historic dam levels. |
| 11 | +Water demands may shift due to choice in cultivated crops, changing water management |
| 12 | +practices, investment into more efficient irrigation systems, choice in crops, or a |
| 13 | +systemic change to water policies. |
| 14 | + |
| 15 | +In essence, this is a cursory investigation into "what might have been" if the regional |
| 16 | +context were different. |
| 17 | + |
| 18 | +For the purpose of this example, the farm water requirements are defined as a volume of |
| 19 | +daily water requirements throughout a growing season. A growing season is the period of time |
| 20 | +over which a crop is cultivated, and is assumed to be between *May* and *February*. |
| 21 | +In practice, water requirements are provided by another model. |
| 22 | + |
| 23 | +```julia |
| 24 | +using CSV, DataFrames, YAML |
| 25 | +using Dates |
| 26 | +using Plots |
| 27 | +using Streamfall |
| 28 | + |
| 29 | +# Load climate data - in this case from a CSV file with data for all nodes. |
| 30 | +# Indicate which columns are precipitation and evaporation data based on partial identifiers |
| 31 | +climate = Climate("../test/data/campaspe/climate/climate.csv", "_rain", "_evap") |
| 32 | + |
| 33 | +# Historic extractions from the dam |
| 34 | +extraction_data = CSV.read("../test/data/campaspe/gauges/dam_extraction.csv", DataFrame; comment="#") |
| 35 | + |
| 36 | +# Load the example network |
| 37 | +sn = load_network("Example Network", "../test/data/campaspe/two_node_network.yml") |
| 38 | + |
| 39 | +# Run the model for the basin to obtain baseline values |
| 40 | +run_basin!(sn, climate; extraction=extraction_data) |
| 41 | +baseline_dam_level = sn[2].level |
| 42 | +baseline_dam_outflow = sn[2].outflow |
| 43 | + |
| 44 | +# Get represented dates for simulation |
| 45 | +sim_dates = Streamfall.timesteps(climate) |
| 46 | + |
| 47 | +# Create DataFrame to use as a template to store water extractions |
| 48 | +extractions = copy(extraction_data) |
| 49 | + |
| 50 | +""" |
| 51 | +Convenience function handling interactions with all "external" models. |
| 52 | +Note: We are using a burn-in period of a year. |
| 53 | +""" |
| 54 | +function run_scenario(sn, climate, extractions, increased_demand; burn_in=366) |
| 55 | + reset!(sn) |
| 56 | + inlets, outlets = find_inlets_and_outlets(sn) |
| 57 | + |
| 58 | + extractions = copy(extractions) |
| 59 | + |
| 60 | + sim_dates = Streamfall.timesteps(climate) |
| 61 | + prep_state!(sn, length(sim_dates)) |
| 62 | + |
| 63 | + for (ts, date) in enumerate(sim_dates) |
| 64 | + if ((month(date) >= 5) || (month(date) <= 2)) |
| 65 | + # Within growing season (May - Feb) |
| 66 | + |
| 67 | + # Additional agricultural water demand (ML/day) from historic conditions |
| 68 | + # This would normally come from another model indicating daily water demands. |
| 69 | + extractions[ts, "406000_releases_[ML]"] += increased_demand |
| 70 | + end |
| 71 | + |
| 72 | + for outlet in outlets |
| 73 | + run_node!(sn, outlet, climate, ts; extraction=extractions) |
| 74 | + end |
| 75 | + end |
| 76 | + |
| 77 | + # Return dam levels for assessment |
| 78 | + return sn[2].level[burn_in:end] |
| 79 | +end |
| 80 | + |
| 81 | +""" |
| 82 | +A hypothetical Critical Threshold index. |
| 83 | +
|
| 84 | +Mean of proportional distance to critical threshold (default 55% of dam capacity). |
| 85 | +Values of 1 indicate the dam is always full, negative values indicate the dam levels are |
| 86 | +below the critical threshold. |
| 87 | +
|
| 88 | +Greater values indicate greater water security but may have trade-offs with regard to |
| 89 | +environmental outcomes and farm profitability. |
| 90 | +""" |
| 91 | +function critical_threshold_index(levels; threshold=0.55) |
| 92 | + max_level = 195 |
| 93 | + min_level = 160 # mAHD |
| 94 | + |
| 95 | + # Essentially proportion of dam capacity |
| 96 | + relative_level = ((levels .- min_level) ./ (max_level - min_level)) |
| 97 | + |
| 98 | + # Distance |
| 99 | + indicator = (relative_level .- threshold) ./ (1.0 - threshold) |
| 100 | + |
| 101 | + return round(mean(indicator), digits=4) |
| 102 | +end |
| 103 | + |
| 104 | + |
| 105 | +ct_index = critical_threshold_index(baseline_dam_level[366:end]) |
| 106 | +f = plot(sim_dates[366:end], baseline_dam_level[366:end]; label="Historic modelled ($(ct_index))") |
| 107 | +results = Vector{Vector{Float64}}(undef, 5) |
| 108 | + |
| 109 | +for (i, daily_demand) in enumerate([-2.0, -0.5, 1.0, 3.0, 5.0]) |
| 110 | + results[i] = run_scenario(sn, climate, extractions, daily_demand) |
| 111 | + |
| 112 | + ct_index = critical_threshold_index(results[i]) |
| 113 | + |
| 114 | + plot!(sim_dates[366:end], results[i]; label = "+$(daily_demand) ML/day ($(ct_index))") |
| 115 | +end |
| 116 | +display(f) |
| 117 | +# savefig("simple_water_demand.png") |
| 118 | + |
| 119 | +Streamfall.temporal_cross_section(sim_dates, calib_data[:, "406000"], sn[2].level) |
| 120 | +``` |
| 121 | + |
| 122 | + |
0 commit comments