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/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/ZDF – NEMO

source: branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/ZDF/zdfbfr.F90 @ 2200

Last change on this file since 2200 was 2104, checked in by cetlod, 14 years ago

update DEV_r2006_merge_TRA_TRC according to review

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 12.0 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   !!----------------------------------------------------------------------
11
12   !!----------------------------------------------------------------------
13   !!   zdf_bfr      : update momentum Kz at the ocean bottom due to the type of bottom friction chosen
14   !!   zdf_bfr_init : read in namelist and control the bottom friction parameters.
15   !!   zdf_bfr_2d   : read in namelist and control the bottom friction parameters.
16   !!----------------------------------------------------------------------
17   USE oce             ! ocean dynamics and tracers variables
18   USE dom_oce         ! ocean space and time domain variables
19   USE zdf_oce         ! ocean vertical physics variables
20   USE in_out_manager  ! I/O manager
21   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
22   USE lib_mpp         ! distributed memory computing
23   USE prtctl          ! Print control
24
25   IMPLICIT NONE
26   PRIVATE
27
28   PUBLIC   zdf_bfr         ! called by step.F90
29   PUBLIC   zdf_bfr_init    ! called by opa.F90
30   
31   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   bfrua , bfrva   !: Bottom friction coefficients set in zdfbfr
32
33   !                                    !!* Namelist nambfr: bottom friction namelist *
34   INTEGER  ::   nn_bfr    = 0           ! = 0/1/2/3 type of bottom friction
35   REAL(wp) ::   rn_bfri1  = 4.0e-4_wp   ! bottom drag coefficient (linear case)
36   REAL(wp) ::   rn_bfri2  = 1.0e-3_wp   ! bottom drag coefficient (non linear case)
37   REAL(wp) ::   rn_bfeb2  = 2.5e-3_wp   ! background bottom turbulent kinetic energy  [m2/s2]
38   REAL(wp) ::   rn_bfrien = 30_wp       ! local factor to enhance coefficient bfri
39   LOGICAL  ::   ln_bfr2d  = .false.     ! logical switch for 2D enhancement
40   
41   REAL(wp), DIMENSION(jpi,jpj) ::   bfrcoef2d = 1.e-3_wp   ! 2D bottom drag coefficient
42
43   !! * Substitutions
44#  include "domzgr_substitute.h90"
45   !!----------------------------------------------------------------------
46   !! NEMO/OPA 3,3 , LOCEAN-IPSL (2010)
47   !! $Id$
48   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
49   !!----------------------------------------------------------------------
50
51CONTAINS
52
53   SUBROUTINE zdf_bfr( kt )
54      !!----------------------------------------------------------------------
55      !!                   ***  ROUTINE zdf_bfr  ***
56      !!                 
57      !! ** Purpose :   compute the bottom friction coefficient.
58      !!
59      !! ** Method  :   Calculate and store part of the momentum trend due   
60      !!              to bottom friction following the chosen friction type
61      !!              (free-slip, linear, or quadratic). The component
62      !!              calculated here is multiplied by the bottom velocity in
63      !!              dyn_bfr to provide the trend term.
64      !!                The coefficients are updated at each time step only
65      !!              in the quadratic case.
66      !!
67      !! ** Action  :   bfrua , bfrva   bottom friction coefficients
68      !!----------------------------------------------------------------------
69      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
70      !!
71      INTEGER  ::   ji, jj         ! dummy loop indices
72      INTEGER  ::   ikbu, ikbum1   ! temporary integers
73      INTEGER  ::   ikbv, ikbvm1   !    -          -
74      REAL(wp) ::   zvu, zuv, zecu, zecv   ! temporary scalars
75      !!----------------------------------------------------------------------
76
77      IF( nn_bfr == 2 ) THEN                 ! quadratic botton friction
78         ! Calculate and store the quadratic bottom friction coefficient bfrua and bfrva
79         ! where bfrUa = C_d*SQRT(u_bot^2 + v_bot^2 + e_b) {U=[u,v]}
80         ! from these the trend due to bottom friction:  -F_h/e3U  can be calculated
81         ! where -F_h/e3U_bot = bfrUa*Ub/e3U_bot {U=[u,v]}
82         !
83# if defined key_vectopt_loop
84         DO jj = 1, 1
85!CDIR NOVERRCHK
86            DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling)
87# else
88!CDIR NOVERRCHK
89         DO jj = 2, jpjm1
90!CDIR NOVERRCHK
91            DO ji = 2, jpim1
92# endif
93               ikbu   = MIN( mbathy(ji+1,jj  ), mbathy(ji,jj) )
94               ikbv   = MIN( mbathy(ji  ,jj+1), mbathy(ji,jj) )
95               ikbum1 = MAX( ikbu-1, 1 )
96               ikbvm1 = MAX( ikbv-1, 1 )
97               !
98               zvu  = 0.25 * (  vn(ji,jj  ,ikbum1) + vn(ji+1,jj  ,ikbum1)     &
99                  &           + vn(ji,jj-1,ikbum1) + vn(ji+1,jj-1,ikbum1)  )
100               zuv  = 0.25 * (  un(ji,jj  ,ikbvm1) + un(ji-1,jj  ,ikbvm1)     &
101                  &           + un(ji,jj+1,ikbvm1) + un(ji-1,jj+1,ikbvm1)  )
102               !
103               zecu = SQRT(  un(ji,jj,ikbum1) * un(ji,jj,ikbum1) + zvu*zvu + rn_bfeb2  )
104               zecv = SQRT(  vn(ji,jj,ikbvm1) * vn(ji,jj,ikbvm1) + zuv*zuv + rn_bfeb2  )
105               !
106               bfrua(ji,jj) = - 0.5 * ( bfrcoef2d(ji,jj) + bfrcoef2d(ji+1,jj  ) ) * zecu 
107               bfrva(ji,jj) = - 0.5 * ( bfrcoef2d(ji,jj) + bfrcoef2d(ji  ,jj+1) ) * zecv
108            END DO
109         END DO
110         !
111         CALL lbc_lnk( bfrua, 'U', 1. )   ;   CALL lbc_lnk( bfrva, 'V', 1. )      ! Lateral boundary condition
112         !
113         IF(ln_ctl)   CALL prt_ctl( tab2d_1=bfrua, clinfo1=' bfr  - u: ', mask1=umask,        &
114            &                       tab2d_2=bfrva, clinfo2=       ' v: ', mask2=vmask,ovlap=1 )
115      ENDIF
116      !
117   END SUBROUTINE zdf_bfr
118
119
120   SUBROUTINE zdf_bfr_init
121      !!----------------------------------------------------------------------
122      !!                  ***  ROUTINE zdf_bfr_init  ***
123      !!                   
124      !! ** Purpose :   Initialization of the bottom friction
125      !!
126      !! ** Method  :   Read the nammbf namelist and check their consistency
127      !!              called at the first timestep (nit000)
128      !!----------------------------------------------------------------------
129      USE iom   ! I/O module for ehanced bottom friction file
130      !!
131      INTEGER ::   inum         ! logical unit for enhanced bottom friction file
132      INTEGER ::   ji, jj       ! dummy loop indexes
133      INTEGER ::   ikbu, ikbv   ! temporary integers
134      INTEGER ::   ictu, ictv   !    -          -
135      REAL(wp) ::  zminbfr, zmaxbfr   ! temporary scalars
136      REAL(wp) ::  zfru, zfrv         !    -         -
137      !!
138      NAMELIST/nambfr/ nn_bfr, rn_bfri1, rn_bfri2, rn_bfeb2, ln_bfr2d, rn_bfrien
139      !!----------------------------------------------------------------------
140
141      REWIND ( numnam )               !* Read Namelist nam_bfr : bottom momentum boundary condition
142      READ   ( numnam, nambfr )
143
144      !                               !* Parameter control and print
145      IF(lwp) WRITE(numout,*)
146      IF(lwp) WRITE(numout,*) 'zdf_bfr : momentum bottom friction'
147      IF(lwp) WRITE(numout,*) '~~~~~~~'
148      IF(lwp) WRITE(numout,*) '   Namelist nam_bfr : set bottom friction parameters'
149
150      SELECT CASE (nn_bfr)
151
152      CASE( 0 )
153         IF(lwp) WRITE(numout,*) '      free-slip '
154         bfrua(:,:) = 0.e0
155         bfrva(:,:) = 0.e0
156         !
157      CASE( 1 )
158         IF(lwp) WRITE(numout,*) '      linear botton friction'
159         IF(lwp) WRITE(numout,*) '      friction coef.   rn_bfri1  = ', rn_bfri1
160         IF( ln_bfr2d ) THEN
161            IF(lwp) WRITE(numout,*) '      read coef. enhancement distribution from file   ln_bfr2d  = ', ln_bfr2d
162            IF(lwp) WRITE(numout,*) '      coef rn_bfri2 enhancement factor                rn_bfrien  = ',rn_bfrien
163         ENDIF
164         !
165         bfrcoef2d(:,:) = rn_bfri1  ! initialize bfrcoef2d to the namelist variable
166         !
167         IF(ln_bfr2d) THEN 
168            ! bfr_coef is a coefficient in [0,1] giving the mask where to apply the bfr enhancement
169            CALL iom_open('bfr_coef.nc',inum)
170            CALL iom_get (inum, jpdom_data, 'bfr_coef',bfrcoef2d,1) ! bfrcoef2d is used as tmp array
171            CALL iom_close(inum)
172            bfrcoef2d(:,:)= rn_bfri1 * ( 1 + rn_bfrien * bfrcoef2d(:,:) )
173         ENDIF
174         bfrua(:,:) = - bfrcoef2d(:,:)
175         bfrva(:,:) = - bfrcoef2d(:,:)
176         !
177      CASE( 2 )
178         IF(lwp) WRITE(numout,*) '      quadratic botton friction'
179         IF(lwp) WRITE(numout,*) '      friction coef.   rn_bfri2  = ', rn_bfri2
180         IF(lwp) WRITE(numout,*) '      background tke   rn_bfeb2  = ', rn_bfeb2
181         IF( ln_bfr2d ) THEN
182            IF(lwp) WRITE(numout,*) '      read coef. enhancement distribution from file   ln_bfr2d  = ', ln_bfr2d
183            IF(lwp) WRITE(numout,*) '      coef rn_bfri2 enhancement factor                rn_bfrien  = ',rn_bfrien
184         ENDIF
185         bfrcoef2d(:,:) = rn_bfri2  ! initialize bfrcoef2d to the namelist variable
186         !
187         IF(ln_bfr2d) THEN 
188            ! bfr_coef is a coefficient in [0,1] giving the mask where to apply the bfr enhancement
189            CALL iom_open('bfr_coef.nc',inum)
190            CALL iom_get (inum, jpdom_data, 'bfr_coef',bfrcoef2d,1) ! bfrcoef2d is used as tmp array
191            CALL iom_close(inum)
192            bfrcoef2d(:,:)= rn_bfri2 * ( 1 + rn_bfrien * bfrcoef2d(:,:) )
193         ENDIF
194         !
195      CASE DEFAULT
196         IF(lwp) WRITE(ctmp1,*) '         bad flag value for nn_bfr = ', nn_bfr
197         CALL ctl_stop( ctmp1 )
198         !
199      END SELECT
200      !
201      ! Basic stability check on bottom friction coefficient
202      !
203      ictu = 0               ! counter for stability criterion breaches at U-pts
204      ictv = 0               ! counter for stability criterion breaches at V-pts
205      zminbfr =  1.e10_wp    ! initialise tracker for minimum of bottom friction coefficient
206      zmaxbfr = -1.e10_wp    ! initialise tracker for maximum of bottom friction coefficient
207      !
208#  if defined key_vectopt_loop
209      DO jj = 1, 1
210!CDIR NOVERRCHK
211         DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling)
212#  else
213!CDIR NOVERRCHK
214      DO jj = 2, jpjm1
215!CDIR NOVERRCHK
216         DO ji = 2, jpim1
217#  endif
218             ikbu = MIN( mbathy(ji+1,jj  ), mbathy(ji,jj) )
219             ikbv = MIN( mbathy(ji  ,jj+1), mbathy(ji,jj) )
220             zfru = 0.5 * fse3u(ji,jj,ikbu) / rdt
221             zfrv = 0.5 * fse3v(ji,jj,ikbv) / rdt
222             IF( ABS( bfrcoef2d(ji,jj) ) > zfru ) THEN
223                IF( ln_ctl ) THEN
224                   WRITE(numout,*) 'BFR ',narea,nimpp+ji,njmpp+jj,ikbu
225                   WRITE(numout,*) 'BFR ',ABS( bfrcoef2d(ji,jj) ), zfru
226                ENDIF
227                ictu = ictu + 1
228             ENDIF
229             IF( ABS( bfrcoef2d(ji,jj) ) > zfrv ) THEN
230                 IF( ln_ctl ) THEN
231                     WRITE(numout,*) 'BFR ', narea, nimpp+ji, njmpp+jj, ikbv
232                     WRITE(numout,*) 'BFR ', bfrcoef2d(ji,jj), zfrv
233                 ENDIF
234                 ictv = ictv + 1
235             ENDIF
236             zminbfr = MIN(  zminbfr, MIN( zfru, ABS( bfrcoef2d(ji,jj) ) )  )
237             zmaxbfr = MAX(  zmaxbfr, MIN( zfrv, ABS( bfrcoef2d(ji,jj) ) )  )
238         END DO
239      END DO
240      IF( lk_mpp ) THEN
241         CALL mpp_sum( ictu )
242         CALL mpp_sum( ictv )
243         CALL mpp_min( zminbfr )
244         CALL mpp_max( zmaxbfr )
245      ENDIF
246      IF( lwp .AND. ictu + ictv .GT. 0 ) THEN
247         WRITE(numout,*) ' Bottom friction stability check failed at ', ictu, ' U-points '
248         WRITE(numout,*) ' Bottom friction stability check failed at ', ictv, ' V-points '
249         WRITE(numout,*) ' Bottom friction coefficient now ranges from: ', zminbfr, ' to ', zmaxbfr
250         WRITE(numout,*) ' Bottom friction coefficient will be reduced where necessary'
251      ENDIF
252      !
253   END SUBROUTINE zdf_bfr_init
254
255   !!======================================================================
256END MODULE zdfbfr
Note: See TracBrowser for help on using the repository browser.