diff --git a/Project.toml b/Project.toml
index d52a0a2e8..9c5b45560 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "Wflow"
 uuid = "d48b7d99-76e7-47ae-b1d5-ff0c1cf9a818"
 authors = ["Deltares and contributors"]
-version = "0.7.2"
+version = "0.7.3"
 
 [deps]
 BasicModelInterface = "59605e27-edc0-445a-b93d-c09a3a50b330"
diff --git a/build/wflow_cli/Setup/Installer/Application.Setup.wixproj b/build/wflow_cli/Setup/Installer/Application.Setup.wixproj
index b74452efd..b8f77e563 100644
--- a/build/wflow_cli/Setup/Installer/Application.Setup.wixproj
+++ b/build/wflow_cli/Setup/Installer/Application.Setup.wixproj
@@ -6,7 +6,7 @@
     <ProductVersion>3.10</ProductVersion>
     <ProjectGuid>{39da8083-e405-42e6-850c-d25685e91f81}</ProjectGuid>
     <SchemaVersion>2.0</SchemaVersion>
-    <ReleaseVersion>0.7.2</ReleaseVersion>
+    <ReleaseVersion>0.7.3</ReleaseVersion>
     <OutputType>Package</OutputType>
     <DefineSolutionProperties>false</DefineSolutionProperties>
   </PropertyGroup>
diff --git a/docs/src/changelog.md b/docs/src/changelog.md
index 5b4ccfaeb..4a2391ed1 100644
--- a/docs/src/changelog.md
+++ b/docs/src/changelog.md
@@ -5,7 +5,7 @@ All notable changes to this project will be documented in this file.
 The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
-## [unreleased]
+## v0.7.3 - 2024-01-12
 
 ### Fixed
 - Documentation: add leakage term to the wflow\_sbm figure, document external input
@@ -15,7 +15,14 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
   (BasicModelInterface.jl), including function arguments, return types and the BMI
   specification that arrays are always flattened (this was not the case for variables stored
   as 2-dimensional arrays or as vector of SVectors).
