Bottom drag discretization corrected

source/model/momentum_mod.f90

@@ -449,7 +449,7 @@ SELECT CASE (Turbpar%par)
       endif
     enddo
     do i = 1, M
-      k2(i) = 0.00632 !max(CE(i)*E1(i)**2/(eps1(i) + ACCk),min_visc) + niu_wat + knum(i) !!
+      k2(i) = max(CE(i)*E1(i)**2/(eps1(i) + ACCk),min_visc) + niu_wat + knum(i) !!
     end do
 
 ! 3. Nickuradze (NICK) formulation: Rodi (1993)
@@ -1232,9 +1232,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.
@@ -1345,7 +1345,7 @@ do k = 1, M+1
 enddo
 
 
-!print*, taux
+!print*, 0.5*(u1(M+1)+u1(M))
 !print*, tauy
 
 !print*, 'Runtime for seiche model: ', totime

source/model/turb_mod.f90

@@ -991,7 +991,7 @@ allocate (rhotemp(1:M), rhosal(1:M))
   endif
 
   do i = 1, M
-    KT(i) = 0.00632*cw_m_row0 !max(CEt(i)*E2(i)**2/(eps2(i) + ACCk),min_diff)*cw_m_row0 !E2, eps2 ! !
+    KT(i) = max(CEt(i)*E2(i)**2/(eps2(i) + ACCk),min_diff)*cw_m_row0 !E2, eps2 ! !
   enddo
 
   ! Diagnostic calculation