New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
zdfbfr.F90 in branches/2012/dev_r3452_NOCL06_LOGLAYER/NEMOGCM/NEMO/OPA_SRC/ZDF – NEMO

source: branches/2012/dev_r3452_NOCL06_LOGLAYER/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfbfr.F90 @ 3487

Last change on this file since 3487 was 3487, checked in by hliu, 11 years ago

Log Layer bottom friction option is added. The two relevant parameters to activate this option has been added in the AMM12 namelist file

  • Property svn:keywords set to Id
File size: 14.5 KB
Line 
1MODULE zdfbfr
2   !!======================================================================
3   !!                       ***  MODULE  zdfbfr  ***
4   !! Ocean physics: Bottom friction
5   !!======================================================================
6   !! History :  OPA  ! 1997-06  (G. Madec, A.-M. Treguier)  Original code
7   !!   NEMO     1.0  ! 2002-06  (G. Madec)  F90: Free form and module
8   !!            3.2  ! 2009-09  (A.C.Coward)  Correction to include barotropic contribution
9   !!            3.3  ! 2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase
10   !!            3.4  ! 2011-11  (H. Liu) implementation of semi-implicit bottom friction option
11   !!                 ! 2012-06  (H. Liu) implementation of Log Layer bottom friction option
12   !!----------------------------------------------------------------------
13
14   !!----------------------------------------------------------------------
15   !!   zdf_bfr      : update momentum Kz at the ocean bottom due to the type of bottom friction chosen
16   !!   zdf_bfr_init : read in namelist and control the bottom friction parameters.
17   !!   zdf_bfr_2d   : read in namelist and control the bottom friction parameters.
18   !!----------------------------------------------------------------------
19   USE oce             ! ocean dynamics and tracers variables
20   USE dom_oce         ! ocean space and time domain variables
21   USE zdf_oce         ! ocean vertical physics variables
22   USE in_out_manager  ! I/O manager
23   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
24   USE lib_mpp         ! distributed memory computing
25   USE prtctl          ! Print control
26   USE timing          ! Timing
27
28   USE phycst, ONLY: vkarmn
29
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC   zdf_bfr         ! called by step.F90
34   PUBLIC   zdf_bfr_init    ! called by opa.F90
35
36   !                                    !!* Namelist nambfr: bottom friction namelist *
37   INTEGER  ::   nn_bfr      = 0           ! = 0/1/2/3 type of bottom friction
38   REAL(wp) ::   rn_bfri1    = 4.0e-4_wp   ! bottom drag coefficient (linear case)
39   REAL(wp) ::   rn_bfri2    = 1.0e-3_wp   ! bottom drag coefficient (non linear case)
40   REAL(wp) ::   rn_bfeb2    = 2.5e-3_wp   ! background bottom turbulent kinetic energy  [m2/s2]
41   REAL(wp) ::   rn_bfrien   = 30._wp      ! local factor to enhance coefficient bfri
42   REAL(wp) ::   rn_bfrz0    = 0.003_wp    ! bottom roughness for loglayer bfr coeff
43   LOGICAL  ::   ln_bfr2d    = .false.     ! logical switch for 2D enhancement
44   LOGICAL  ::   ln_loglayer = .false.     ! switch for log layer bfr coeff.
45   LOGICAL , PUBLIC                            ::  ln_bfrimp   = .false.  ! switch for implicit bottom friction
46   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::  bfrcoef2d            ! 2D bottom drag coefficient
47
48   !! * Substitutions
49#  include "vectopt_loop_substitute.h90"
50#  include "domzgr_substitute.h90"
51   !!----------------------------------------------------------------------
52   !! NEMO/OPA 4.0 , NEMO Consortium (2011)
53   !! $Id$
54   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
55   !!----------------------------------------------------------------------
56CONTAINS
57
58   INTEGER FUNCTION zdf_bfr_alloc()
59      !!----------------------------------------------------------------------
60      !!                ***  FUNCTION zdf_bfr_alloc  ***
61      !!----------------------------------------------------------------------
62      ALLOCATE( bfrcoef2d(jpi,jpj), STAT=zdf_bfr_alloc )
63      !
64      IF( lk_mpp             )   CALL mpp_sum ( zdf_bfr_alloc )
65      IF( zdf_bfr_alloc /= 0 )   CALL ctl_warn('zdf_bfr_alloc: failed to allocate arrays.')
66   END FUNCTION zdf_bfr_alloc
67
68
69   SUBROUTINE zdf_bfr( kt )
70      !!----------------------------------------------------------------------
71      !!                   ***  ROUTINE zdf_bfr  ***
72      !!                 
73      !! ** Purpose :   compute the bottom friction coefficient.
74      !!
75      !! ** Method  :   Calculate and store part of the momentum trend due   
76      !!              to bottom friction following the chosen friction type
77      !!              (free-slip, linear, or quadratic). The component
78      !!              calculated here is multiplied by the bottom velocity in
79      !!              dyn_bfr to provide the trend term.
80      !!                The coefficients are updated at each time step only
81      !!              in the quadratic case.
82      !!
83      !! ** Action  :   bfrua , bfrva   bottom friction coefficients
84      !!----------------------------------------------------------------------
85      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
86      !!
87      INTEGER  ::   ji, jj       ! dummy loop indices
88      INTEGER  ::   ikbu, ikbv   ! local integers
89      REAL(wp) ::   zvu, zuv, zecu, zecv   ! temporary scalars
90      REAL(wp) ::   ztmp         ! temporary scalars
91      !!----------------------------------------------------------------------
92      !
93      IF( nn_timing == 1 )  CALL timing_start('zdf_bfr')
94      !
95      IF( nn_bfr == 2 ) THEN                 ! quadratic botton friction
96         ! Calculate and store the quadratic bottom friction coefficient bfrua and bfrva
97         ! where bfrUa = C_d*SQRT(u_bot^2 + v_bot^2 + e_b) {U=[u,v]}
98         ! from these the trend due to bottom friction:  -F_h/e3U  can be calculated
99         ! where -F_h/e3U_bot = bfrUa*Ub/e3U_bot {U=[u,v]}
100         !
101
102         IF(ln_loglayer) THEN       ! "log layer" bottom friction coefficient
103#  if defined key_vectopt_loop
104           DO jj = 1, 1
105             DO ji = 1, jpij   ! vector opt. (forced unrolling)
106#  else
107           DO jj = 1, jpj
108             DO ji = 1, jpi
109#  endif
110                ztmp = 0.5_wp * fse3t(ji,jj,mbkt(ji,jj))
111                ztmp = max(ztmp, rn_bfrz0)
112                bfrcoef2d(ji,jj) = ( log( ztmp / rn_bfrz0 ) / vkarmn ) ** (-2)
113             END DO
114           END DO
115         ENDIF
116
117# if defined key_vectopt_loop
118         DO jj = 1, 1
119!CDIR NOVERRCHK
120            DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling)
121# else
122!CDIR NOVERRCHK
123         DO jj = 2, jpjm1
124!CDIR NOVERRCHK
125            DO ji = 2, jpim1
126# endif
127               ikbu = mbku(ji,jj)         ! ocean bottom level at u- and v-points
128               ikbv = mbkv(ji,jj)         ! (deepest ocean u- and v-points)
129               !
130               zvu  = 0.25 * (  vn(ji,jj  ,ikbu) + vn(ji+1,jj  ,ikbu)     &
131                  &           + vn(ji,jj-1,ikbu) + vn(ji+1,jj-1,ikbu)  )
132               zuv  = 0.25 * (  un(ji,jj  ,ikbv) + un(ji-1,jj  ,ikbv)     &
133                  &           + un(ji,jj+1,ikbv) + un(ji-1,jj+1,ikbv)  )
134               !
135               zecu = SQRT(  un(ji,jj,ikbu) * un(ji,jj,ikbu) + zvu*zvu + rn_bfeb2  )
136               zecv = SQRT(  vn(ji,jj,ikbv) * vn(ji,jj,ikbv) + zuv*zuv + rn_bfeb2  )
137               !
138               bfrua(ji,jj) = - 0.5_wp * ( bfrcoef2d(ji,jj) + bfrcoef2d(ji+1,jj  ) ) * zecu 
139               bfrva(ji,jj) = - 0.5_wp * ( bfrcoef2d(ji,jj) + bfrcoef2d(ji  ,jj+1) ) * zecv
140            END DO
141         END DO
142
143
144         !
145         CALL lbc_lnk( bfrua, 'U', 1. )   ;   CALL lbc_lnk( bfrva, 'V', 1. )      ! Lateral boundary condition
146         !
147         IF(ln_ctl)   CALL prt_ctl( tab2d_1=bfrua, clinfo1=' bfr  - u: ', mask1=umask,        &
148            &                       tab2d_2=bfrva, clinfo2=       ' v: ', mask2=vmask,ovlap=1 )
149      ENDIF
150
151      !
152      IF( nn_timing == 1 )  CALL timing_stop('zdf_bfr')
153      !
154   END SUBROUTINE zdf_bfr
155
156
157   SUBROUTINE zdf_bfr_init
158      !!----------------------------------------------------------------------
159      !!                  ***  ROUTINE zdf_bfr_init  ***
160      !!                   
161      !! ** Purpose :   Initialization of the bottom friction
162      !!
163      !! ** Method  :   Read the nammbf namelist and check their consistency
164      !!              called at the first timestep (nit000)
165      !!----------------------------------------------------------------------
166      USE iom   ! I/O module for ehanced bottom friction file
167      !!
168      INTEGER   ::   inum         ! logical unit for enhanced bottom friction file
169      INTEGER   ::   ji, jj       ! dummy loop indexes
170      INTEGER   ::   ikbu, ikbv   ! temporary integers
171      INTEGER   ::   ictu, ictv   !    -          -
172      REAL(wp)  ::   zminbfr, zmaxbfr   ! temporary scalars
173      REAL(wp)  ::   zfru, zfrv         !    -         -
174      !!
175      NAMELIST/nambfr/ nn_bfr, rn_bfri1, rn_bfri2, rn_bfeb2, rn_bfrz0, ln_bfr2d, &
176                    &  rn_bfrien, ln_bfrimp, ln_loglayer
177      !!----------------------------------------------------------------------
178      !
179      IF( nn_timing == 1 )  CALL timing_start('zdf_bfr_init')
180      !
181      REWIND ( numnam )               !* Read Namelist nam_bfr : bottom momentum boundary condition
182      READ   ( numnam, nambfr )
183
184      !                               !* Parameter control and print
185      IF(lwp) WRITE(numout,*)
186      IF(lwp) WRITE(numout,*) 'zdf_bfr : momentum bottom friction'
187      IF(lwp) WRITE(numout,*) '~~~~~~~'
188      IF(lwp) WRITE(numout,*) '   Namelist nam_bfr : set bottom friction parameters'
189
190      !                              ! allocate zdfbfr arrays
191      IF( zdf_bfr_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'zdf_bfr_init : unable to allocate arrays' )
192
193      !                              ! Make sure ln_zdfexp=.false. when use implicit bfr
194      IF( ln_bfrimp .AND. ln_zdfexp ) THEN
195         IF(lwp) THEN
196            WRITE(numout,*)
197            WRITE(numout,*) 'Implicit bottom friction can only be used when ln_zdfexp=.false.'
198            WRITE(numout,*) '         but you set: ln_bfrimp=.true. and ln_zdfexp=.true.!!!!'
199            WRITE(ctmp1,*)  '         set either ln_zdfexp = .false or ln_bfrimp = .false.'
200            CALL ctl_stop( ctmp1 )
201         END IF
202      END IF
203
204      SELECT CASE (nn_bfr)
205      !
206      CASE( 0 )
207         IF(lwp) WRITE(numout,*) '      free-slip '
208         bfrua(:,:) = 0.e0
209         bfrva(:,:) = 0.e0
210         !
211      CASE( 1 )
212         IF(lwp) WRITE(numout,*) '      linear botton friction'
213         IF(lwp) WRITE(numout,*) '      friction coef.   rn_bfri1  = ', rn_bfri1
214         IF( ln_bfr2d ) THEN
215            IF(lwp) WRITE(numout,*) '      read coef. enhancement distribution from file   ln_bfr2d  = ', ln_bfr2d
216            IF(lwp) WRITE(numout,*) '      coef rn_bfri2 enhancement factor                rn_bfrien  = ',rn_bfrien
217         ENDIF
218         !
219         bfrcoef2d(:,:) = rn_bfri1  ! initialize bfrcoef2d to the namelist variable
220         !
221         IF(ln_bfr2d) THEN 
222            ! bfr_coef is a coefficient in [0,1] giving the mask where to apply the bfr enhancement
223            CALL iom_open('bfr_coef.nc',inum)
224            CALL iom_get (inum, jpdom_data, 'bfr_coef',bfrcoef2d,1) ! bfrcoef2d is used as tmp array
225            CALL iom_close(inum)
226            bfrcoef2d(:,:) = rn_bfri1 * ( 1 + rn_bfrien * bfrcoef2d(:,:) )
227         ENDIF
228         bfrua(:,:) = - bfrcoef2d(:,:)
229         bfrva(:,:) = - bfrcoef2d(:,:)
230         !
231      CASE( 2 )
232         IF(lwp) WRITE(numout,*) '      quadratic botton friction'
233         IF(lwp) WRITE(numout,*) '      friction coef.   rn_bfri2  = ', rn_bfri2
234         IF(lwp) WRITE(numout,*) '      background tke   rn_bfeb2  = ', rn_bfeb2
235         IF( ln_bfr2d ) THEN
236            IF(lwp) WRITE(numout,*) '      read coef. enhancement distribution from file   ln_bfr2d  = ', ln_bfr2d
237            IF(lwp) WRITE(numout,*) '      coef rn_bfri2 enhancement factor                rn_bfrien  = ',rn_bfrien
238         ENDIF
239         bfrcoef2d(:,:) = rn_bfri2  ! initialize bfrcoef2d to the namelist variable
240
241         !
242         IF(ln_bfr2d) THEN 
243            ! bfr_coef is a coefficient in [0,1] giving the mask where to apply the bfr enhancement
244            CALL iom_open('bfr_coef.nc',inum)
245            CALL iom_get (inum, jpdom_data, 'bfr_coef',bfrcoef2d,1) ! bfrcoef2d is used as tmp array
246            CALL iom_close(inum)
247            bfrcoef2d(:,:)= rn_bfri2 * ( 1 + rn_bfrien * bfrcoef2d(:,:) )
248         ENDIF
249         !
250      CASE DEFAULT
251         IF(lwp) WRITE(ctmp1,*) '         bad flag value for nn_bfr = ', nn_bfr
252         CALL ctl_stop( ctmp1 )
253         !
254      END SELECT
255      IF(lwp) WRITE(numout,*) '      implicit bottom friction switch                ln_bfrimp  = ', ln_bfrimp
256      !
257      ! Basic stability check on bottom friction coefficient
258      !
259      ictu = 0               ! counter for stability criterion breaches at U-pts
260      ictv = 0               ! counter for stability criterion breaches at V-pts
261      zminbfr =  1.e10_wp    ! initialise tracker for minimum of bottom friction coefficient
262      zmaxbfr = -1.e10_wp    ! initialise tracker for maximum of bottom friction coefficient
263      !
264#  if defined key_vectopt_loop
265      DO jj = 1, 1
266!CDIR NOVERRCHK
267         DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling)
268#  else
269!CDIR NOVERRCHK
270      DO jj = 2, jpjm1
271!CDIR NOVERRCHK
272         DO ji = 2, jpim1
273#  endif
274             ikbu = mbku(ji,jj)       ! deepest ocean level at u- and v-points
275             ikbv = mbkv(ji,jj)
276             zfru = 0.5 * fse3u(ji,jj,ikbu) / rdt
277             zfrv = 0.5 * fse3v(ji,jj,ikbv) / rdt
278             IF( ABS( bfrcoef2d(ji,jj) ) > zfru ) THEN
279                IF( ln_ctl ) THEN
280                   WRITE(numout,*) 'BFR ', narea, nimpp+ji, njmpp+jj, ikbu
281                   WRITE(numout,*) 'BFR ', ABS( bfrcoef2d(ji,jj) ), zfru
282                ENDIF
283                ictu = ictu + 1
284             ENDIF
285             IF( ABS( bfrcoef2d(ji,jj) ) > zfrv ) THEN
286                 IF( ln_ctl ) THEN
287                     WRITE(numout,*) 'BFR ', narea, nimpp+ji, njmpp+jj, ikbv
288                     WRITE(numout,*) 'BFR ', bfrcoef2d(ji,jj), zfrv
289                 ENDIF
290                 ictv = ictv + 1
291             ENDIF
292             zminbfr = MIN(  zminbfr, MIN( zfru, ABS( bfrcoef2d(ji,jj) ) )  )
293             zmaxbfr = MAX(  zmaxbfr, MIN( zfrv, ABS( bfrcoef2d(ji,jj) ) )  )
294         END DO
295      END DO
296      IF( lk_mpp ) THEN
297         CALL mpp_sum( ictu )
298         CALL mpp_sum( ictv )
299         CALL mpp_min( zminbfr )
300         CALL mpp_max( zmaxbfr )
301      ENDIF
302      IF( .NOT.ln_bfrimp) THEN
303      IF( lwp .AND. ictu + ictv > 0 ) THEN
304         WRITE(numout,*) ' Bottom friction stability check failed at ', ictu, ' U-points '
305         WRITE(numout,*) ' Bottom friction stability check failed at ', ictv, ' V-points '
306         WRITE(numout,*) ' Bottom friction coefficient now ranges from: ', zminbfr, ' to ', zmaxbfr
307         WRITE(numout,*) ' Bottom friction coefficient will be reduced where necessary'
308      ENDIF
309      ENDIF
310      !
311      IF( nn_timing == 1 )  CALL timing_stop('zdf_bfr_init')
312      !
313   END SUBROUTINE zdf_bfr_init
314
315   !!======================================================================
316END MODULE zdfbfr
Note: See TracBrowser for help on using the repository browser.