-- Bump compat for NCDatasets to 0.13.
+- Bump compat for NCDatasets to 0.13, 0.14.
+- The solution for lake outflow as part of the Modified Puls Approach. The inflow and
+  outflow variables are defined for period `Δt`, and not at `t1` and `t2` (instantaneous) as
+  in the original mass balance equation of the Modified Puls Approach. Because of this, the
+  terms of the quadratic equation (and solution) were fixed.
+- Use `kvfrac` for the computation of vertical saturated hydraulic conductivity at the
+  bottom of the soil layer, since `kvfrac` is also used for the computation of vertical
+  unsaturated flow. 
 
 ### Changed
 - For cyclic parameters different cyclic time inputs are supported (only one common cyclic
diff --git a/docs/src/model_docs/lateral/waterbodies.md b/docs/src/model_docs/lateral/waterbodies.md
index 87331d458..335ab2867 100644
--- a/docs/src/model_docs/lateral/waterbodies.md
+++ b/docs/src/model_docs/lateral/waterbodies.md
@@ -91,8 +91,8 @@ Inserting this equation in the mass balance gives:
 The solution for Q is then:
 
 ```math
-    Q = { \left( -LF + \sqrt{LF^{2} + 2  \left( SI - \dfrac{A*H_{0}}{\Delta t} \right)
-    } \right) }^{2} \text{for } SI > \dfrac{A H_{0}}{\Delta t}  \text{ and where}\\
+    Q = { \left( \dfrac{-LF + \sqrt{LF^{2} + 4  \left( SI - \dfrac{A*H_{0}}{\Delta t} \right)}}
+    {2} \right) }^{2} \text{for } SI > \dfrac{A H_{0}}{\Delta t}  \text{ and where}\\
     LF = \dfrac{A}{\Delta t \sqrt{\alpha}} \\~\\
     Q = 0 \text{ for } SI \leq \dfrac{A*H_{0}}{\Delta t}
 ```
diff --git a/docs/src/model_docs/params_lateral.md b/docs/src/model_docs/params_lateral.md
index 56cca1277..50326452e 100644
--- a/docs/src/model_docs/params_lateral.md
+++ b/docs/src/model_docs/params_lateral.md
@@ -191,7 +191,7 @@ scales (hillslope) in reality, not represented by the model resolution.
 | `dw` | drain width | m | - |
 | `zi` | pseudo-water table depth (top of the saturated zone) | m | - |
 | `exfiltwater` | exfiltration (groundwater above surface level, saturated excess conditions) | m Δt⁻¹ | - |
-| `recharge` | net recharge to saturated store  | m Δt⁻¹ | - |
+| `recharge` | net recharge to saturated store  | m``^2`` Δt⁻¹ | - |
 | `ssf` | subsurface flow  | m``^3`` d``{-1}``  | - |
 | `ssfin` | inflow from upstream cells | m``^3`` d``{-1}``  | - |
 | `ssfmax` | maximum subsurface flow | m``^2`` d``{-1}``  | - |
diff --git a/src/flow.jl b/src/flow.jl
index 4aaa05e1c..af77daa89 100644
--- a/src/flow.jl
+++ b/src/flow.jl
@@ -398,7 +398,7 @@ end
     dw::Vector{T} | "m"                    # Flow width [m]
     zi::Vector{T} | "m"                    # Pseudo-water table depth [m] (top of the saturated zone)
     exfiltwater::Vector{T} | "m Δt-1"      # Exfiltration [m Δt⁻¹] (groundwater above surface level, saturated excess conditions)
-    recharge::Vector{T} | "m Δt-1"         # Net recharge to saturated store [m Δt⁻¹]
+    recharge::Vector{T} | "m2 Δt-1"        # Net recharge to saturated store [m² Δt⁻¹]
     ssf::Vector{T} | "m3 d-1"              # Subsurface flow [m³ d⁻¹]
     ssfin::Vector{T} | "m3 d-1"            # Inflow from upstream cells [m³ d⁻¹]
     ssfmax::Vector{T} | "m2 d-1"           # Maximum subsurface flow [m² d⁻¹]
diff --git a/src/reservoir_lake.jl b/src/reservoir_lake.jl
index 410745eef..5abd0506a 100644
--- a/src/reservoir_lake.jl
+++ b/src/reservoir_lake.jl
@@ -513,7 +513,8 @@ function update(lake::Lake, i, inflow, doy, timestepsecs)
 
         if si_factor_adj > 0.0
             outflow = pow(
-                -lakefactor + pow((pow(lakefactor, 2.0) + 2.0 * si_factor_adj), 0.5),
+                0.5 *
+                (-lakefactor + pow((pow(lakefactor, 2.0) + 4.0 * si_factor_adj), 0.5)),
                 2.0,
             )
         else
diff --git a/src/sbm.jl b/src/sbm.jl
index e9ae822e3..954b0ab51 100644
--- a/src/sbm.jl
+++ b/src/sbm.jl
@@ -928,7 +928,7 @@ function update_until_recharge(sbm::SBM, config)
                 actcapflux = actcapflux + toadd
             end
         end
-        deepksat = sbm.kv₀[i] * exp(-sbm.f[i] * sbm.soilthickness[i])
+        deepksat = sbm.kvfrac[i][end] * sbm.kv₀[i] * exp(-sbm.f[i] * sbm.soilthickness[i])
         deeptransfer = min(satwaterdepth, deepksat)
         actleakage = max(0.0, min(sbm.maxleakage[i], deeptransfer))
 
diff --git a/test/reservoir_lake.jl b/test/reservoir_lake.jl
index c05a1f7fc..66d251fe6 100644
--- a/test/reservoir_lake.jl
+++ b/test/reservoir_lake.jl
@@ -50,10 +50,10 @@ end
     )
 
     Wflow.update(lake, 1, 2500.0, 181, 86400.0)
-    @test lake.outflow[1] ≈ 85.31903276150577
-    @test lake.totaloutflow[1] ≈ 7.371564430594098e6
-    @test lake.storage[1] ≈ 3.551103576940606e9
-    @test lake.waterlevel[1] ≈ 19.672569557695734
+    @test lake.outflow[1] ≈ 85.14292808113598
+    @test lake.totaloutflow[1] ≈ 7.356348986210149e6
+    @test lake.storage[1] ≈ 3.55111879238499e9
+    @test lake.waterlevel[1] ≈ 19.672653848925634
     @test lake.precipitation[1] ≈ 20.0
     @test lake.evaporation[1] ≈ 3.2
 end