Merge pull request #152 from lanl-ansi/make_si_units

ADD: Functionality for transforming network and solution data between SI and per-unit systems

CHANGELOG.md

@@ -4,6 +4,9 @@ WaterModels.jl Change Log
 ### Staged
 - nothing
 
+### v0.9.3
+- Add a `make_si_units!` function for transforming network data.
+
 ### v0.9.2
 - Add support for Interpolations v0.14.
 

Project.toml

@@ -2,7 +2,7 @@ name = "WaterModels"
 uuid = "7c60b362-08f4-5b14-8680-cd67a3e18348"
 authors = ["Byron Tasseff"]
 repo = "https://github.com/lanl-ansi/WaterModels.jl"
-version = "0.9.2"
+version = "0.9.3"
 
 [deps]
 HiGHS = "87dc4568-4c63-4d18-b0c0-bb2238e4078b"

docs/src/examples.md

@@ -31,8 +31,12 @@ set_flow_partitions_si!(data_mn, 1.0, 1.0e-4)
 highs = JuMP.optimizer_with_attributes(HiGHS.Optimizer, "time_limit" => 60.0)
 result = solve_mn_owf(data_mn, LRDWaterModel, highs)
 ```
-Note that results are presented in an automatically-applied per-unit system.
-The instance is also challenging, and only a feasible solution is returned within the time limit for the script above.
+The instance is challenging, and only a feasible solution is returned within the time limit for the script above.
+Also note that results are presented in an automatically-applied per-unit system.
+To convert the solution to SI units, the following can be executed:
+```julia
+make_si_units!(result["solution"])
+```
 
 ## References
 [1] Alperovits, E., & Shamir, U. (1977). Design of optimal water distribution systems. _Water Resources Research_, _13_(6), 885-900.

docs/src/index.md

@@ -57,12 +57,15 @@ After solving the problem, results can then be analyzed, e.g.,
 # The termination status of the optimization solver.
 result["termination_status"]
 
+# Transform solution data to SI units.
+make_si_units!(result["solution"])
+
 # The flow along pipe 4 in cubic meters per second.
-result["solution"]["pipe"]["4"]["q"] * result["solution"]["base_flow"]
+result["solution"]["pipe"]["4"]["q"]
 
 # The total hydraulic head at node 2 in meters.
-result["solution"]["node"]["2"]["h"] * result["solution"]["base_head"]
+result["solution"]["node"]["2"]["h"]
 
 # The pressure head at node 2 in meters.
-result["solution"]["node"]["2"]["p"] * result["solution"]["base_head"]
+result["solution"]["node"]["2"]["p"]
 ```

docs/src/quickguide.md

@@ -94,24 +94,24 @@ result = solve_des(data, LRDWaterModel, HiGHS.Optimizer)
 ```
 
 For example, the algorithm's runtime and final objective value can be accessed with
-```
+```julia
 result["solve_time"] # Total solve time required (seconds).
 result["objective"] # Final objective value (in units of the objective).
 ```
 
 The `"solution"` field contains detailed information about the solution produced by the `solve` method.
 For example, the following dictionary comprehension can be used to inspect the flows in the solution:
-```
+```julia
 flows = Dict(name => data["q"] for (name, data) in result["solution"]["des_pipe"])
 ```
 
 To determine the design pipes that were selected via the optimization, the following can be used:
-```
+```julia
 pipes_selected = filter(x -> x.second["status"] == 1, result["solution"]["des_pipe"])
 ```
 
 To retrieve the subset of the original pipe dataset, the following can be used:
-```
+```julia
 pipes_subset = filter(x -> x.first in keys(pipes_selected), data["des_pipe"])
 ```
 
@@ -151,4 +151,11 @@ wm = instantiate_model(data, LRDWaterModel, WaterModels.build_des);
 println(wm.model)
 
 result = optimize_model!(wm, optimizer = HiGHS.Optimizer)
+```
+
+## Solution Unit Conversion
+The default behavior of WaterModels produces solution results in non-dimensionalized units.
+To recover solutions in SI units, the following function can be used:
+```julia
+make_si_units!(result["solution"])
 ```

docs/src/result-data.md

@@ -219,4 +219,8 @@ Note also that per-unit base quantities are also reported in the solution data d
 ```julia
 # Should return `true`.
 result["solution"]["base_flow"] == data["base_flow"]
+```
+For convenience, solution data can be transformed to SI units using
+```julia
+make_si_units!(result["solution"])
 ```

src/core/constraint_template.jl

@@ -396,9 +396,9 @@ function constraint_pipe_head_loss(
 
     wm_data = get_wm_data(wm.data)
     head_loss, viscosity = wm_data["head_loss"], wm_data["viscosity"]
-    base_length = get(wm_data, "base_length", 1.0)
-    base_mass = get(wm_data, "base_mass", 1.0)
-    base_time = get(wm_data, "base_time", 1.0)
+    base_length = wm_data["per_unit"] ? wm_data["base_length"] : 1.0
+    base_mass = wm_data["per_unit"] ? wm_data["base_mass"] : 1.0
+    base_time = wm_data["per_unit"] ? wm_data["base_time"] : 1.0
 
     r = _calc_pipe_resistance(
         ref(wm, nw, :pipe, a),
@@ -650,9 +650,9 @@ function constraint_on_off_des_pipe_head_loss(
     exponent = _get_exponent_from_head_loss_form(head_loss)
 
     wm_data = get_wm_data(wm.data)
-    base_length = get(wm_data, "base_length", 1.0)
-    base_mass = get(wm_data, "base_mass", 1.0)
-    base_time = get(wm_data, "base_time", 1.0)
+    base_length = wm_data["per_unit"] ? wm_data["base_length"] : 1.0
+    base_mass = wm_data["per_unit"] ? wm_data["base_mass"] : 1.0
+    base_time = wm_data["per_unit"] ? wm_data["base_time"] : 1.0
 
     r = _calc_pipe_resistance(des_pipe, head_loss, viscosity, base_length, base_mass, base_time)
     q_max_reverse = min(get(des_pipe, "flow_max_reverse", 0.0), 0.0)