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/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/ZDF – NEMO

source: branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfbfr.F90 @ 3211

Last change on this file since 3211 was 3211, checked in by spickles2, 12 years ago

Stephen Pickles, 11 Dec 2011

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