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

source: trunk/NEMO/OPA_SRC/ZDF/zdfbfr.F90 @ 1671

Last change on this file since 1671 was 1671, checked in by rblod, 15 years ago

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