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.
zdfdrg.F90 in NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/src/OCE/ZDF – NEMO

source: NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/src/OCE/ZDF/zdfdrg.F90 @ 11671

Last change on this file since 11671 was 11671, checked in by acc, 5 years ago

Branch 2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles. Final, non-substantive changes to complete this branch. These changes remove all REWIND statements on the old namelist fortran units (now character variables for internal files). These changes have been left until last since they are easily repeated via a script and it may be preferable to use the previous revision for merge purposes and reapply these last changes separately. This branch has been fully SETTE tested.
  • Property svn:keywords set to Id
File size: 24.6 KB
Line 
1MODULE zdfdrg
2   !!======================================================================
3   !!                       ***  MODULE  zdfdrg  ***
4   !! Ocean physics: top and/or 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   !!            4.0  ! 2017-05  (G. Madec) zdfbfr becomes zdfdrg + variable names change
13   !!                                     + drag defined at t-point + new user interface + top drag (ocean cavities)
14   !!----------------------------------------------------------------------
15
16   !!----------------------------------------------------------------------
17   !!   zdf_drg       : update bottom friction coefficient (non-linear bottom friction only)
18   !!   zdf_drg_exp   : compute the top & bottom friction in explicit case
19   !!   zdf_drg_init  : read in namdrg namelist and control the bottom friction parameters.
20   !!       drg_init  :
21   !!----------------------------------------------------------------------
22   USE oce            ! ocean dynamics and tracers variables
23   USE phycst  , ONLY : vkarmn
24   USE dom_oce        ! ocean space and time domain variables
25   USE zdf_oce        ! ocean vertical physics variables
26   USE trd_oce        ! trends: ocean variables
27   USE trddyn         ! trend manager: dynamics
28   !
29   USE in_out_manager ! I/O manager
30   USE iom            ! I/O module
31   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
32   USE lib_mpp        ! distributed memory computing
33   USE prtctl         ! Print control
34
35   IMPLICIT NONE
36   PRIVATE
37
38   PUBLIC   zdf_drg         ! called by zdf_phy
39   PUBLIC   zdf_drg_exp     ! called by dyn_zdf
40   PUBLIC   zdf_drg_init    ! called by zdf_phy_init
41
42   !                                 !!* Namelist namdrg: nature of drag coefficient namelist *
43   LOGICAL          ::   ln_OFF       ! free-slip       : Cd = 0
44   LOGICAL          ::   ln_lin       !     linear  drag: Cd = Cd0_lin
45   LOGICAL          ::   ln_non_lin   ! non-linear  drag: Cd = Cd0_nl |U|
46   LOGICAL          ::   ln_loglayer  ! logarithmic drag: Cd = vkarmn/log(z/z0)
47   LOGICAL , PUBLIC ::   ln_drgimp    ! implicit top/bottom friction flag
48
49   !                                 !!* Namelist namdrg_top & _bot: TOP or BOTTOM coefficient namelist *
50   REAL(wp)         ::   rn_Cd0       !: drag coefficient                                           [ - ]
51   REAL(wp)         ::   rn_Uc0       !: characteristic velocity (linear case: tau=rho*Cd0*Uc0*u)   [m/s]
52   REAL(wp)         ::   rn_Cdmax     !: drag value maximum (ln_loglayer=T)                         [ - ]
53   REAL(wp)         ::   rn_z0        !: roughness          (ln_loglayer=T)                         [ m ]
54   REAL(wp)         ::   rn_ke0       !: background kinetic energy (non-linear case)                [m2/s2]
55   LOGICAL          ::   ln_boost     !: =T regional boost of Cd0 ; =F Cd0 horizontally uniform
56   REAL(wp)         ::   rn_boost     !: local boost factor                                         [ - ]
57
58   REAL(wp), PUBLIC ::   r_Cdmin_top, r_Cdmax_top, r_z0_top, r_ke0_top   ! set from namdrg_top namelist values
59   REAL(wp), PUBLIC ::   r_Cdmin_bot, r_Cdmax_bot, r_z0_bot, r_ke0_bot   !  -    -  namdrg_bot    -       -
60
61   INTEGER ::              ndrg       ! choice of the type of drag coefficient
62   !                                  ! associated indices:
63   INTEGER, PARAMETER ::   np_OFF      = 0   ! free-slip: drag set to zero
64   INTEGER, PARAMETER ::   np_lin      = 1   !     linear drag: Cd = Cd0_lin
65   INTEGER, PARAMETER ::   np_non_lin  = 2   ! non-linear drag: Cd = Cd0_nl |U|
66   INTEGER, PARAMETER ::   np_loglayer = 3   ! non linear drag (logarithmic formulation): Cd = vkarmn/log(z/z0)
67
68   LOGICAL , PUBLIC ::   l_zdfdrg           !: flag to update at each time step the top/bottom Cd
69   LOGICAL          ::   l_log_not_linssh   !: flag to update at each time step the position ot the velocity point
70   !
71   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC ::   rCd0_top, rCd0_bot   !: precomputed top/bottom drag coeff. at t-point (>0)
72   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC ::   rCdU_top, rCdU_bot   !: top/bottom drag coeff. at t-point (<0)  [m/s]
73
74   !! * Substitutions
75#  include "vectopt_loop_substitute.h90"
76   !!----------------------------------------------------------------------
77   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
78   !! $Id$
79   !! Software governed by the CeCILL license (see ./LICENSE)
80   !!----------------------------------------------------------------------
81CONTAINS
82
83   SUBROUTINE zdf_drg( kt, k_mk, pCdmin, pCdmax, pz0, pke0, pCd0,   &   ! <<== in
84      &                                                     pCdU )      ! ==>> out : bottom drag [m/s]
85      !!----------------------------------------------------------------------
86      !!                   ***  ROUTINE zdf_drg  ***
87      !!
88      !! ** Purpose :   update the top/bottom drag coefficient (non-linear case only)
89      !!
90      !! ** Method  :   In non linear friction case, the drag coeficient is
91      !!              a function of the velocity:
92      !!                          Cd = cd0 * |U+Ut|   
93      !!              where U is the top or bottom velocity and
94      !!                    Ut a tidal velocity (Ut^2 = Tidal kinetic energy
95      !!                       assumed here here to be constant)
96      !!              Depending on the input variable, the top- or bottom drag is compted
97      !!
98      !! ** Action  :   p_Cd   drag coefficient at t-point
99      !!----------------------------------------------------------------------
100      INTEGER                 , INTENT(in   ) ::   kt       ! ocean time-step index
101      !                       !               !!         !==  top or bottom variables  ==!
102      INTEGER , DIMENSION(:,:), INTENT(in   ) ::   k_mk     ! wet level (1st or last)
103      REAL(wp)                , INTENT(in   ) ::   pCdmin   ! min drag value
104      REAL(wp)                , INTENT(in   ) ::   pCdmax   ! max drag value
105      REAL(wp)                , INTENT(in   ) ::   pz0      ! roughness
106      REAL(wp)                , INTENT(in   ) ::   pke0     ! background tidal KE
107      REAL(wp), DIMENSION(:,:), INTENT(in   ) ::   pCd0     ! masked precomputed part of Cd0
108      REAL(wp), DIMENSION(:,:), INTENT(  out) ::   pCdU     ! = - Cd*|U|   (t-points) [m/s]
109      !!
110      INTEGER ::   ji, jj   ! dummy loop indices
111      INTEGER ::   imk      ! local integers
112      REAL(wp)::   zzz, zut, zvt, zcd   ! local scalars
113      !!----------------------------------------------------------------------
114      !
115      IF( l_log_not_linssh ) THEN     !==  "log layer"  ==!   compute Cd and -Cd*|U|
116         DO jj = 2, jpjm1
117            DO ji = 2, jpim1
118               imk = k_mk(ji,jj)          ! ocean bottom level at t-points
119               zut = un(ji,jj,imk) + un(ji-1,jj,imk)     ! 2 x velocity at t-point
120               zvt = vn(ji,jj,imk) + vn(ji,jj-1,imk)
121               zzz = 0.5_wp * e3t_n(ji,jj,imk)           ! altitude below/above (top/bottom) the boundary
122               !
123!!JC: possible WAD implementation should modify line below if layers vanish
124               zcd = (  vkarmn / LOG( zzz / pz0 )  )**2
125               zcd = pCd0(ji,jj) * MIN(  MAX( pCdmin , zcd ) , pCdmax  )   ! here pCd0 = mask*boost
126               pCdU(ji,jj) = - zcd * SQRT(  0.25 * ( zut*zut + zvt*zvt ) + pke0  )
127            END DO
128         END DO
129      ELSE                                            !==  standard Cd  ==!
130         DO jj = 2, jpjm1
131            DO ji = 2, jpim1
132               imk = k_mk(ji,jj)    ! ocean bottom level at t-points
133               zut = un(ji,jj,imk) + un(ji-1,jj,imk)     ! 2 x velocity at t-point
134               zvt = vn(ji,jj,imk) + vn(ji,jj-1,imk)
135               !                                                           ! here pCd0 = mask*boost * drag
136               pCdU(ji,jj) = - pCd0(ji,jj) * SQRT(  0.25 * ( zut*zut + zvt*zvt ) + pke0  )
137            END DO
138         END DO
139      ENDIF
140      !
141      IF(ln_ctl)   CALL prt_ctl( tab2d_1=pCdU, clinfo1=' Cd*U ')
142      !
143   END SUBROUTINE zdf_drg
144
145
146   SUBROUTINE zdf_drg_exp( kt, pub, pvb, pua, pva )
147      !!----------------------------------------------------------------------
148      !!                  ***  ROUTINE zdf_drg_exp  ***
149      !!
150      !! ** Purpose :   compute and add the explicit top and bottom frictions.
151      !!
152      !! ** Method  :   in explicit case,
153      !!
154      !!              NB: in implicit case the calculation is performed in dynzdf.F90
155      !!
156      !! ** Action  :   (pua,pva)   momentum trend increased by top & bottom friction trend
157      !!---------------------------------------------------------------------
158      INTEGER                         , INTENT(in   ) ::   kt         ! ocean time-step index
159      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pub, pvb   ! the two components of the before velocity
160      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pua, pva   ! the two components of the velocity tendency
161      !!
162      INTEGER  ::   ji, jj       ! dummy loop indexes
163      INTEGER  ::   ikbu, ikbv   ! local integers
164      REAL(wp) ::   zm1_2dt      ! local scalar
165      REAL(wp) ::   zCdu, zCdv   !   -      -
166      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrdu, ztrdv
167      !!---------------------------------------------------------------------
168      !
169!!gm bug : time step is only rdt (not 2 rdt if euler start !)
170      zm1_2dt = - 1._wp / ( 2._wp * rdt )
171
172      IF( l_trddyn ) THEN      ! trends: store the input trends
173         ALLOCATE( ztrdu(jpi,jpj,jpk) , ztrdv(jpi,jpj,jpk) )
174         ztrdu(:,:,:) = pua(:,:,:)
175         ztrdv(:,:,:) = pva(:,:,:)
176      ENDIF
177
178      DO jj = 2, jpjm1
179         DO ji = 2, jpim1
180            ikbu = mbku(ji,jj)          ! deepest wet ocean u- & v-levels
181            ikbv = mbkv(ji,jj)
182            !
183            ! Apply stability criteria on absolute value  : abs(bfr/e3) < 1/(2dt) => bfr/e3 > -1/(2dt)
184            zCdu = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / e3u_n(ji,jj,ikbu)
185            zCdv = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / e3v_n(ji,jj,ikbv)
186            !
187            pua(ji,jj,ikbu) = pua(ji,jj,ikbu) + MAX(  zCdu , zm1_2dt  ) * pub(ji,jj,ikbu)
188            pva(ji,jj,ikbv) = pva(ji,jj,ikbv) + MAX(  zCdv , zm1_2dt  ) * pvb(ji,jj,ikbv)
189         END DO
190      END DO
191      !
192      IF( ln_isfcav ) THEN        ! ocean cavities
193         DO jj = 2, jpjm1
194            DO ji = 2, jpim1
195               ikbu = miku(ji,jj)          ! first wet ocean u- & v-levels
196               ikbv = mikv(ji,jj)
197               !
198               ! Apply stability criteria on absolute value  : abs(bfr/e3) < 1/(2dt) => bfr/e3 > -1/(2dt)
199               zCdu = 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / e3u_n(ji,jj,ikbu)    ! NB: Cdtop masked
200               zCdv = 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / e3v_n(ji,jj,ikbv)
201               !
202               pua(ji,jj,ikbu) = pua(ji,jj,ikbu) + MAX(  zCdu , zm1_2dt  ) * pub(ji,jj,ikbu)
203               pva(ji,jj,ikbv) = pva(ji,jj,ikbv) + MAX(  zCdv , zm1_2dt  ) * pvb(ji,jj,ikbv)
204           END DO
205         END DO
206      ENDIF
207      !
208      IF( l_trddyn ) THEN      ! trends: send trends to trddyn for further diagnostics
209         ztrdu(:,:,:) = pua(:,:,:) - ztrdu(:,:,:)
210         ztrdv(:,:,:) = pva(:,:,:) - ztrdv(:,:,:)
211         CALL trd_dyn( ztrdu(:,:,:), ztrdv(:,:,:), jpdyn_bfr, kt )
212         DEALLOCATE( ztrdu, ztrdv )
213      ENDIF
214      !                                          ! print mean trends (used for debugging)
215      IF(ln_ctl)   CALL prt_ctl( tab3d_1=pua, clinfo1=' bfr  - Ua: ', mask1=umask,               &
216         &                       tab3d_2=pva, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
217      !
218   END SUBROUTINE zdf_drg_exp
219
220
221   SUBROUTINE zdf_drg_init
222      !!----------------------------------------------------------------------
223      !!                  ***  ROUTINE zdf_brg_init  ***
224      !!
225      !! ** Purpose :   Initialization of the bottom friction
226      !!
227      !! ** Method  :   Read the namdrg namelist and check their consistency
228      !!                called at the first timestep (nit000)
229      !!----------------------------------------------------------------------
230      INTEGER   ::   ji, jj      ! dummy loop indexes
231      INTEGER   ::   ios, ioptio   ! local integers
232      !!
233      NAMELIST/namdrg/ ln_OFF, ln_lin, ln_non_lin, ln_loglayer, ln_drgimp
234      !!----------------------------------------------------------------------
235      !
236      !                     !==  drag nature  ==!
237      !
238      READ  ( numnam_ref, namdrg, IOSTAT = ios, ERR = 901)
239901   IF( ios /= 0 )   CALL ctl_nam( ios , 'namdrg in reference namelist' )
240      READ  ( numnam_cfg, namdrg, IOSTAT = ios, ERR = 902 )
241902   IF( ios >  0 )   CALL ctl_nam( ios , 'namdrg in configuration namelist' )
242      IF(lwm) WRITE ( numond, namdrg )
243      !
244      IF(lwp) THEN
245         WRITE(numout,*)
246         WRITE(numout,*) 'zdf_drg_init : top and/or bottom drag setting'
247         WRITE(numout,*) '~~~~~~~~~~~~'
248         WRITE(numout,*) '   Namelist namdrg : top/bottom friction choices'
249         WRITE(numout,*) '      free-slip       : Cd = 0                  ln_OFF      = ', ln_OFF 
250         WRITE(numout,*) '      linear  drag    : Cd = Cd0                ln_lin      = ', ln_lin
251         WRITE(numout,*) '      non-linear  drag: Cd = Cd0_nl |U|         ln_non_lin  = ', ln_non_lin
252         WRITE(numout,*) '      logarithmic drag: Cd = vkarmn/log(z/z0)   ln_loglayer = ', ln_loglayer
253         WRITE(numout,*) '      implicit friction                         ln_drgimp   = ', ln_drgimp
254      ENDIF
255      !
256      ioptio = 0                       ! set ndrg and control check
257      IF( ln_OFF      ) THEN   ;   ndrg = np_OFF        ;   ioptio = ioptio + 1   ;   ENDIF
258      IF( ln_lin      ) THEN   ;   ndrg = np_lin        ;   ioptio = ioptio + 1   ;   ENDIF
259      IF( ln_non_lin  ) THEN   ;   ndrg = np_non_lin    ;   ioptio = ioptio + 1   ;   ENDIF
260      IF( ln_loglayer ) THEN   ;   ndrg = np_loglayer   ;   ioptio = ioptio + 1   ;   ENDIF
261      !
262      IF( ioptio /= 1 )   CALL ctl_stop( 'zdf_drg_init: Choose ONE type of drag coef in namdrg' )
263      !
264      !
265      !                     !==  BOTTOM drag setting  ==!   (applied at seafloor)
266      !
267      ALLOCATE( rCd0_bot(jpi,jpj), rCdU_bot(jpi,jpj) )
268      CALL drg_init( 'BOTTOM'   , mbkt       ,                                         &   ! <== in
269         &           r_Cdmin_bot, r_Cdmax_bot, r_z0_bot, r_ke0_bot, rCd0_bot, rCdU_bot )   ! ==> out
270
271      !
272      !                     !==  TOP drag setting  ==!   (applied at the top of ocean cavities)
273      !
274      IF( ln_isfcav ) THEN              ! Ocean cavities: top friction setting
275         ALLOCATE( rCd0_top(jpi,jpj), rCdU_top(jpi,jpj) )
276         CALL drg_init( 'TOP   '   , mikt       ,                                         &   ! <== in
277            &           r_Cdmin_top, r_Cdmax_top, r_z0_top, r_ke0_top, rCd0_top, rCdU_top )   ! ==> out
278      ENDIF
279      !
280   END SUBROUTINE zdf_drg_init
281
282
283   SUBROUTINE drg_init( cd_topbot, k_mk,  &
284      &                 pCdmin, pCdmax, pz0, pke0, pCd0, pCdU ) 
285      !!----------------------------------------------------------------------
286      !!                  ***  ROUTINE drg_init  ***
287      !!
288      !! ** Purpose :   Initialization of the top/bottom friction CdO and Cd
289      !!              from namelist parameters
290      !!----------------------------------------------------------------------
291      CHARACTER(len=6)        , INTENT(in   ) ::   cd_topbot       ! top/ bot indicator
292      INTEGER , DIMENSION(:,:), INTENT(in   ) ::   k_mk            ! 1st/last  wet level
293      REAL(wp)                , INTENT(  out) ::   pCdmin, pCdmax  ! min and max drag coef. [-]
294      REAL(wp)                , INTENT(  out) ::   pz0             ! roughness              [m]
295      REAL(wp)                , INTENT(  out) ::   pke0            ! background KE          [m2/s2]
296      REAL(wp), DIMENSION(:,:), INTENT(  out) ::   pCd0            ! masked precomputed part of the non-linear drag coefficient
297      REAL(wp), DIMENSION(:,:), INTENT(  out) ::   pCdU            ! minus linear drag*|U| at t-points  [m/s]
298      !!
299      CHARACTER(len=40) ::   cl_namdrg, cl_file, cl_varname, cl_namref, cl_namcfg  ! local names
300      INTEGER ::   ji, jj              ! dummy loop indexes
301      LOGICAL ::   ll_top, ll_bot      ! local logical
302      INTEGER ::   ios, inum, imk      ! local integers
303      REAL(wp)::   zmsk, zzz, zcd      ! local scalars
304      REAL(wp), DIMENSION(jpi,jpj) ::   zmsk_boost   ! 2D workspace
305      !!
306      NAMELIST/namdrg_top/ rn_Cd0, rn_Uc0, rn_Cdmax, rn_ke0, rn_z0, ln_boost, rn_boost
307      NAMELIST/namdrg_bot/ rn_Cd0, rn_Uc0, rn_Cdmax, rn_ke0, rn_z0, ln_boost, rn_boost
308      !!----------------------------------------------------------------------
309      !
310      !                          !==  set TOP / BOTTOM specificities  ==!
311      ll_top = .FALSE.
312      ll_bot = .FALSE.
313      !
314      SELECT CASE (cd_topbot)
315      CASE( 'TOP   ' )
316         ll_top = .TRUE.
317         cl_namdrg  = 'namdrg_top'
318         cl_namref  = 'namdrg_top in reference     namelist'
319         cl_namcfg  = 'namdrg_top in configuration namelist'
320         cl_file    = 'tfr_coef.nc'
321         cl_varname = 'tfr_coef'
322      CASE( 'BOTTOM' )
323         ll_bot = .TRUE.
324         cl_namdrg  = 'namdrg_bot'
325         cl_namref  = 'namdrg_bot  in reference     namelist'
326         cl_namcfg  = 'namdrg_bot  in configuration namelist'
327         cl_file    = 'bfr_coef.nc'
328         cl_varname = 'bfr_coef'
329      CASE DEFAULT
330         CALL ctl_stop( 'drg_init: bad value for cd_topbot ' )
331      END SELECT
332      !
333      !                          !==  read namlist  ==!
334      !
335      IF(ll_top)   READ  ( numnam_ref, namdrg_top, IOSTAT = ios, ERR = 901)
336      IF(ll_bot)   READ  ( numnam_ref, namdrg_bot, IOSTAT = ios, ERR = 901)
337901   IF( ios /= 0 )   CALL ctl_nam( ios , TRIM(cl_namref) )
338      IF(ll_top)   READ  ( numnam_cfg, namdrg_top, IOSTAT = ios, ERR = 902 )
339      IF(ll_bot)   READ  ( numnam_cfg, namdrg_bot, IOSTAT = ios, ERR = 902 )
340902   IF( ios >  0 )   CALL ctl_nam( ios , TRIM(cl_namcfg) )
341      IF(lwm .AND. ll_top)   WRITE ( numond, namdrg_top )
342      IF(lwm .AND. ll_bot)   WRITE ( numond, namdrg_bot )
343      !
344      IF(lwp) THEN
345         WRITE(numout,*)
346         WRITE(numout,*) '   Namelist ',TRIM(cl_namdrg),' : set ',TRIM(cd_topbot),' friction parameters'
347         WRITE(numout,*) '      drag coefficient                        rn_Cd0   = ', rn_Cd0
348         WRITE(numout,*) '      characteristic velocity (linear case)   rn_Uc0   = ', rn_Uc0, ' m/s'
349         WRITE(numout,*) '      non-linear drag maximum                 rn_Cdmax = ', rn_Cdmax
350         WRITE(numout,*) '      background kinetic energy  (n-l case)   rn_ke0   = ', rn_ke0
351         WRITE(numout,*) '      bottom roughness           (n-l case)   rn_z0    = ', rn_z0
352         WRITE(numout,*) '      set a regional boost of Cd0             ln_boost = ', ln_boost
353         WRITE(numout,*) '         associated boost factor              rn_boost = ', rn_boost
354      ENDIF
355      !
356      !                          !==  return some namelist parametres  ==!   (used in non_lin and loglayer cases)
357      pCdmin = rn_Cd0
358      pCdmax = rn_Cdmax
359      pz0    = rn_z0
360      pke0   = rn_ke0
361      !
362      !                          !==  mask * boost factor  ==!
363      !
364      IF( ln_boost ) THEN           !* regional boost:   boost factor = 1 + regional boost
365         IF(lwp) WRITE(numout,*)
366         IF(lwp) WRITE(numout,*) '   ==>>>   use a regional boost read in ', TRIM(cl_file), ' file'
367         IF(lwp) WRITE(numout,*) '           using enhancement factor of ', rn_boost
368         ! cl_varname is a coefficient in [0,1] giving where to apply the regional boost
369         CALL iom_open ( TRIM(cl_file), inum )
370         CALL iom_get  ( inum, jpdom_data, TRIM(cl_varname), zmsk_boost, 1 )
371         CALL iom_close( inum)
372         zmsk_boost(:,:) = 1._wp + rn_boost * zmsk_boost(:,:)
373         !
374      ELSE                          !* no boost:   boost factor = 1
375         zmsk_boost(:,:) = 1._wp
376      ENDIF
377      !                             !* mask outside ocean cavities area (top) or land area (bot)
378      IF(ll_top)   zmsk_boost(:,:) = zmsk_boost(:,:) * ssmask(:,:) * (1. - tmask(:,:,1) )  ! none zero in ocean cavities only
379      IF(ll_bot)   zmsk_boost(:,:) = zmsk_boost(:,:) * ssmask(:,:)                         ! x seafloor mask
380      !
381      !
382      SELECT CASE( ndrg )
383      !
384      CASE( np_OFF  )            !==  No top/bottom friction  ==!   (pCdU = 0)
385         IF(lwp) WRITE(numout,*)
386         IF(lwp) WRITE(numout,*) '   ==>>>   ',TRIM(cd_topbot),' free-slip, friction set to zero'
387         !
388         l_zdfdrg = .FALSE.         ! no time variation of the drag: set it one for all
389         !
390         pCdU(:,:) = 0._wp
391         pCd0(:,:) = 0._wp
392         !
393      CASE( np_lin )             !==  linear friction  ==!   (pCdU = Cd0 * Uc0)
394         IF(lwp) WRITE(numout,*)
395         IF(lwp) WRITE(numout,*) '   ==>>>   linear ',TRIM(cd_topbot),' friction (constant coef = Cd0*Uc0 = ', rn_Cd0*rn_Uc0, ')'
396         !
397         l_zdfdrg = .FALSE.         ! no time variation of the Cd*|U| : set it one for all
398         !                     
399         pCd0(:,:) = rn_Cd0 * zmsk_boost(:,:)  !* constant in time drag coefficient (= mask (and boost) Cd0)
400         pCdU(:,:) = - pCd0(:,:) * rn_Uc0      !  using a constant velocity
401         !
402      CASE( np_non_lin )         !== non-linear friction  ==!   (pCd0 = Cd0 )
403         IF(lwp) WRITE(numout,*)
404         IF(lwp) WRITE(numout,*) '   ==>>>   quadratic ',TRIM(cd_topbot),' friction (propotional to module of the velocity)'
405         IF(lwp) WRITE(numout,*) '   with    a drag coefficient Cd0 = ', rn_Cd0, ', and'
406         IF(lwp) WRITE(numout,*) '           a background velocity module of (rn_ke0)^1/2 = ', SQRT(rn_ke0), 'm/s)'
407         !
408         l_zdfdrg = .TRUE.          !* Cd*|U| updated at each time-step (it depends on ocean velocity)
409         !
410         pCd0(:,:) = rn_Cd0 * zmsk_boost(:,:)  !* constant in time proportionality coefficient (= mask (and boost) Cd0)
411         pCdU(:,:) = 0._wp                     ! 
412         !
413      CASE( np_loglayer )       !== logarithmic layer formulation of friction  ==!   (CdU = (vkarman log(z/z0))^2 |U| )
414         IF(lwp) WRITE(numout,*)
415         IF(lwp) WRITE(numout,*) '   ==>>>   quadratic ',TRIM(cd_topbot),' drag (propotional to module of the velocity)'
416         IF(lwp) WRITE(numout,*) '   with   a logarithmic Cd0 formulation Cd0 = ( vkarman log(z/z0) )^2 ,'
417         IF(lwp) WRITE(numout,*) '          a background velocity module of (rn_ke0)^1/2 = ', SQRT(pke0), 'm/s), '
418         IF(lwp) WRITE(numout,*) '          a logarithmic formulation: a roughness of ', pz0, ' meters,   and '
419         IF(lwp) WRITE(numout,*) '          a proportionality factor bounded by min/max values of ', pCdmin, pCdmax
420         !
421         l_zdfdrg = .TRUE.          !* Cd*|U| updated at each time-step (it depends on ocean velocity)
422         !
423         IF( ln_linssh ) THEN       !* pCd0 = (v log(z/z0))^2   as velocity points have a fixed z position
424            IF(lwp) WRITE(numout,*)
425            IF(lwp) WRITE(numout,*) '   N.B.   linear free surface case, Cd0 computed one for all'
426            !
427            l_log_not_linssh = .FALSE.    !- don't update Cd at each time step
428            !
429            DO jj = 1, jpj                   ! pCd0 = mask (and boosted) logarithmic drag coef.
430               DO ji = 1, jpi
431                  zzz =  0.5_wp * e3t_0(ji,jj,k_mk(ji,jj))
432                  zcd = (  vkarmn / LOG( zzz / rn_z0 )  )**2
433                  pCd0(ji,jj) = zmsk_boost(ji,jj) * MIN(  MAX( rn_Cd0 , zcd ) , rn_Cdmax  )  ! rn_Cd0 < Cd0 < rn_Cdmax
434               END DO
435            END DO
436         ELSE                       !* Cd updated at each time-step ==> pCd0 = mask * boost
437            IF(lwp) WRITE(numout,*)
438            IF(lwp) WRITE(numout,*) '   N.B.   non-linear free surface case, Cd0 updated at each time-step '
439            !
440            l_log_not_linssh = .TRUE.     ! compute the drag coef. at each time-step
441            !
442            pCd0(:,:) = zmsk_boost(:,:)
443         ENDIF
444         pCdU(:,:) = 0._wp          ! initialisation to zero (will be updated at each time step)
445         !
446      CASE DEFAULT
447         CALL ctl_stop( 'drg_init: bad flag value for ndrg ' )
448      END SELECT
449      !
450   END SUBROUTINE drg_init
451
452   !!======================================================================
453END MODULE zdfdrg
Note: See TracBrowser for help on using the repository browser.