Fix required states sbm_gwf model type (#424)

* Fix required states `sbm_gwf` model type
Added `h_av` of river and land domain for vertical `SBM` concept.

* Update changelog.md

* Fix ldd map of test model sbm-groundwater

* Extend state checking test

build/create_binaries/download_test_data.jl

@@ -45,9 +45,9 @@ forcing_sbm_gw_path = testdata(
 )
 forcing_meuse_path = testdata(v"0.2.8", "forcing_meuse.nc", "forcing_meuse.nc")
 staticmaps_sbm_gw_path =
-    testdata(v"0.2.2", "staticmaps-sbm-groundwater.nc", "staticmaps-sbm-groundwater.nc")
+    testdata(v"0.2.3", "staticmaps-sbm-groundwater.nc", "staticmaps-sbm-groundwater.nc")
 instates_sbm_gw_path =
-    testdata(v"0.2.2", "instates-example-sbm-gwf.nc", "instates-example-sbm-gwf.nc")
+    testdata(v"0.2.3", "instates-example-sbm-gwf.nc", "instates-example-sbm-gwf.nc")
 lake_sh_1_path = testdata(v"0.2.1", "lake_sh_1.csv", "lake_sh_1.csv")
 lake_sh_2_path = testdata(v"0.2.1", "lake_sh_2.csv", "lake_sh_2.csv")
 lake_hq_2_path = testdata(v"0.2.1", "lake_hq_2.csv", "lake_hq_2.csv")

docs/src/changelog.md

@@ -35,6 +35,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
   zero when `LoopVectorization` is applied (`@tturbo`) to the for loop of these functions.
   This issue occured on a virtual machine, Windows 10 Enterprise, with Intel(R) Xeon(R) Gold
   6144 CPU (2 processors). This has been fixed by replacing `@tturbo` with `reduction` of
+  `Polyester.jl`.
+- Fixed required states of the model type `sbm_gwf`: added `h_av` for the river and land
+  domain.
 
 ### Changed
 - Stop exposing scalar variables through BMI. The `BMI.get_value_ptr` function was not

src/states.jl

@@ -122,7 +122,7 @@ function extract_required_states(config::Config)
         if land_routing == "local-inertial"
             land_states = (:qx, :qy, :h, :h_av)
         elseif land_routing == "kinematic-wave"
-            if model_type == "sbm"
+            if model_type == "sbm" || model_type == "sbm_gwf"
                 land_states = (:q, :h, :h_av)
             else
                 land_states = (:q, :h)
@@ -152,7 +152,7 @@ function extract_required_states(config::Config)
             :outlagg,
             :outgrav,
         )
-    elseif model_type == "sbm"
+    elseif model_type == "sbm" || model_type == "sbm_gwf"
         river_states = (:q, :h, :h_av)
     else
         river_states = (:q, :h)

test/io.jl

@@ -445,7 +445,6 @@ end
 @testset "State checking" begin
     tomlpath = joinpath(@__DIR__, "sbm_config.toml")
     config = Wflow.Config(tomlpath)
-    config = Wflow.Config(tomlpath)
 
     # Extracting required states and test if some are covered (not all are tested!)
     required_states = Wflow.extract_required_states(config)
@@ -454,6 +453,8 @@ end
     @test (:vertical, :canopystorage) in required_states
     @test (:lateral, :subsurface, :ssf) in required_states
     @test (:lateral, :river, :q) in required_states
+    @test (:lateral, :river, :h_av) in required_states
+    @test (:lateral, :land, :h_av) in required_states
     @test !((:lateral, :river, :lake, :waterlevel) in required_states)
 
     # Adding an unused state the see if the right warning message is thrown
@@ -470,4 +471,16 @@ end
     delete!(config.state["vertical"], "additional_state")
     delete!(config.state["vertical"], "snow")
     @test_throws ArgumentError Wflow.check_states(config)
+
+    # Extracting required states for model type sbm_gwf and test if some are covered
+    tomlpath = joinpath(@__DIR__, "sbm_gwf_config.toml")
+    config = Wflow.Config(tomlpath)
+    required_states = Wflow.extract_required_states(config)
+    @test (:vertical, :satwaterdepth) in required_states
+    @test (:vertical, :ustorelayerdepth) in required_states
+    @test (:vertical, :canopystorage) in required_states
+    @test (:lateral, :subsurface, :flow, :aquifer, :head) in required_states
+    @test (:lateral, :river, :q) in required_states
+    @test (:lateral, :river, :h_av) in required_states
+    @test (:lateral, :land, :h_av) in required_states
 end

