Discretization of the bottom logarithmic layer drag is corrected

source/model/momentum_mod.f90

@@ -171,6 +171,7 @@ real(kind=4), allocatable :: Amatrix(:,:), rhs(:) !to be passed to lapack routin
 real(kind=ireals) :: tau_ix, tau_iy
 real(kind=ireals) :: Lxi, Lxj, Lyi, Lyj
 real(kind=ireals) :: tau_grx, tau_gry
+real(kind=ireals) :: velx_log, vely_log
 real(kind=ireals) :: time1, time2, totime
 real(kind=ireals) :: scross, hydr_radius, disch, speedav 
 real(kind=ireals), allocatable :: row0_(:)
@@ -330,8 +331,10 @@ endif
  !coef2 = g*Cz_gr**(-2)*sqrt(u1(M+1)**2+v1(M+1)**2) !*0.005 !*row0
 
  !Bottom momentum exchange coefficient from logarithmic profile
- coef2 = sqrt(u1(M+1)**2+v1(M+1)**2)*kappa**2/log(h1*0.25*ddz(M)/z0_bot)**2
- velfrict_bot = sqrt(u1(M+1)**2+v1(M+1)**2)*kappa/log(h1*0.25*ddz(M)/z0_bot)
+ velx_log = 0.5*(u1(M+1)+u1(M))
+ vely_log = 0.5*(v1(M+1)+v1(M))
+ coef2 = sqrt(velx_log**2+vely_log**2)*kappa**2/log(h1*0.5*ddz(M)/z0_bot)**2
+ velfrict_bot = sqrt(velx_log**2+vely_log**2)*kappa/log(h1*0.5*ddz(M)/z0_bot)
  tau_gr = row0*velfrict_bot**2
 
 ! Friction at the water - upper ice interface
@@ -1227,9 +1230,9 @@ bcbot : if (cuette%par == 0 .or. cuette%par == 11) then
   
   do i = 1, 2
     xx3(i) = horvisc_ON*0.5*pi*pi*horvisc%par*dt/LxLyCDL(M+1,i)**2 !horizontal viscosity for 1-st harmonic mode
-    cm_(M+1,i,i) = (xx + yy*(1 + xx3(i)) + coef2*ddz(M)*h1/k2(M) + &
+    cm_(M+1,i,i) = (xx + yy*(1 + xx3(i)) + 0.5*coef2*ddz(M)*h1/k2(M) + &
     & 0.5*yy*dt*a_veg*c_veg*sqrt(um(M+1)**2 + vm(M+1)**2)*veg_sw(M+1))
-    am_(M+1,i,i) = xx !- coef2*ddz(M)*h1/k2(M)
+    am_(M+1,i,i) = xx - 0.5*coef2*ddz(M)*h1/k2(M)
   enddo
   
   am_(M+1,1,2) = 0.