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.
trabbl.F90 in branches/dev_001_GM/NEMO/OPA_SRC/TRA – NEMO

source: branches/dev_001_GM/NEMO/OPA_SRC/TRA/trabbl.F90 @ 786

Last change on this file since 786 was 786, checked in by gm, 16 years ago

dev_001_GM - merge TRC-TRA on OPA only, trabbl & zpshde not done and trdmld not OK - compilation OK

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 23.0 KB
Line 
1MODULE trabbl
2   !!==============================================================================
3   !!                       ***  MODULE  trabbl  ***
4   !! Ocean physics :  advective and/or diffusive bottom boundary layer scheme
5   !!==============================================================================
6   !! History :  8.0  !  96-06  (L. Mortier)  Original code
7   !!            8.0  !  97-11  (G. Madec)  Optimization
8   !!            8.5  !  02-08  (G. Madec)  free form + modules
9   !!----------------------------------------------------------------------
10#if   defined key_trabbl_dif   ||   defined key_trabbl_adv   || defined key_esopa
11   !!----------------------------------------------------------------------
12   !!   'key_trabbl_dif'   or            diffusive bottom boundary layer
13   !!   'key_trabbl_adv'                 advective bottom boundary layer
14   !!----------------------------------------------------------------------
15   !!----------------------------------------------------------------------
16   !!   tra_bbl_dif  : update the active tracer trends due to the bottom
17   !!                  boundary layer (diffusive only)
18   !!   tra_bbl_adv  : update the active tracer trends due to the bottom
19   !!                  boundary layer (advective and/or diffusive)
20   !!   tra_bbl_init : initialization, namlist read, parameters control
21   !!----------------------------------------------------------------------
22   USE oce                ! ocean dynamics and active tracers
23   USE dom_oce            ! ocean space and time domain
24   USE trdmod             ! ocean active tracers trends
25   USE trdmod_oce         ! ocean variables trends
26   USE in_out_manager     ! I/O manager
27   USE lbclnk             ! ocean lateral boundary conditions
28   USE prtctl             ! Print control
29
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC tra_bbl_dif    ! routine called by step.F90
34   PUBLIC tra_bbl_adv    ! routine called by step.F90
35
36   !!* Namelist nambbl: bottom boundary layer
37   REAL(wp), PUBLIC ::   atrbbl = 1.e+3   !: lateral coeff. for bottom boundary layer scheme (m2/s)
38
39# if defined key_trabbl_dif
40   LOGICAL, PUBLIC, PARAMETER ::   lk_trabbl_dif = .TRUE.          !: diffusive bottom boundary layer flag
41# else
42   LOGICAL, PUBLIC, PARAMETER ::   lk_trabbl_dif = .FALSE.         !: diffusive bottom boundary layer flag
43# endif
44
45# if defined key_trabbl_adv
46   LOGICAL, PUBLIC, PARAMETER ::   lk_trabbl_adv = .TRUE.   !: advective bottom boundary layer flag
47   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   u_bbl      !: 3 components of the velocity
48   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   v_bbl      !: associated with advective BBL
49   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   w_bbl      !: (only affect tracer)
50# else
51   LOGICAL, PUBLIC, PARAMETER ::   lk_trabbl_adv = .FALSE.  !: advective bottom boundary layer flag
52# endif
53
54   INTEGER, DIMENSION(jpi,jpj) ::   mbkt          ! vertical index of the bottom ocean T-level
55   INTEGER, DIMENSION(jpi,jpj) ::   mbku, mbkv    ! vertical index of the bottom ocean U/V-level
56
57   !! * Substitutions
58#  include "domzgr_substitute.h90"
59#  include "vectopt_loop_substitute.h90"
60   !!----------------------------------------------------------------------
61   !!   OPA 9.0 , LOCEAN-IPSL (2006)
62   !! $Header$
63   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
64   !!----------------------------------------------------------------------
65
66CONTAINS
67
68   SUBROUTINE tra_bbl_dif( kt )
69      !!----------------------------------------------------------------------
70      !!                   ***  ROUTINE tra_bbl_dif  ***
71      !!
72      !! ** Purpose :   Compute the before tracer (t & s) trend associated
73      !!      with the bottom boundary layer and add it to the general trend
74      !!      of tracer equations. The bottom boundary layer is supposed to be
75      !!      a purely diffusive bottom boundary layer.
76      !!
77      !! ** Method  :   When the product grad( rho) * grad(h) < 0 (where grad
78      !!      is an along bottom slope gradient) an additional lateral diffu-
79      !!      sive trend along the bottom slope is added to the general tracer
80      !!      trend, otherwise nothing is done.
81      !!      Second order operator (laplacian type) with variable coefficient
82      !!      computed as follow for temperature (idem on s):
83      !!         difft = 1/(e1t*e2t*e3t) { di-1[ ahbt e2u*e3u/e1u di[ztb] ]
84      !!                                 + dj-1[ ahbt e1v*e3v/e2v dj[ztb] ] }
85      !!      where ztb is a 2D array: the bottom ocean temperature and ahtb
86      !!      is a time and space varying diffusive coefficient defined by:
87      !!         ahbt = zahbp    if grad(rho).grad(h) < 0
88      !!              = 0.       otherwise.
89      !!      Note that grad(.) is the along bottom slope gradient. grad(rho)
90      !!      is evaluated using the local density (i.e. referenced at the
91      !!      local depth). Typical value of ahbt is 2000 m2/s (equivalent to
92      !!      a downslope velocity of 20 cm/s if the condition for slope
93      !!      convection is satified)
94      !!      Add this before trend to the general trend (ta,sa) of the
95      !!      botton ocean tracer point:
96      !!         ta = ta + difft
97      !!
98      !! ** Action  : - update (ta,sa) at the bottom level with the bottom
99      !!                boundary layer trend
100      !!              - save the trends in ztrdt/ztrds ('key_trdtra')
101      !!
102      !! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591.
103      !!----------------------------------------------------------------------
104      USE oce, ONLY :   ztrdt => ua   ! use ua as 3D workspace   
105      USE oce, ONLY :   ztrds => va   ! use va as 3D workspace   
106      USE eosbn2                      ! equation of state
107      !!
108      INTEGER, INTENT( in ) ::   kt   ! ocean time-step
109      !!
110      INTEGER  ::   ji, jj                   ! dummy loop indices
111      INTEGER  ::   ik
112      INTEGER  ::   ii0, ii1, ij0, ij1       ! temporary integers
113      INTEGER  ::   iku1, iku2, ikv1,ikv2   ! temporary intergers
114      REAL(wp) ::   ze3u, ze3v              ! temporary scalars
115      INTEGER  ::   iku, ikv
116      REAL(wp) ::   zsign, zt, zs, zh, zalbet   ! temporary scalars
117      REAL(wp) ::   zgdrho, zbtr, zta, zsa
118      REAL(wp), DIMENSION(jpi,jpj) ::   zki, zkj, zkw, zkx, zky, zkz   ! 2D workspace arrays
119      REAL(wp), DIMENSION(jpi,jpj) ::   ztnb, zsnb, zdep
120      REAL(wp), DIMENSION(jpi,jpj) ::   ztbb, zsbb, zahu, zahv
121      REAL(wp) ::    fsalbt, pft, pfs, pfh   ! statement function
122      !!----------------------------------------------------------------------
123      ! ratio alpha/beta
124      ! ================
125      !  fsalbt: ratio of thermal over saline expension coefficients
126      !       pft :  potential temperature in degrees celcius
127      !       pfs :  salinity anomaly (s-35) in psu
128      !       pfh :  depth in meters
129
130      fsalbt( pft, pfs, pfh ) =                                              &
131         ( ( ( -0.255019e-07 * pft + 0.298357e-05 ) * pft                    &
132                                   - 0.203814e-03 ) * pft                    &
133                                   + 0.170907e-01 ) * pft                    &
134                                   + 0.665157e-01                            &
135         +(-0.678662e-05 * pfs - 0.846960e-04 * pft + 0.378110e-02 ) * pfs   &
136         +  ( ( - 0.302285e-13 * pfh                                         &
137                - 0.251520e-11 * pfs                                         &
138                + 0.512857e-12 * pft * pft          ) * pfh                  &
139                                     - 0.164759e-06   * pfs                  &
140             +(   0.791325e-08 * pft - 0.933746e-06 ) * pft                  &
141                                     + 0.380374e-04 ) * pfh   
142      !!----------------------------------------------------------------------
143
144      IF( kt == nit000 )   CALL tra_bbl_init
145
146      IF( l_trdtra )   THEN         ! Save ta and sa trends
147         ztrdt(:,:,:) = ta(:,:,:) 
148         ztrds(:,:,:) = sa(:,:,:) 
149      ENDIF
150
151      ! 0. 2D fields of bottom temperature and salinity, and bottom slope
152      ! -----------------------------------------------------------------
153      ! mbathy= number of w-level, minimum value=1 (cf dommsk.F)
154#  if defined key_vectopt_loop   &&   ! defined key_mpp_omp
155      jj = 1
156      DO ji = 1, jpij   ! vector opt. (forced unrolling)
157#  else
158      DO jj = 1, jpj
159         DO ji = 1, jpi
160#  endif
161            ik = mbkt(ji,jj)                              ! index of the bottom ocean T-level
162            ztnb(ji,jj) = tn(ji,jj,ik) * tmask(ji,jj,1)   ! masked now T and S at ocean bottom
163            zsnb(ji,jj) = sn(ji,jj,ik) * tmask(ji,jj,1)
164            ztbb(ji,jj) = tb(ji,jj,ik) * tmask(ji,jj,1)   ! masked before T and S at ocean bottom
165            zsbb(ji,jj) = sb(ji,jj,ik) * tmask(ji,jj,1)
166            zdep(ji,jj) = fsdept(ji,jj,ik)                ! depth of the ocean bottom T-level
167#  if ! defined key_vectopt_loop   ||   defined key_mpp_omp
168         END DO
169#  endif
170      END DO
171
172      IF( ln_zps ) THEN      ! partial steps correction
173# if defined key_vectopt_loop   &&   ! defined key_mpp_omp
174         jj = 1
175         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
176# else
177         DO jj = 1, jpjm1
178            DO ji = 1, jpim1
179# endif
180               iku1 = MAX( mbathy(ji+1,jj  )-1, 1 )
181               iku2 = MAX( mbathy(ji  ,jj  )-1, 1 )
182               ikv1 = MAX( mbathy(ji  ,jj+1)-1, 1 )
183               ikv2 = MAX( mbathy(ji  ,jj  )-1, 1 )
184               ze3u = MIN( fse3u(ji,jj,iku1), fse3u(ji,jj,iku2) ) 
185               ze3v = MIN( fse3v(ji,jj,ikv1), fse3v(ji,jj,ikv2) ) 
186               zahu(ji,jj) = atrbbl * e2u(ji,jj) * ze3u / e1u(ji,jj) * umask(ji,jj,1)
187               zahv(ji,jj) = atrbbl * e1v(ji,jj) * ze3v / e2v(ji,jj) * vmask(ji,jj,1)
188# if ! defined key_vectopt_loop   ||   defined key_mpp_omp
189            END DO
190# endif
191         END DO
192      ELSE                    ! z-coordinate - full steps or s-coordinate
193#   if defined key_vectopt_loop   &&   ! defined key_mpp_omp
194         jj = 1
195         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
196#   else
197         DO jj = 1, jpjm1
198            DO ji = 1, jpim1
199#   endif
200               iku = mbku(ji,jj)
201               ikv = mbkv(ji,jj)
202               zahu(ji,jj) = atrbbl * e2u(ji,jj) * fse3u(ji,jj,iku) / e1u(ji,jj) * umask(ji,jj,1)
203               zahv(ji,jj) = atrbbl * e1v(ji,jj) * fse3v(ji,jj,ikv) / e2v(ji,jj) * vmask(ji,jj,1)
204#   if ! defined key_vectopt_loop   ||   defined key_mpp_omp
205            END DO
206#   endif
207         END DO
208      ENDIF
209
210      ! 1. Criteria of additional bottom diffusivity: grad(rho).grad(h)<0
211      ! --------------------------------------------
212      ! Sign of the local density gradient along the i- and j-slopes
213      ! multiplied by the slope of the ocean bottom
214
215      SELECT CASE ( neos )
216
217      CASE ( 0 )                 ! Jackett and McDougall (1994) formulation
218
219#  if defined key_vectopt_loop   &&   ! defined key_mpp_omp
220      jj = 1
221      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
222#  else
223      DO jj = 1, jpjm1
224         DO ji = 1, jpim1
225#  endif
226            ! temperature, salinity anomalie and depth
227            zt = 0.5 * ( ztnb(ji,jj) + ztnb(ji+1,jj) )
228            zs = 0.5 * ( zsnb(ji,jj) + zsnb(ji+1,jj) ) - 35.0
229            zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) )
230            ! masked ratio alpha/beta
231            zalbet = fsalbt( zt, zs, zh )*umask(ji,jj,1)
232            ! local density gradient along i-bathymetric slope
233            zgdrho = zalbet * ( ztnb(ji+1,jj) - ztnb(ji,jj) )   &
234                   -          ( zsnb(ji+1,jj) - zsnb(ji,jj) )
235            ! sign of local i-gradient of density multiplied by the i-slope
236            zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) )
237            zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj)
238#  if ! defined key_vectopt_loop   ||   defined key_mpp_omp
239         END DO
240#  endif
241      END DO
242
243#  if defined key_vectopt_loop   &&   ! defined key_mpp_omp
244      jj = 1
245      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
246#  else
247      DO jj = 1, jpjm1
248         DO ji = 1, jpim1
249#  endif
250            ! temperature, salinity anomalie and depth
251            zt = 0.5 * ( ztnb(ji,jj+1) + ztnb(ji,jj) )
252            zs = 0.5 * ( zsnb(ji,jj+1) + zsnb(ji,jj) ) - 35.0
253            zh = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) )
254            ! masked ratio alpha/beta
255            zalbet = fsalbt( zt, zs, zh )*vmask(ji,jj,1)
256            ! local density gradient along j-bathymetric slope
257            zgdrho = zalbet * ( ztnb(ji,jj+1) - ztnb(ji,jj) )   &
258               &   -          ( zsnb(ji,jj+1) - zsnb(ji,jj) )
259            ! sign of local j-gradient of density multiplied by the j-slope
260            zsign = sign( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) )
261            zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj)
262#  if ! defined key_vectopt_loop   ||   defined key_mpp_omp
263         END DO
264#  endif
265      END DO
266
267      CASE ( 1 )               ! Linear formulation function of temperature only
268                               !
269#  if defined key_vectopt_loop   &&   ! defined key_mpp_omp
270      jj = 1
271      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
272#  else
273      DO jj = 1, jpjm1
274         DO ji = 1, jpim1
275#  endif
276            ! local 'density/temperature' gradient along i-bathymetric slope
277            zgdrho =  ztnb(ji+1,jj) - ztnb(ji,jj) 
278            ! sign of local i-gradient of density multiplied by the i-slope
279            zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) )
280            zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj)
281#  if ! defined key_vectopt_loop   ||   defined key_mpp_omp
282         END DO
283#  endif
284      END DO
285
286#  if defined key_vectopt_loop   &&   ! defined key_mpp_omp
287      jj = 1
288      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
289#  else
290      DO jj = 1, jpjm1
291         DO ji = 1, jpim1
292#  endif
293            ! local density gradient along j-bathymetric slope
294            zgdrho =  ztnb(ji,jj+1) - ztnb(ji,jj) 
295            ! sign of local j-gradient of density multiplied by the j-slope
296            zsign = sign( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) )
297            zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj)
298#  if ! defined key_vectopt_loop   ||   defined key_mpp_omp
299         END DO
300#  endif
301      END DO
302
303      CASE ( 2 )               ! Linear formulation function of temperature and salinity
304
305#  if defined key_vectopt_loop   &&   ! defined key_mpp_omp
306      jj = 1
307      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
308#  else
309      DO jj = 1, jpjm1
310         DO ji = 1, jpim1
311#  endif     
312            ! local density gradient along i-bathymetric slope
313            zgdrho = - ( rbeta*( zsnb(ji+1,jj) - zsnb(ji,jj) )   &
314                     -  ralpha*( ztnb(ji+1,jj) - ztnb(ji,jj) ) )
315            ! sign of local i-gradient of density multiplied by the i-slope
316            zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) )
317            zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj)
318#  if ! defined key_vectopt_loop   ||   defined key_mpp_omp
319         END DO
320#  endif
321      END DO
322
323#  if defined key_vectopt_loop   &&   ! defined key_mpp_omp
324      jj = 1
325      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
326#  else
327      DO jj = 1, jpjm1
328         DO ji = 1, jpim1
329#  endif     
330            ! local density gradient along j-bathymetric slope
331            zgdrho = - ( rbeta*( zsnb(ji,jj+1) - zsnb(ji,jj) )   &
332                     -  ralpha*( ztnb(ji,jj+1) - ztnb(ji,jj) ) )   
333            ! sign of local j-gradient of density multiplied by the j-slope
334            zsign = sign( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) )
335            zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj)
336#  if ! defined key_vectopt_loop   ||   defined key_mpp_omp
337         END DO
338#  endif
339      END DO
340     
341      CASE DEFAULT
342
343         WRITE(ctmp1,*) '          bad flag value for neos = ', neos
344         CALL ctl_stop(ctmp1)
345
346      END SELECT
347
348      ! 2. Additional second order diffusive trends
349      ! -------------------------------------------
350
351      ! first derivative (gradient)
352#  if defined key_vectopt_loop   &&   ! defined key_mpp_omp
353      jj = 1
354      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
355#  else
356      DO jj = 1, jpjm1
357         DO ji = 1, jpim1
358#  endif
359            zkx(ji,jj) = zki(ji,jj) * ( ztbb(ji+1,jj) - ztbb(ji,jj) )
360            zkz(ji,jj) = zki(ji,jj) * ( zsbb(ji+1,jj) - zsbb(ji,jj) )
361
362            zky(ji,jj) = zkj(ji,jj) * ( ztbb(ji,jj+1) - ztbb(ji,jj) )
363            zkw(ji,jj) = zkj(ji,jj) * ( zsbb(ji,jj+1) - zsbb(ji,jj) )
364#  if ! defined key_vectopt_loop   ||   defined key_mpp_omp
365         END DO
366#  endif
367      END DO
368
369      IF( cp_cfg == "orca" ) THEN
370         !
371         SELECT CASE ( jp_cfg )
372         !                                           ! =======================
373         CASE ( 2 )                                  !  ORCA_R2 configuration
374            !                                        ! =======================
375            ! Gibraltar enhancement of BBL
376            ij0 = 102   ;   ij1 = 102
377            ii0 = 139   ;   ii1 = 140 
378            zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) )
379            zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) )
380            !
381            ! Red Sea enhancement of BBL
382            ij0 =  88   ;   ij1 =  88
383            ii0 = 161   ;   ii1 = 162
384            zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 10.e0 * zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) )
385            zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 10.e0 * zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) )
386            !
387            !                                        ! =======================
388         CASE ( 4 )                                  !  ORCA_R4 configuration
389            !                                        ! =======================
390            ! Gibraltar enhancement of BBL
391            ij0 =  52   ;   ij1 =  52
392            ii0 =  70   ;   ii1 =  71 
393            zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) )
394            zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) )
395            !
396         END SELECT
397      !
398      ENDIF
399
400
401      ! second derivative (divergence) and add to the general tracer trend
402#  if defined key_vectopt_loop   &&   ! defined key_mpp_omp
403      jj = 1
404      DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling)
405#  else
406      DO jj = 2, jpjm1
407         DO ji = 2, jpim1
408#  endif
409            ik = max( mbathy(ji,jj)-1, 1 )
410            zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,ik) )
411            zta = (  zkx(ji,jj) - zkx(ji-1,jj  )    &
412               &   + zky(ji,jj) - zky(ji  ,jj-1)  ) * zbtr
413            zsa = (  zkz(ji,jj) - zkz(ji-1,jj  )    &
414               &   + zkw(ji,jj) - zkw(ji  ,jj-1)  ) * zbtr
415            ta(ji,jj,ik) = ta(ji,jj,ik) + zta
416            sa(ji,jj,ik) = sa(ji,jj,ik) + zsa
417#  if ! defined key_vectopt_loop   ||   defined key_mpp_omp
418         END DO
419#  endif
420      END DO
421
422      IF( l_trdtra ) THEN      ! save the BBL lateral diffusion trends for diagnostic
423         ztrdt(:,:,:) = ta(:,:,:) - ztrdt(:,:,:)
424         ztrds(:,:,:) = sa(:,:,:) - ztrds(:,:,:)
425         CALL trd_tra( kt, jp_tem, jpt_trd_bbl, 'TRA', ptrd3d=ztrdt)
426         CALL trd_tra( kt, jp_sal, jpt_trd_bbl, 'TRA', ptrd3d=ztrds)
427      ENDIF
428
429      IF(ln_ctl)   CALL prt_ctl( tab3d_1=ta, clinfo1=' bbl  - Ta: ', mask1=tmask,   &
430         &                       tab3d_2=sa, clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
431      !
432   END SUBROUTINE tra_bbl_dif
433
434# if defined key_trabbl_adv
435   !!----------------------------------------------------------------------
436   !!   'key_trabbl_adv'                    advective bottom boundary layer
437   !!----------------------------------------------------------------------
438#  include "trabbl_adv.h90"
439# else
440   !!----------------------------------------------------------------------
441   !!   Default option :                 NO advective bottom boundary layer
442   !!----------------------------------------------------------------------
443   SUBROUTINE tra_bbl_adv (kt )              ! Empty routine
444      INTEGER, INTENT(in) :: kt
445      WRITE(*,*) 'tra_bbl_adv: You should not have seen this print! error?', kt
446   END SUBROUTINE tra_bbl_adv
447# endif
448
449   SUBROUTINE tra_bbl_init
450      !!----------------------------------------------------------------------
451      !!                  ***  ROUTINE tra_bbl_init  ***
452      !!
453      !! ** Purpose :   Initialization for the bottom boundary layer scheme.
454      !!
455      !! ** Method  :   Read the nambbl namelist and check the parameters
456      !!      called by tra_bbl at the first timestep (nit000)
457      !!----------------------------------------------------------------------
458      INTEGER ::   ji, jj      ! dummy loop indices
459      REAL(wp),  DIMENSION(jpi,jpj) :: zmbk 
460
461      NAMELIST/nambbl/ atrbbl
462      !!----------------------------------------------------------------------
463
464      REWIND ( numnam )              ! Read Namelist nambbl : bottom boundary layer scheme
465      READ   ( numnam, nambbl )
466
467      IF(lwp) THEN                   ! Parameter control and print
468         WRITE(numout,*)
469         WRITE(numout,*) 'tra_bbl_init : '
470         WRITE(numout,*) '~~~~~~~~~~~~'
471         IF (lk_trabbl_dif )   WRITE(numout,*) '               * Diffusive Bottom Boundary Layer'
472         IF( lk_trabbl_adv )   WRITE(numout,*) '               * Advective Bottom Boundary Layer'
473         WRITE(numout,*) '       Namelist nambbl : set bbl parameters'
474         WRITE(numout,*) '          bottom boundary layer coef.    atrbbl = ', atrbbl
475      ENDIF
476 
477      DO jj = 1, jpj
478         DO ji = 1, jpi
479            mbkt(ji,jj) = MAX( mbathy(ji,jj) - 1, 1 )   ! vertical index of the bottom ocean T-level
480         END DO
481      END DO
482      DO jj = 1, jpjm1
483         DO ji = 1, jpim1
484            mbku(ji,jj) = MAX( MIN( mbathy(ji+1,jj  ), mbathy(ji,jj) ) - 1, 1 )
485            mbkv(ji,jj) = MAX( MIN( mbathy(ji  ,jj+1), mbathy(ji,jj) ) - 1, 1 )
486         END DO
487      END DO
488
489      zmbk(:,:) = FLOAT( mbku (:,:) )   
490      CALL lbc_lnk(zmbk,'U',1.)
491      mbku(:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
492
493      zmbk(:,:) = FLOAT( mbkv (:,:) )   
494      CALL lbc_lnk(zmbk,'V',1.)
495      mbkv(:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
496
497# if defined key_trabbl_adv
498      w_bbl(:,:,:) = 0.e0          ! initialisation of w_bbl to zero
499# endif
500      !
501   END SUBROUTINE tra_bbl_init
502
503#else
504   !!----------------------------------------------------------------------
505   !!   Dummy module :                      No bottom boundary layer scheme
506   !!----------------------------------------------------------------------
507   LOGICAL, PUBLIC, PARAMETER ::   lk_trabbl_dif = .FALSE.   !: diff bbl flag
508   LOGICAL, PUBLIC, PARAMETER ::   lk_trabbl_adv = .FALSE.   !: adv  bbl flag
509CONTAINS
510   SUBROUTINE tra_bbl_dif( kt )              ! Empty routine
511      WRITE(*,*) 'tra_bbl_dif: You should not have seen this print! error?', kt
512   END SUBROUTINE tra_bbl_dif
513   SUBROUTINE tra_bbl_adv( kt )              ! Empty routine
514      WRITE(*,*) 'tra_bbl_adv: You should not have seen this print! error?', kt
515   END SUBROUTINE tra_bbl_adv
516#endif
517
518   !!======================================================================
519END MODULE trabbl
Note: See TracBrowser for help on using the repository browser.