test/run_sbm_gwf.jl

@@ -152,37 +152,37 @@ model = Wflow.run_timestep(model)
 @testset "second timestep warm start" begin
     sbm = model.vertical
     @test sbm.runoff[1] == 0.0
-    @test sbm.soilevap[1] == 0.2927279656884887
+    @test sbm.soilevap[1] ≈ 0.2926821507238362f0
     @test sbm.transpiration[1] ≈ 1.0122634204681036f0
 end
 
 @testset "overland flow warm start (kinematic wave)" begin
     q = model.lateral.land.q_av
-    @test sum(q) ≈ 1.4411427003142072f-5
+    @test sum(q) ≈ 1.4589771292158736f-5
 end
 
 @testset "river domain warm start (kinematic wave)" begin
     q = model.lateral.river.q_av
     river = model.lateral.river
-    @test sum(q) ≈ 0.011317441400219936f0
-    @test q[6] ≈ 0.002255753266287542f0
-    @test river.volume[6] ≈ 2.1212499727096956f0
-    @test river.inwater[6] ≈ -6.888767545965421f-5
-    @test q[13] ≈ 8.664553314598283f-5
-    @test q[network.river.order[end]] ≈ 0.002258255913217909f0
+    @test sum(q) ≈ 0.011909898793648826f0
+    @test q[6] ≈ 0.0024332762253494573f0
+    @test river.volume[6] ≈ 2.226607059761764f0
+    @test river.inwater[6] ≈ -1.3378627340716941f-5
+    @test q[13] ≈ 7.318205070445492f-5
+    @test q[network.river.order[end]] ≈ 0.002470511395424463f0
 end
 
 @testset "groundwater warm start" begin
     gw = model.lateral.subsurface
     @test gw.river.stage[1] ≈ 1.2031171676781156f0
     @test gw.flow.aquifer.head[17:21] ≈ [
-        1.2195537323067487f0,
-        1.2741313108104333f0,
+        1.227685769857181f0,
+        1.286975535107666f0,
         1.7999999523162842f0,
-        1.5862473868459186f0,
-        1.202268433263572f0,
+        1.5892917098077068f0,
+        1.2093841308150408f0,
     ]
-    @test gw.river.flux[1] ≈ -6.044849112655228f0
+    @test gw.river.flux[1] ≈ -6.680320482009705f0
     @test gw.drain.flux[1] ≈ 0.0
     @test gw.recharge.rate[19] ≈ -0.0014241196552847502f0
 end

test/runtests.jl

@@ -62,9 +62,9 @@ forcing_sbm_gw_path = testdata(
 )
 forcing_meuse_path = testdata(v"0.2.8", "forcing_meuse.nc", "forcing_meuse.nc")
 staticmaps_sbm_gw_path =
-    testdata(v"0.2.2", "staticmaps-sbm-groundwater.nc", "staticmaps-sbm-groundwater.nc")
+    testdata(v"0.2.3", "staticmaps-sbm-groundwater.nc", "staticmaps-sbm-groundwater.nc")
 instates_sbm_gw_path =
-    testdata(v"0.2.2", "instates-example-sbm-gwf.nc", "instates-example-sbm-gwf.nc")
+    testdata(v"0.2.3", "instates-example-sbm-gwf.nc", "instates-example-sbm-gwf.nc")
 lake_sh_1_path = testdata(v"0.2.1", "lake_sh_1.csv", "lake_sh_1.csv")
 lake_sh_2_path = testdata(v"0.2.1", "lake_sh_2.csv", "lake_sh_2.csv")
 lake_hq_2_path = testdata(v"0.2.1", "lake_hq_2.csv", "lake_hq_2.csv")

test/sbm_gwf_config.toml

@@ -25,10 +25,12 @@ ustorelayerdepth = "ustorelayerdepth"
 
 [state.lateral.river]
 h = "h_river"
+h_av = "h_av_river"
 q = "q_river"
 
 [state.lateral.land]
 h = "h_land"
+h_av = "h_av_land"
 q = "q_land"
 
 [state.lateral.subsurface.flow.aquifer]