- fixed bug in Ks_tmp calculation for soilDB == 0

src/MPR/mo_mpr_runoff.f90

@@ -135,7 +135,7 @@ contains
     !-----------------------------
     ! FAST INTERFLOW
     !-----------------------------
-    ! HL1 = f(soil properties; No refrence found)
+    ! HL1 = f(soil properties; No reference found)
     ! Based on the saturation deficit from the field capacity status
     ! seems more reasonable and intutative.
     ! NOTE: This value for the sandy soils will have higher value of HL1, as compared to 
@@ -144,9 +144,9 @@ contains
     L1_HL1 = upscale_arithmetic_mean(nL0_in_L1, Upp_row_L1, Low_row_L1, &
             Lef_col_L1, Rig_col_L1, cell_id0, mask0, nodata_dp, tmp)
 
-    ! 1/K0 = f(terrian slope) [Booij, et. al.(2005), JoH]
-    ! Steeper slopes resists (1/K0) fast water flows lesser as 
-    ! compared to that on the flater slope areas. 
+    ! 1/K0 = f(terrain slope) [Booij, et. al.(2005), JoH]
+    ! Steeper slopes resists (1/K0) fast water flows lesser as
+    ! compared to that on the flater slope areas.
     ! Assuming that above relationship holds for all kind of land cover classes
 
     ! In the forested area surface resistance to fast interflow is higher as compared
@@ -173,16 +173,16 @@ contains
     L1_K1 = upscale_arithmetic_mean(nL0_in_L1, Upp_row_L1, Low_row_L1, &
             Lef_col_L1, Rig_col_L1, cell_id0, mask0, nodata_dp, tmp)
 
-    ! minimum value of K1 is 1-day 
+    ! minimum value of K1 is 1-day
     L1_K1 = merge(2.0_dp, L1_K1, L1_K1 .lt. 2.0_dp)
 
 
     ! alpha = f(soil type; variabitity of Ks)
-    ! Lower the alpha (exponent of slow interflow) means lower amount of 
-    ! water released from the storage to contribute for slow interflow. 
+    ! Lower the alpha (exponent of slow interflow) means lower amount of
+    ! water released from the storage to contribute for slow interflow.
     ! For instance sandy soils will have lower value of alpha as comapred to
     ! the clayey soils.
-    ! This assumption is quite realistic in physical sense... 
+    ! This assumption is quite realistic in physical sense...
     tmp = merge(param(5) * (1.0_dp / KsVar_H0) * (1.0_dp / (1.0_dp + SMs_FC0)), &
             nodata_dp, cell_id0 .ne. nodata_i4)
     L1_alpha = upscale_arithmetic_mean(nL0_in_L1, Upp_row_L1, Low_row_L1, Lef_col_L1, &

src/MPR/mo_mpr_soilmoist.f90

@@ -247,7 +247,6 @@ contains
     Ks = 0.0_dp
     Db = 0.0_dp
     if(allocated(Ks_non_till)) Ks_non_till = 0.0_dp
-    ! TODO:  is it necessary to loop over 1:maxval, isn't it required to only loop over existing classes?!
     max_LCover = maxval(LCOVER0)
     ! select case according to a given soil database flag
     SELECT CASE(iFlag_soilDB)
@@ -284,13 +283,12 @@ contains
               ! Effect of organic matter content
               ! This is taken into account in a simplified form by using
               ! the ratio of(Bd / BdOM)
-              Ks_tmp = Ks_tmp * (DbM(i, j) / Db(i, j, L))
-              Ks(i, j, L) = Ks_tmp
+              Ks(i, j, L) = Ks(i, j, L) * (DbM(i, j) / Db(i, j, L))
               ! estimated SMs_till & van Genuchten's shape parameter (n)
               call Genuchten(thetaS_till(i, j, L), Genu_Mual_n, Genu_Mual_alpha, &
                       param(4 : 9), sand(i, j), clay(i, j), Db(i, j, L))
               ! estimating field capacity
-              call field_cap(thetaFC_till(i, j, L), Ks_tmp, thetaS_till(i, j, L), Genu_Mual_n)
+              call field_cap(thetaFC_till(i, j, L), Ks(i, j, L), thetaS_till(i, j, L), Genu_Mual_n)
               ! estimating permanent wilting point
               call PWP(Genu_Mual_n, Genu_Mual_alpha, thetaS_till(i, j, L), thetaPW_till(i, j, L))
             end do
@@ -301,7 +299,7 @@ contains
                     param(4 : 9), sand(i, j), clay(i, j), DbM(i, j))
             ! estimate field capacity
             call field_cap(thetaFC(i, j - tmp_minSoilHorizon), &
-                    Ks_tmp, thetaS(i, j - tmp_minSoilHorizon), Genu_Mual_n)
+                    Ks(i, j, 1), thetaS(i, j - tmp_minSoilHorizon), Genu_Mual_n)
             ! estimate permanent wilting point
             call PWP(Genu_Mual_n, Genu_Mual_alpha, thetaS(i, j - tmp_minSoilHorizon), &
                     thetaPW(i, j - tmp_minSoilHorizon))