From 4becf7b1fd3dd5842ec865dbd5fbf3401039744a Mon Sep 17 00:00:00 2001 From: dalmijn <92092029+dalmijn@users.noreply.github.com> Date: Tue, 11 Jun 2024 15:42:21 +0200 Subject: [PATCH] Quick boundary fix (#430) * Quick boundary fix * Unpack `aquifer` and `constanthead` for readability * Update changelog --------- Co-authored-by: Willem van Verseveld --- docs/src/changelog.md | 3 +++ src/sbm_gwf_model.jl | 19 ++++++++++--------- 2 files changed, 13 insertions(+), 9 deletions(-) diff --git a/docs/src/changelog.md b/docs/src/changelog.md index c37e167a6..a0bd2349c 100644 --- a/docs/src/changelog.md +++ b/docs/src/changelog.md @@ -38,6 +38,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 `Polyester.jl`. - Fixed required states of the model type `sbm_gwf`: added `h_av` for the river and land domain. +- For the `constanthead` boundary of `GroundwaterFlow` the `head` should not be changed (was + set to `top` elevation of the `aquifer` if `head` > `top`), and `exfiltwater` should be 0 + for these boundary cells. ### Changed - Stop exposing scalar variables through BMI. The `BMI.get_value_ptr` function was not diff --git a/src/sbm_gwf_model.jl b/src/sbm_gwf_model.jl index 6a76b5bdf..d715399d2 100644 --- a/src/sbm_gwf_model.jl +++ b/src/sbm_gwf_model.jl @@ -483,6 +483,8 @@ function update(model::Model{N,L,V,R,W,T}) where {N,L,V,R,W,T<:SbmGwfModel} @unpack lateral, vertical, network, clock, config = model inds_riv = network.index_river + aquifer = lateral.subsurface.flow.aquifer + constanthead = lateral.subsurface.flow.constanthead # extract water levels h_av [m] from the land and river domains # this is used to limit open water evaporation @@ -512,7 +514,7 @@ function update(model::Model{N,L,V,R,W,T}) where {N,L,V,R,W,T<:SbmGwfModel} # determine stable time step for groundwater flow conductivity_profile = get(config.input.lateral.subsurface, "conductivity_profile", "uniform") - dt_gw = stable_timestep(lateral.subsurface.flow.aquifer, conductivity_profile) # time step in day (Float64) + dt_gw = stable_timestep(aquifer, conductivity_profile) # time step in day (Float64) dt_sbm = (vertical.dt / tosecond(basetimestep)) # vertical.dt is in seconds (Float64) if dt_gw < dt_sbm @warn( @@ -529,18 +531,17 @@ function update(model::Model{N,L,V,R,W,T}) where {N,L,V,R,W,T<:SbmGwfModel} # determine excess water depth [m] (exfiltwater) in groundwater domain (head > surface) # and reset head - exfiltwater = - ( - lateral.subsurface.flow.aquifer.head .- - min.(lateral.subsurface.flow.aquifer.head, lateral.subsurface.flow.aquifer.top) - ) .* storativity(lateral.subsurface.flow.aquifer) - lateral.subsurface.flow.aquifer.head .= - min.(lateral.subsurface.flow.aquifer.head, lateral.subsurface.flow.aquifer.top) + exfiltwater = (aquifer.head .- min.(aquifer.head, aquifer.top)) .* storativity(aquifer) + aquifer.head .= min.(aquifer.head, aquifer.top) + + # Adjust for constant head boundary of groundwater domain + exfiltwater[constanthead.index] .= 0 + aquifer.head[constanthead.index] .= constanthead.head # update vertical sbm concept (runoff, ustorelayerdepth and satwaterdepth) update_after_subsurfaceflow( vertical, - (network.land.altitude .- lateral.subsurface.flow.aquifer.head) .* 1000.0, # zi [mm] in vertical concept SBM + (network.land.altitude .- aquifer.head) .* 1000.0, # zi [mm] in vertical concept SBM exfiltwater .* 1000.0, )