Changeset 2528 for trunk/NEMOGCM/NEMO/OPA_SRC/TRA/trabbl.F90
- Timestamp:
- 2010-12-27T18:33:53+01:00 (13 years ago)
- File:
-
- 1 edited
-
trunk/NEMOGCM/NEMO/OPA_SRC/TRA/trabbl.F90 (modified) (3 diffs, 1 prop)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/TRA/trabbl.F90
- Property svn:eol-style deleted
r1601 r2528 4 4 !! Ocean physics : advective and/or diffusive bottom boundary layer scheme 5 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 6 !! History : OPA ! 1996-06 (L. Mortier) Original code 7 !! 8.0 ! 1997-11 (G. Madec) Optimization 8 !! NEMO 1.0 ! 2002-08 (G. Madec) free form + modules 9 !! - ! 2004-01 (A. de Miranda, G. Madec, J.M. Molines ) add advective bbl 10 !! 3.3 ! 2009-11 (G. Madec) merge trabbl and trabbl_adv + style + optimization 11 !! - ! 2010-04 (G. Madec) Campin & Goosse advective bbl 12 !! - ! 2010-06 (C. Ethe, G. Madec) merge TRA-TRC 13 !! - ! 2010-11 (G. Madec) add mbk. arrays associated to the deepest ocean level 14 !!---------------------------------------------------------------------- 15 #if defined key_trabbl || defined key_esopa 16 !!---------------------------------------------------------------------- 17 !! 'key_trabbl' or bottom boundary layer 18 !!---------------------------------------------------------------------- 19 !! tra_bbl : update the tracer trends due to the bottom boundary layer (advective and/or diffusive) 20 !! tra_bbl_dif : generic routine to compute bbl diffusive trend 21 !! tra_bbl_adv : generic routine to compute bbl advective trend 22 !! bbl : computation of bbl diffu. flux coef. & transport in bottom boundary layer 23 !! tra_bbl_init : initialization, namelist read, parameters control 24 !!---------------------------------------------------------------------- 25 USE oce ! ocean dynamics and active tracers 26 USE dom_oce ! ocean space and time domain 27 USE phycst ! physical constant 28 USE eosbn2 ! equation of state 29 USE trdmod_oce ! trends: ocean variables 30 USE trdtra ! trends: active tracers 31 USE iom ! IOM server 32 USE in_out_manager ! I/O manager 33 USE lbclnk ! ocean lateral boundary conditions 34 USE prtctl ! Print control 29 35 30 36 IMPLICIT NONE 31 37 PRIVATE 32 38 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 :: rn_ahtbbl = 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 39 PUBLIC tra_bbl ! routine called by step.F90 40 PUBLIC tra_bbl_init ! routine called by opa.F90 41 PUBLIC tra_bbl_dif ! routine called by trcbbl.F90 42 PUBLIC tra_bbl_adv ! - - - - 43 PUBLIC bbl ! routine called by trcbbl.F90 and dtadyn.F90 44 45 LOGICAL, PUBLIC, PARAMETER :: lk_trabbl = .TRUE. !: bottom boundary layer flag 46 47 ! !!* Namelist nambbl * 48 INTEGER , PUBLIC :: nn_bbl_ldf = 0 !: =1 : diffusive bbl or not (=0) 49 INTEGER , PUBLIC :: nn_bbl_adv = 0 !: =1/2 : advective bbl or not (=0) 50 ! ! =1 : advective bbl using the bottom ocean velocity 51 ! ! =2 : - - using utr_bbl proportional to grad(rho) 52 REAL(wp), PUBLIC :: rn_ahtbbl = 1.e3_wp !: along slope bbl diffusive coefficient [m2/s] 53 REAL(wp), PUBLIC :: rn_gambbl = 10.0_wp !: lateral coeff. for bottom boundary layer scheme [s] 54 55 REAL(wp), DIMENSION(jpi,jpj), PUBLIC :: utr_bbl , vtr_bbl ! u- (v-) transport in the bottom boundary layer 56 REAL(wp), DIMENSION(jpi,jpj), PUBLIC :: ahu_bbl , ahv_bbl ! masked diffusive bbl coefficients at u and v-points 57 58 INTEGER , DIMENSION(jpi,jpj) :: mbku_d , mbkv_d ! vertical index of the "lower" bottom ocean U/V-level 59 INTEGER , DIMENSION(jpi,jpj) :: mgrhu , mgrhv ! = +/-1, sign of grad(H) in u-(v-)direction 60 REAL(wp), DIMENSION(jpi,jpj) :: ahu_bbl_0, ahv_bbl_0 ! diffusive bbl flux coefficients at u and v-points 61 REAL(wp), DIMENSION(jpi,jpj) :: e3u_bbl_0, e3v_bbl_0 ! thichness of the bbl (e3) at u and v-points 62 REAL(wp), DIMENSION(jpi,jpj) :: e1e2t_r ! thichness of the bbl (e3) at u and v-points 63 LOGICAL, PUBLIC :: l_bbl !: flag to compute bbl diffu. flux coef and transport 56 64 57 65 !! * Substitutions … … 59 67 # include "vectopt_loop_substitute.h90" 60 68 !!---------------------------------------------------------------------- 61 !! OPA 9.0 , LOCEAN-IPSL (2006) 62 !! $Id$ 63 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 64 !!---------------------------------------------------------------------- 65 69 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 70 !! $Id$ 71 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 72 !!---------------------------------------------------------------------- 66 73 CONTAINS 67 74 68 SUBROUTINE tra_bbl _dif( kt )69 !!---------------------------------------------------------------------- 70 !! *** ROUTINE tra_bbl_dif***71 !! 75 SUBROUTINE tra_bbl( kt ) 76 !!---------------------------------------------------------------------- 77 !! *** ROUTINE bbl *** 78 !! 72 79 !! ** 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 80 !! with the bottom boundary layer and add it to the general 81 !! trend of tracer equations. 82 !! 83 !! ** Method : Depending on namtra_bbl namelist parameters the bbl 84 !! diffusive and/or advective contribution to the tracer trend 85 !! is added to the general tracer trend 86 !!---------------------------------------------------------------------- 87 INTEGER, INTENT( in ) :: kt ! ocean time-step 88 !! 89 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ztrdt, ztrds 90 !!---------------------------------------------------------------------- 91 92 IF( l_trdtra ) THEN !* Save ta and sa trends 93 ALLOCATE( ztrdt(jpi,jpj,jpk) ) ; ztrdt(:,:,:) = tsa(:,:,:,jp_tem) 94 ALLOCATE( ztrds(jpi,jpj,jpk) ) ; ztrds(:,:,:) = tsa(:,:,:,jp_sal) 95 ENDIF 96 97 IF( l_bbl ) CALL bbl( kt, 'TRA' ) !* bbl coef. and transport (only if not already done in trcbbl) 98 99 100 IF( nn_bbl_ldf == 1 ) THEN !* Diffusive bbl 101 CALL tra_bbl_dif( tsb, tsa, jpts ) 102 IF( ln_ctl ) & 103 CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' bbl_ldf - Ta: ', mask1=tmask, & 104 & tab3d_2=tsa(:,:,:,jp_sal), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 105 ! lateral boundary conditions ; just need for outputs 106 CALL lbc_lnk( ahu_bbl, 'U', 1. ) ; CALL lbc_lnk( ahv_bbl, 'V', 1. ) 107 CALL iom_put( "ahu_bbl", ahu_bbl ) ! bbl diffusive flux i-coef 108 CALL iom_put( "ahv_bbl", ahv_bbl ) ! bbl diffusive flux j-coef 109 END IF 110 111 IF( nn_bbl_adv /= 0 ) THEN !* Advective bbl 112 CALL tra_bbl_adv( tsb, tsa, jpts ) 113 IF(ln_ctl) & 114 CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' bbl_adv - Ta: ', mask1=tmask, & 115 & tab3d_2=tsa(:,:,:,jp_sal), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 116 ! lateral boundary conditions ; just need for outputs 117 CALL lbc_lnk( utr_bbl, 'U', 1. ) ; CALL lbc_lnk( vtr_bbl, 'V', 1. ) 118 CALL iom_put( "uoce_bbl", utr_bbl ) ! bbl i-transport 119 CALL iom_put( "voce_bbl", vtr_bbl ) ! bbl j-transport 120 END IF 121 122 IF( l_trdtra ) THEN ! save the horizontal diffusive trends for further diagnostics 123 ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:) 124 ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:) 125 CALL trd_tra( kt, 'TRA', jp_tem, jptra_trd_bbl, ztrdt ) 126 CALL trd_tra( kt, 'TRA', jp_sal, jptra_trd_bbl, ztrds ) 127 DEALLOCATE( ztrdt ) ; DEALLOCATE( ztrds ) 128 ENDIF 129 ! 130 END SUBROUTINE tra_bbl 131 132 133 SUBROUTINE tra_bbl_dif( ptb, pta, kjpt ) 134 !!---------------------------------------------------------------------- 135 !! *** ROUTINE tra_bbl_dif *** 136 !! 137 !! ** Purpose : Computes the bottom boundary horizontal and vertical 138 !! advection terms. 139 !! 140 !! ** Method : 141 !! * diffusive bbl (nn_bbl_ldf=1) : 142 !! When the product grad( rho) * grad(h) < 0 (where grad is an 143 !! along bottom slope gradient) an additional lateral 2nd order 144 !! diffusion along the bottom slope is added to the general 145 !! tracer trend, otherwise the additional trend is set to 0. 146 !! A typical value of ahbt is 2000 m2/s (equivalent to 92 147 !! a downslope velocity of 20 cm/s if the condition for slope 93 148 !! 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') 149 !! 150 !! ** Action : pta increased by the bbl diffusive trend 101 151 !! 102 152 !! 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 119 REAL(wp), DIMENSION(jpi,jpj) :: ztnb, zsnb, zdep, ztbb, zsbb, zahu, zahv 120 !! 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 ) = & 153 !! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430. 154 !!---------------------------------------------------------------------- 155 INTEGER , INTENT(in ) :: kjpt ! number of tracers 156 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields 157 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 158 !! 159 INTEGER :: ji, jj, jn ! dummy loop indices 160 INTEGER :: ik ! local integers 161 REAL(wp) :: zbtr ! local scalars 162 REAL(wp), DIMENSION(jpi,jpj) :: zptb ! tracer trend 163 !!---------------------------------------------------------------------- 164 ! 165 DO jn = 1, kjpt ! tracer loop 166 ! ! =========== 167 # if defined key_vectopt_loop 168 DO jj = 1, 1 ! vector opt. (forced unrolling) 169 DO ji = 1, jpij 170 #else 171 DO jj = 1, jpj 172 DO ji = 1, jpi 173 #endif 174 ik = mbkt(ji,jj) ! bottom T-level index 175 zptb(ji,jj) = ptb(ji,jj,ik,jn) ! bottom before T and S 176 END DO 177 END DO 178 ! ! Compute the trend 179 # if defined key_vectopt_loop 180 DO jj = 1, 1 ! vector opt. (forced unrolling) 181 DO ji = jpi+1, jpij-jpi-1 182 # else 183 DO jj = 2, jpjm1 184 DO ji = 2, jpim1 185 # endif 186 ik = mbkt(ji,jj) ! bottom T-level index 187 zbtr = e1e2t_r(ji,jj) / fse3t(ji,jj,ik) 188 pta(ji,jj,ik,jn) = pta(ji,jj,ik,jn) & 189 & + ( ahu_bbl(ji ,jj ) * ( zptb(ji+1,jj ) - zptb(ji ,jj ) ) & 190 & - ahu_bbl(ji-1,jj ) * ( zptb(ji ,jj ) - zptb(ji-1,jj ) ) & 191 & + ahv_bbl(ji ,jj ) * ( zptb(ji ,jj+1) - zptb(ji ,jj ) ) & 192 & - ahv_bbl(ji ,jj-1) * ( zptb(ji ,jj ) - zptb(ji ,jj-1) ) ) * zbtr 193 END DO 194 END DO 195 ! ! =========== 196 END DO ! end tracer 197 ! ! =========== 198 END SUBROUTINE tra_bbl_dif 199 200 201 SUBROUTINE tra_bbl_adv( ptb, pta, kjpt ) 202 !!---------------------------------------------------------------------- 203 !! *** ROUTINE trc_bbl *** 204 !! 205 !! ** Purpose : Compute the before passive tracer trend associated 206 !! with the bottom boundary layer and add it to the general trend 207 !! of tracer equations. 208 !! ** Method : advective bbl (nn_bbl_adv = 1 or 2) : 209 !! nn_bbl_adv = 1 use of the ocean near bottom velocity as bbl velocity 210 !! nn_bbl_adv = 2 follow Campin and Goosse (1999) implentation i.e. 211 !! transport proportional to the along-slope density gradient 212 !! 213 !! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591. 214 !! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430. 215 !!---------------------------------------------------------------------- 216 INTEGER , INTENT(in ) :: kjpt ! number of tracers 217 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields 218 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 219 !! 220 INTEGER :: ji, jj, jk, jn ! dummy loop indices 221 INTEGER :: iis , iid , ijs , ijd ! local integers 222 INTEGER :: ikus, ikud, ikvs, ikvd ! - - 223 REAL(wp) :: zbtr, ztra ! local scalars 224 REAL(wp) :: zu_bbl, zv_bbl ! - - 225 !!---------------------------------------------------------------------- 226 ! 227 ! ! =========== 228 DO jn = 1, kjpt ! tracer loop 229 ! ! =========== 230 # if defined key_vectopt_loop 231 DO jj = 1, 1 232 DO ji = 1, jpij-jpi-1 ! vector opt. (forced unrolling) 233 # else 234 DO jj = 1, jpjm1 235 DO ji = 1, jpim1 ! CAUTION start from i=1 to update i=2 when cyclic east-west 236 # endif 237 IF( utr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero i-direction bbl advection 238 ! down-slope i/k-indices (deep) & up-slope i/k indices (shelf) 239 iid = ji + MAX( 0, mgrhu(ji,jj) ) ; iis = ji + 1 - MAX( 0, mgrhu(ji,jj) ) 240 ikud = mbku_d(ji,jj) ; ikus = mbku(ji,jj) 241 zu_bbl = ABS( utr_bbl(ji,jj) ) 242 ! 243 ! ! up -slope T-point (shelf bottom point) 244 zbtr = e1e2t_r(iis,jj) / fse3t(iis,jj,ikus) 245 ztra = zu_bbl * ( ptb(iid,jj,ikus,jn) - ptb(iis,jj,ikus,jn) ) * zbtr 246 pta(iis,jj,ikus,jn) = pta(iis,jj,ikus,jn) + ztra 247 ! 248 DO jk = ikus, ikud-1 ! down-slope upper to down T-point (deep column) 249 zbtr = e1e2t_r(iid,jj) / fse3t(iid,jj,jk) 250 ztra = zu_bbl * ( ptb(iid,jj,jk+1,jn) - ptb(iid,jj,jk,jn) ) * zbtr 251 pta(iid,jj,jk,jn) = pta(iid,jj,jk,jn) + ztra 252 END DO 253 ! 254 zbtr = e1e2t_r(iid,jj) / fse3t(iid,jj,ikud) 255 ztra = zu_bbl * ( ptb(iis,jj,ikus,jn) - ptb(iid,jj,ikud,jn) ) * zbtr 256 pta(iid,jj,ikud,jn) = pta(iid,jj,ikud,jn) + ztra 257 ENDIF 258 ! 259 IF( vtr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero j-direction bbl advection 260 ! down-slope j/k-indices (deep) & up-slope j/k indices (shelf) 261 ijd = jj + MAX( 0, mgrhv(ji,jj) ) ; ijs = jj + 1 - MAX( 0, mgrhv(ji,jj) ) 262 ikvd = mbkv_d(ji,jj) ; ikvs = mbkv(ji,jj) 263 zv_bbl = ABS( vtr_bbl(ji,jj) ) 264 ! 265 ! up -slope T-point (shelf bottom point) 266 zbtr = e1e2t_r(ji,ijs) / fse3t(ji,ijs,ikvs) 267 ztra = zv_bbl * ( ptb(ji,ijd,ikvs,jn) - ptb(ji,ijs,ikvs,jn) ) * zbtr 268 pta(ji,ijs,ikvs,jn) = pta(ji,ijs,ikvs,jn) + ztra 269 ! 270 DO jk = ikvs, ikvd-1 ! down-slope upper to down T-point (deep column) 271 zbtr = e1e2t_r(ji,ijd) / fse3t(ji,ijd,jk) 272 ztra = zv_bbl * ( ptb(ji,ijd,jk+1,jn) - ptb(ji,ijd,jk,jn) ) * zbtr 273 pta(ji,ijd,jk,jn) = pta(ji,ijd,jk,jn) + ztra 274 END DO 275 ! ! down-slope T-point (deep bottom point) 276 zbtr = e1e2t_r(ji,ijd) / fse3t(ji,ijd,ikvd) 277 ztra = zv_bbl * ( ptb(ji,ijs,ikvs,jn) - ptb(ji,ijd,ikvd,jn) ) * zbtr 278 pta(ji,ijd,ikvd,jn) = pta(ji,ijd,ikvd,jn) + ztra 279 ENDIF 280 END DO 281 ! 282 END DO 283 ! ! =========== 284 END DO ! end tracer 285 ! ! =========== 286 END SUBROUTINE tra_bbl_adv 287 288 289 SUBROUTINE bbl( kt, cdtype ) 290 !!---------------------------------------------------------------------- 291 !! *** ROUTINE bbl *** 292 !! 293 !! ** Purpose : Computes the bottom boundary horizontal and vertical 294 !! advection terms. 295 !! 296 !! ** Method : 297 !! * diffusive bbl (nn_bbl_ldf=1) : 298 !! When the product grad( rho) * grad(h) < 0 (where grad is an 299 !! along bottom slope gradient) an additional lateral 2nd order 300 !! diffusion along the bottom slope is added to the general 301 !! tracer trend, otherwise the additional trend is set to 0. 302 !! A typical value of ahbt is 2000 m2/s (equivalent to 303 !! a downslope velocity of 20 cm/s if the condition for slope 304 !! convection is satified) 305 !! * advective bbl (nn_bbl_adv=1 or 2) : 306 !! nn_bbl_adv = 1 use of the ocean velocity as bbl velocity 307 !! nn_bbl_adv = 2 follow Campin and Goosse (1999) implentation 308 !! i.e. transport proportional to the along-slope density gradient 309 !! 310 !! NB: the along slope density gradient is evaluated using the 311 !! local density (i.e. referenced at a common local depth). 312 !! 313 !! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591. 314 !! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430. 315 !!---------------------------------------------------------------------- 316 INTEGER , INTENT(in ) :: kt ! ocean time-step index 317 CHARACTER(len=3), INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 318 !! 319 INTEGER :: ji, jj ! dummy loop indices 320 INTEGER :: ik ! local integers 321 INTEGER :: iis , iid , ijs , ijd ! - - 322 INTEGER :: ikus, ikud, ikvs, ikvd ! - - 323 REAL(wp) :: zsign, zsigna, zgbbl ! local scalars 324 REAL(wp) :: zgdrho, zt, zs, zh ! - - 325 REAL(wp), DIMENSION(jpi,jpj) :: zub, zvb, ztb, zsb, zdep ! 2D workspace 326 !! 327 REAL(wp) :: fsalbt, fsbeta, pft, pfs, pfh ! statement function 328 !!----------------------- zv_bbl----------------------------------------------- 329 ! ratio alpha/beta = fsalbt : ratio of thermal over saline expension coefficients 330 ! ================ pft : potential temperature in degrees celcius 331 ! pfs : salinity anomaly (s-35) in psu 332 ! pfh : depth in meters 333 ! nn_eos = 0 (Jackett and McDougall 1994 formulation) 334 fsalbt( pft, pfs, pfh ) = & ! alpha/beta 131 335 ( ( ( -0.255019e-07 * pft + 0.298357e-05 ) * pft & 132 &- 0.203814e-03 ) * pft &133 &+ 0.170907e-01 ) * pft &134 &+ 0.665157e-01 &336 - 0.203814e-03 ) * pft & 337 + 0.170907e-01 ) * pft & 338 + 0.665157e-01 & 135 339 +(-0.678662e-05 * pfs - 0.846960e-04 * pft + 0.378110e-02 ) * pfs & 136 340 + ( ( - 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(:,:,:) 341 - 0.251520e-11 * pfs & 342 + 0.512857e-12 * pft * pft ) * pfh & 343 - 0.164759e-06 * pfs & 344 +( 0.791325e-08 * pft - 0.933746e-06 ) * pft & 345 + 0.380374e-04 ) * pfh 346 fsbeta( pft, pfs, pfh ) = & ! beta 347 ( ( -0.415613e-09 * pft + 0.555579e-07 ) * pft & 348 - 0.301985e-05 ) * pft & 349 + 0.785567e-03 & 350 + ( 0.515032e-08 * pfs & 351 + 0.788212e-08 * pft - 0.356603e-06 ) * pfs & 352 +( ( 0.121551e-17 * pfh & 353 - 0.602281e-15 * pfs & 354 - 0.175379e-14 * pft + 0.176621e-12 ) * pfh & 355 + 0.408195e-10 * pfs & 356 + ( - 0.213127e-11 * pft + 0.192867e-09 ) * pft & 357 - 0.121555e-07 ) * pfh 358 !!---------------------------------------------------------------------- 359 360 IF( kt == nit000 ) THEN 361 IF(lwp) WRITE(numout,*) 362 IF(lwp) WRITE(numout,*) 'trabbl:bbl : Compute bbl velocities and diffusive coefficients in ', cdtype 363 IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 149 364 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 155 DO jj = 1, 1 156 DO ji = 1, jpij ! vector opt. (forced unrolling) 157 # else 365 366 ! !* bottom temperature, salinity, velocity and depth 367 #if defined key_vectopt_loop 368 DO jj = 1, 1 ! vector opt. (forced unrolling) 369 DO ji = 1, jpij 370 #else 158 371 DO jj = 1, jpj 159 372 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 373 #endif 374 ik = mbkt(ji,jj) ! bottom T-level index 375 ztb (ji,jj) = tsb(ji,jj,ik,jp_tem) * tmask(ji,jj,1) ! bottom before T and S 376 zsb (ji,jj) = tsb(ji,jj,ik,jp_sal) * tmask(ji,jj,1) 377 zdep(ji,jj) = fsdept_0(ji,jj,ik) ! bottom T-level reference depth 378 ! 379 zub(ji,jj) = un(ji,jj,mbku(ji,jj)) ! bottom velocity 380 zvb(ji,jj) = vn(ji,jj,mbkv(ji,jj)) 167 381 END DO 168 382 END DO 169 170 IF( ln_zps ) THEN ! partial steps correction 171 # if defined key_vectopt_loop 172 DO jj = 1, 1 173 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 174 # else 175 DO jj = 1, jpjm1 176 DO ji = 1, jpim1 177 # endif 178 iku1 = MAX( mbathy(ji+1,jj )-1, 1 ) 179 iku2 = MAX( mbathy(ji ,jj )-1, 1 ) 180 ikv1 = MAX( mbathy(ji ,jj+1)-1, 1 ) 181 ikv2 = MAX( mbathy(ji ,jj )-1, 1 ) 182 ze3u = MIN( fse3u(ji,jj,iku1), fse3u(ji,jj,iku2) ) 183 ze3v = MIN( fse3v(ji,jj,ikv1), fse3v(ji,jj,ikv2) ) 184 zahu(ji,jj) = rn_ahtbbl * e2u(ji,jj) * ze3u / e1u(ji,jj) * umask(ji,jj,1) 185 zahv(ji,jj) = rn_ahtbbl * e1v(ji,jj) * ze3v / e2v(ji,jj) * vmask(ji,jj,1) 186 END DO 187 END DO 188 ELSE ! z-coordinate - full steps or s-coordinate 189 # if defined key_vectopt_loop 190 DO jj = 1, 1 191 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 192 # else 193 DO jj = 1, jpjm1 194 DO ji = 1, jpim1 195 # endif 196 iku = mbku(ji,jj) 197 ikv = mbkv(ji,jj) 198 zahu(ji,jj) = rn_ahtbbl * e2u(ji,jj) * fse3u(ji,jj,iku) / e1u(ji,jj) * umask(ji,jj,1) 199 zahv(ji,jj) = rn_ahtbbl * e1v(ji,jj) * fse3v(ji,jj,ikv) / e2v(ji,jj) * vmask(ji,jj,1) 200 END DO 201 END DO 202 ENDIF 203 204 ! 1. Criteria of additional bottom diffusivity: grad(rho).grad(h)<0 205 ! -------------------------------------------- 206 ! Sign of the local density gradient along the i- and j-slopes 207 ! multiplied by the slope of the ocean bottom 208 209 SELECT CASE ( nn_eos ) 210 ! 211 CASE ( 0 ) !== Jackett and McDougall (1994) formulation ==! 212 # if defined key_vectopt_loop 213 DO jj = 1, 1 214 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 215 # else 216 DO jj = 1, jpjm1 217 DO ji = 1, jpim1 218 # endif 219 ! temperature, salinity anomalie and depth 220 zt = 0.5 * ( ztnb(ji,jj) + ztnb(ji+1,jj) ) 221 zs = 0.5 * ( zsnb(ji,jj) + zsnb(ji+1,jj) ) - 35.0 383 384 ! !-------------------! 385 IF( nn_bbl_ldf == 1 ) THEN ! diffusive bbl ! 386 ! !-------------------! 387 DO jj = 1, jpjm1 ! (criteria for non zero flux: grad(rho).grad(h) < 0 ) 388 DO ji = 1, jpim1 389 ! ! i-direction 390 zt = 0.5 * ( ztb (ji,jj) + ztb (ji+1,jj) ) ! T, S anomalie, and depth 391 zs = 0.5 * ( zsb (ji,jj) + zsb (ji+1,jj) ) - 35.0 222 392 zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) ) 223 ! masked ratio alpha/beta 224 zalbet = fsalbt( zt, zs, zh )*umask(ji,jj,1) 225 ! local density gradient along i-bathymetric slope 226 zgdrho = zalbet * ( ztnb(ji+1,jj) - ztnb(ji,jj) ) & 227 - ( zsnb(ji+1,jj) - zsnb(ji,jj) ) 228 ! sign of local i-gradient of density multiplied by the i-slope 229 zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) ) 230 zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj) 393 ! ! masked bbl i-gradient of density 394 zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji+1,jj) - ztb(ji,jj) ) & 395 & - ( zsb(ji+1,jj) - zsb(ji,jj) ) ) * umask(ji,jj,1) 396 ! 397 zsign = SIGN( 0.5, - zgdrho * REAL( mgrhu(ji,jj) ) ) ! sign of ( i-gradient * i-slope ) 398 ahu_bbl(ji,jj) = ( 0.5 - zsign ) * ahu_bbl_0(ji,jj) ! masked diffusive flux coeff. 231 399 ! 232 ! temperature, salinity anomalie and depth233 zt = 0.5 * ( zt nb(ji,jj+1) + ztnb(ji,jj) )234 zs = 0.5 * ( zs nb(ji,jj+1) + zsnb(ji,jj) ) - 35.0400 ! ! j-direction 401 zt = 0.5 * ( ztb (ji,jj+1) + ztb (ji,jj) ) ! T, S anomalie, and depth 402 zs = 0.5 * ( zsb (ji,jj+1) + zsb (ji,jj) ) - 35.0 235 403 zh = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) ) 236 ! masked ratio alpha/beta 237 zalbet = fsalbt( zt, zs, zh )*vmask(ji,jj,1) 238 ! local density gradient along j-bathymetric slope 239 zgdrho = zalbet * ( ztnb(ji,jj+1) - ztnb(ji,jj) ) & 240 - ( zsnb(ji,jj+1) - zsnb(ji,jj) ) 241 ! sign of local j-gradient of density multiplied by the j-slope 242 zsign = sign( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) ) 243 zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj) 404 ! ! masked bbl j-gradient of density 405 zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji,jj+1) - ztb(ji,jj) ) & 406 & - ( zsb(ji,jj+1) - zsb(ji,jj) ) ) * vmask(ji,jj,1) 407 ! 408 zsign = SIGN( 0.5, -zgdrho * REAL( mgrhv(ji,jj) ) ) ! sign of ( j-gradient * j-slope ) 409 ahv_bbl(ji,jj) = ( 0.5 - zsign ) * ahv_bbl_0(ji,jj) 410 ! 244 411 END DO 245 412 END DO 246 413 ! 247 CASE ( 1 ) !== Linear formulation function of temperature only ==! 248 # if defined key_vectopt_loop 249 DO jj = 1, 1 250 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 251 # else 252 DO jj = 1, jpjm1 253 DO ji = 1, jpim1 254 # endif 255 ! local 'density/temperature' gradient along i-bathymetric slope 256 zgdrho = ztnb(ji+1,jj) - ztnb(ji,jj) 257 ! sign of local i-gradient of density multiplied by the i-slope 258 zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) ) 259 zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj) 260 ! 261 ! local density gradient along j-bathymetric slope 262 zgdrho = ztnb(ji,jj+1) - ztnb(ji,jj) 263 ! sign of local j-gradient of density multiplied by the j-slope 264 zsign = sign( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) ) 265 zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj) 414 ENDIF 415 416 ! !-------------------! 417 IF( nn_bbl_adv /= 0 ) THEN ! advective bbl ! 418 ! !-------------------! 419 SELECT CASE ( nn_bbl_adv ) !* bbl transport type 420 ! 421 CASE( 1 ) != use of upper velocity 422 DO jj = 1, jpjm1 ! criteria: grad(rho).grad(h)<0 and grad(rho).grad(h)<0 423 DO ji = 1, fs_jpim1 ! vector opt. 424 ! ! i-direction 425 zt = 0.5 * ( ztb (ji,jj) + ztb (ji+1,jj) ) ! T, S anomalie, and depth 426 zs = 0.5 * ( zsb (ji,jj) + zsb (ji+1,jj) ) - 35.0 427 zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) ) 428 ! ! masked bbl i-gradient of density 429 zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji+1,jj) - ztb(ji,jj) ) & 430 & - ( zsb(ji+1,jj) - zsb(ji,jj) ) ) * umask(ji,jj,1) 431 ! 432 zsign = SIGN( 0.5, - zgdrho * REAL( mgrhu(ji,jj) ) ) ! sign of i-gradient * i-slope 433 zsigna= SIGN( 0.5, zub(ji,jj) * REAL( mgrhu(ji,jj) ) ) ! sign of u * i-slope 434 ! 435 ! ! bbl velocity 436 utr_bbl(ji,jj) = ( 0.5 + zsigna ) * ( 0.5 - zsign ) * e2u(ji,jj) * e3u_bbl_0(ji,jj) * zub(ji,jj) 437 ! 438 ! ! j-direction 439 zt = 0.5 * ( ztb (ji,jj+1) + ztb (ji,jj) ) ! T, S anomalie, and depth 440 zs = 0.5 * ( zsb (ji,jj+1) + zsb (ji,jj) ) - 35.0 441 zh = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) ) 442 ! ! masked bbl j-gradient of density 443 zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji,jj+1) - ztb(ji,jj) ) & 444 & - ( zsb(ji,jj+1) - zsb(ji,jj) ) ) * vmask(ji,jj,1) 445 zsign = SIGN( 0.5, - zgdrho * REAL( mgrhv(ji,jj) ) ) ! sign of j-gradient * j-slope 446 zsigna= SIGN( 0.5, zvb(ji,jj) * REAL( mgrhv(ji,jj) ) ) ! sign of u * i-slope 447 ! 448 ! ! bbl velocity 449 vtr_bbl(ji,jj) = ( 0.5 + zsigna ) * ( 0.5 - zsign ) * e1v(ji,jj) * e3v_bbl_0(ji,jj) * zvb(ji,jj) 450 END DO 266 451 END DO 452 ! 453 CASE( 2 ) != bbl velocity = F( delta rho ) 454 zgbbl = grav * rn_gambbl 455 DO jj = 1, jpjm1 ! criteria: rho_up > rho_down 456 DO ji = 1, fs_jpim1 ! vector opt. 457 ! ! i-direction 458 ! down-slope T-point i/k-index (deep) & up-slope T-point i/k-index (shelf) 459 iid = ji + MAX( 0, mgrhu(ji,jj) ) ; iis = ji + 1 - MAX( 0, mgrhu(ji,jj) ) 460 ikud = mbku_d(ji,jj) ; ikus = mbku(ji,jj) 461 ! 462 ! ! mid-depth density anomalie (up-slope minus down-slope) 463 zt = 0.5 * ( ztb (ji,jj) + ztb (ji+1,jj) ) ! mid slope depth of T, S, and depth 464 zs = 0.5 * ( zsb (ji,jj) + zsb (ji+1,jj) ) - 35.0 465 zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) ) 466 zgdrho = fsbeta( zt, zs, zh ) & 467 & * ( fsalbt( zt, zs, zh ) * ( ztb(iid,jj) - ztb(iis,jj) ) & 468 & - ( zsb(iid,jj) - zsb(iis,jj) ) ) * umask(ji,jj,1) 469 zgdrho = MAX( 0.e0, zgdrho ) ! only if shelf is denser than deep 470 ! 471 ! ! bbl transport (down-slope direction) 472 utr_bbl(ji,jj) = e2u(ji,jj) * e3u_bbl_0(ji,jj) * zgbbl * zgdrho * REAL( mgrhu(ji,jj) ) 473 ! 474 ! ! j-direction 475 ! down-slope T-point j/k-index (deep) & of the up -slope T-point j/k-index (shelf) 476 ijd = jj + MAX( 0, mgrhv(ji,jj) ) ; ijs = jj + 1 - MAX( 0, mgrhv(ji,jj) ) 477 ikvd = mbkv_d(ji,jj) ; ikvs = mbkv(ji,jj) 478 ! 479 ! ! mid-depth density anomalie (up-slope minus down-slope) 480 zt = 0.5 * ( ztb (ji,jj) + ztb (ji,jj+1) ) ! mid slope depth of T, S, and depth 481 zs = 0.5 * ( zsb (ji,jj) + zsb (ji,jj+1) ) - 35.0 482 zh = 0.5 * ( zdep(ji,jj) + zdep(ji,jj+1) ) 483 zgdrho = fsbeta( zt, zs, zh ) & 484 & * ( fsalbt( zt, zs, zh ) * ( ztb(ji,ijd) - ztb(ji,ijs) ) & 485 & - ( zsb(ji,ijd) - zsb(ji,ijs) ) ) * vmask(ji,jj,1) 486 zgdrho = MAX( 0.e0, zgdrho ) ! only if shelf is denser than deep 487 ! 488 ! ! bbl transport (down-slope direction) 489 vtr_bbl(ji,jj) = e1v(ji,jj) * e3v_bbl_0(ji,jj) * zgbbl * zgdrho * REAL( mgrhv(ji,jj) ) 490 END DO 491 END DO 492 END SELECT 493 ! 494 ENDIF 495 ! 496 END SUBROUTINE bbl 497 498 499 SUBROUTINE tra_bbl_init 500 !!---------------------------------------------------------------------- 501 !! *** ROUTINE tra_bbl_init *** 502 !! 503 !! ** Purpose : Initialization for the bottom boundary layer scheme. 504 !! 505 !! ** Method : Read the nambbl namelist and check the parameters 506 !! called by tra_bbl at the first timestep (nit000) 507 !!---------------------------------------------------------------------- 508 INTEGER :: ji, jj ! dummy loop indices 509 INTEGER :: ii0, ii1, ij0, ij1 ! temporary integer 510 REAL(wp), DIMENSION(jpi,jpj) :: zmbk ! 2D workspace 511 !! 512 NAMELIST/nambbl/ nn_bbl_ldf, nn_bbl_adv, rn_ahtbbl, rn_gambbl 513 !!---------------------------------------------------------------------- 514 515 REWIND ( numnam ) !* Read Namelist nambbl : bottom boundary layer scheme 516 READ ( numnam, nambbl ) 517 ! 518 l_bbl = .TRUE. !* flag to compute bbl coef and transport 519 ! 520 IF(lwp) THEN !* Parameter control and print 521 WRITE(numout,*) 522 WRITE(numout,*) 'tra_bbl_init : bottom boundary layer initialisation' 523 WRITE(numout,*) '~~~~~~~~~~~~' 524 WRITE(numout,*) ' Namelist nambbl : set bbl parameters' 525 WRITE(numout,*) ' diffusive bbl (=1) or not (=0) nn_bbl_ldf = ', nn_bbl_ldf 526 WRITE(numout,*) ' advective bbl (=1/2) or not (=0) nn_bbl_adv = ', nn_bbl_adv 527 WRITE(numout,*) ' diffusive bbl coefficient rn_ahtbbl = ', rn_ahtbbl, ' m2/s' 528 WRITE(numout,*) ' advective bbl coefficient rn_gambbl = ', rn_gambbl, ' s' 529 ENDIF 530 531 IF( nn_bbl_adv == 1 ) WRITE(numout,*) ' * Advective BBL using upper velocity' 532 IF( nn_bbl_adv == 2 ) WRITE(numout,*) ' * Advective BBL using velocity = F( delta rho)' 533 534 IF( nn_eos /= 0 ) CALL ctl_stop ( ' bbl parameterisation requires eos = 0. We stop.' ) 535 536 537 ! !* inverse of surface of T-cells 538 e1e2t_r(:,:) = 1.0 / ( e1t(:,:) * e2t(:,:) ) 539 540 ! !* vertical index of "deep" bottom u- and v-points 541 DO jj = 1, jpjm1 ! (the "shelf" bottom k-indices are mbku and mbkv) 542 DO ji = 1, jpim1 543 mbku_d(ji,jj) = MAX( mbkt(ji+1,jj ) , mbkt(ji,jj) ) ! >= 1 as mbkt=1 over land 544 mbkv_d(ji,jj) = MAX( mbkt(ji ,jj+1) , mbkt(ji,jj) ) 267 545 END DO 268 ! 269 CASE ( 2 ) !== Linear formulation function of temperature and salinity ==! 270 # if defined key_vectopt_loop 271 DO jj = 1, 1 272 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 273 # else 274 DO jj = 1, jpjm1 275 DO ji = 1, jpim1 276 # endif 277 ! local density gradient along i-bathymetric slope 278 zgdrho = - ( rn_beta *( zsnb(ji+1,jj) - zsnb(ji,jj) ) & 279 & - rn_alpha*( ztnb(ji+1,jj) - ztnb(ji,jj) ) ) 280 ! sign of local i-gradient of density multiplied by the i-slope 281 zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) ) 282 zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj) 283 ! 284 ! local density gradient along j-bathymetric slope 285 zgdrho = - ( rn_beta *( zsnb(ji,jj+1) - zsnb(ji,jj) ) & 286 & - rn_alpha*( ztnb(ji,jj+1) - ztnb(ji,jj) ) ) 287 ! sign of local j-gradient of density multiplied by the j-slope 288 zsign = sign( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) ) 289 zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj) 290 END DO 546 END DO 547 ! converte into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk 548 zmbk(:,:) = REAL( mbku_d(:,:), wp ) ; CALL lbc_lnk(zmbk,'U',1.) ; mbku_d(:,:) = MAX( INT( zmbk(:,:) ), 1 ) 549 zmbk(:,:) = REAL( mbkv_d(:,:), wp ) ; CALL lbc_lnk(zmbk,'V',1.) ; mbkv_d(:,:) = MAX( INT( zmbk(:,:) ), 1 ) 550 551 !* sign of grad(H) at u- and v-points 552 mgrhu(jpi,:) = 0. ; mgrhu(:,jpj) = 0. ; mgrhv(jpi,:) = 0. ; mgrhv(:,jpj) = 0. 553 DO jj = 1, jpjm1 554 DO ji = 1, jpim1 555 mgrhu(ji,jj) = INT( SIGN( 1.e0, fsdept_0(ji+1,jj,mbkt(ji+1,jj)) - fsdept_0(ji,jj,mbkt(ji,jj)) ) ) 556 mgrhv(ji,jj) = INT( SIGN( 1.e0, fsdept_0(ji,jj+1,mbkt(ji,jj+1)) - fsdept_0(ji,jj,mbkt(ji,jj)) ) ) 291 557 END DO 292 !293 END SELECT294 295 ! 2. Additional second order diffusive trends296 ! -------------------------------------------297 298 ! first derivative (gradient)299 # if defined key_vectopt_loop300 jj = 1301 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling)302 # else303 DO jj = 1, jpjm1304 DO ji = 1, jpim1305 # endif306 zkx(ji,jj) = zki(ji,jj) * ( ztbb(ji+1,jj) - ztbb(ji,jj) )307 zkz(ji,jj) = zki(ji,jj) * ( zsbb(ji+1,jj) - zsbb(ji,jj) )308 309 zky(ji,jj) = zkj(ji,jj) * ( ztbb(ji,jj+1) - ztbb(ji,jj) )310 zkw(ji,jj) = zkj(ji,jj) * ( zsbb(ji,jj+1) - zsbb(ji,jj) )311 # if ! defined key_vectopt_loop312 END DO313 # endif314 558 END DO 315 559 316 IF( cp_cfg == "orca" ) THEN 560 DO jj = 1, jpjm1 !* bbl thickness at u- (v-) point 561 DO ji = 1, jpim1 ! minimum of top & bottom e3u_0 (e3v_0) 562 e3u_bbl_0(ji,jj) = MIN( fse3u_0(ji,jj,mbkt(ji+1,jj )), fse3u_0(ji,jj,mbkt(ji,jj)) ) 563 e3v_bbl_0(ji,jj) = MIN( fse3v_0(ji,jj,mbkt(ji ,jj+1)), fse3v_0(ji,jj,mbkt(ji,jj)) ) 564 END DO 565 END DO 566 CALL lbc_lnk( e3u_bbl_0, 'U', 1. ) ; CALL lbc_lnk( e3v_bbl_0, 'V', 1. ) ! lateral boundary conditions 567 568 ! !* masked diffusive flux coefficients 569 ahu_bbl_0(:,:) = rn_ahtbbl * e2u(:,:) * e3u_bbl_0(:,:) / e1u(:,:) * umask(:,:,1) 570 ahv_bbl_0(:,:) = rn_ahtbbl * e1v(:,:) * e3v_bbl_0(:,:) / e2v(:,:) * vmask(:,:,1) 571 572 573 IF( cp_cfg == "orca" ) THEN !* ORCA configuration : regional enhancement of ah_bbl 317 574 ! 318 575 SELECT CASE ( jp_cfg ) 319 ! ! ======================= 320 CASE ( 2 ) ! ORCA_R2 configuration 321 ! ! ======================= 322 ! Gibraltar enhancement of BBL 323 ij0 = 102 ; ij1 = 102 576 CASE ( 2 ) ! ORCA_R2 577 ij0 = 102 ; ij1 = 102 ! Gibraltar enhancement of BBL 324 578 ii0 = 139 ; ii1 = 140 325 zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1))326 zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1))579 ahu_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) = 4.e0*ahu_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) 580 ahv_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) = 4.e0*ahv_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) 327 581 ! 328 ! Red Sea enhancement of BBL 329 ij0 = 88 ; ij1 = 88 582 ij0 = 88 ; ij1 = 88 ! Red Sea enhancement of BBL 330 583 ii0 = 161 ; ii1 = 162 331 zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 10.e0 * zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1))332 zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 10.e0 * zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1))584 ahu_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) = 10.e0*ahu_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) 585 ahv_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) = 10.e0*ahv_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) 333 586 ! 334 ! ! ======================= 335 CASE ( 4 ) ! ORCA_R4 configuration 336 ! ! ======================= 337 ! Gibraltar enhancement of BBL 338 ij0 = 52 ; ij1 = 52 587 CASE ( 4 ) ! ORCA_R4 588 ij0 = 52 ; ij1 = 52 ! Gibraltar enhancement of BBL 339 589 ii0 = 70 ; ii1 = 71 340 zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) 341 zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) 342 ! 590 ahu_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) = 4.e0*ahu_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) 591 ahv_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) = 4.e0*ahv_bbl_0(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1)) 343 592 END SELECT 344 !593 ! 345 594 ENDIF 346 347 348 ! second derivative (divergence) and add to the general tracer trend349 # if defined key_vectopt_loop350 DO jj = 1, 1351 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling)352 # else353 DO jj = 2, jpjm1354 DO ji = 2, jpim1355 # endif356 ik = max( mbathy(ji,jj)-1, 1 )357 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,ik) )358 zta = ( zkx(ji,jj) - zkx(ji-1,jj ) &359 + zky(ji,jj) - zky(ji ,jj-1) ) * zbtr360 zsa = ( zkz(ji,jj) - zkz(ji-1,jj ) &361 + zkw(ji,jj) - zkw(ji ,jj-1) ) * zbtr362 ta(ji,jj,ik) = ta(ji,jj,ik) + zta363 sa(ji,jj,ik) = sa(ji,jj,ik) + zsa364 END DO365 END DO366 367 IF( l_trdtra ) THEN ! save the BBL lateral diffusion trends for diagnostic368 ztrdt(:,:,:) = ta(:,:,:) - ztrdt(:,:,:)369 ztrds(:,:,:) = sa(:,:,:) - ztrds(:,:,:)370 CALL trd_mod(ztrdt, ztrds, jptra_trd_bbl, 'TRA', kt)371 ENDIF372 373 IF(ln_ctl) CALL prt_ctl( tab3d_1=ta, clinfo1=' bbl - Ta: ', mask1=tmask, &374 & tab3d_2=sa, clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' )375 !376 END SUBROUTINE tra_bbl_dif377 378 # if defined key_trabbl_adv379 !!----------------------------------------------------------------------380 !! 'key_trabbl_adv' advective bottom boundary layer381 !!----------------------------------------------------------------------382 # include "trabbl_adv.h90"383 # else384 !!----------------------------------------------------------------------385 !! Default option : NO advective bottom boundary layer386 !!----------------------------------------------------------------------387 SUBROUTINE tra_bbl_adv (kt ) ! Empty routine388 INTEGER, INTENT(in) :: kt389 WRITE(*,*) 'tra_bbl_adv: You should not have seen this print! error?', kt390 END SUBROUTINE tra_bbl_adv391 # endif392 393 SUBROUTINE tra_bbl_init394 !!----------------------------------------------------------------------395 !! *** ROUTINE tra_bbl_init ***396 !!397 !! ** Purpose : Initialization for the bottom boundary layer scheme.398 !!399 !! ** Method : Read the nambbl namelist and check the parameters400 !! called by tra_bbl at the first timestep (nit000)401 !!----------------------------------------------------------------------402 INTEGER :: ji, jj ! dummy loop indices403 REAL(wp), DIMENSION(jpi,jpj) :: zmbk404 405 NAMELIST/nambbl/ rn_ahtbbl406 !!----------------------------------------------------------------------407 408 REWIND ( numnam ) ! Read Namelist nambbl : bottom boundary layer scheme409 READ ( numnam, nambbl )410 411 IF(lwp) THEN ! Parameter control and print412 WRITE(numout,*)413 WRITE(numout,*) 'tra_bbl_init : '414 WRITE(numout,*) '~~~~~~~~~~~~'415 IF( lk_trabbl_dif ) WRITE(numout,*) ' * Diffusive Bottom Boundary Layer'416 IF( lk_trabbl_adv ) WRITE(numout,*) ' * Advective Bottom Boundary Layer'417 WRITE(numout,*) ' Namelist nambbl : set bbl parameters'418 WRITE(numout,*) ' bottom boundary layer coef. rn_ahtbbl = ', rn_ahtbbl419 ENDIF420 421 DO jj = 1, jpj422 DO ji = 1, jpi423 mbkt(ji,jj) = MAX( mbathy(ji,jj) - 1, 1 ) ! vertical index of the bottom ocean T-level424 END DO425 END DO426 DO jj = 1, jpjm1427 DO ji = 1, jpim1428 mbku(ji,jj) = MAX( MIN( mbathy(ji+1,jj ), mbathy(ji,jj) ) - 1, 1 )429 mbkv(ji,jj) = MAX( MIN( mbathy(ji ,jj+1), mbathy(ji,jj) ) - 1, 1 )430 END DO431 END DO432 433 zmbk(:,:) = FLOAT( mbku (:,:) )434 CALL lbc_lnk(zmbk,'U',1.)435 mbku(:,:) = MAX( INT( zmbk(:,:) ), 1 )436 437 zmbk(:,:) = FLOAT( mbkv (:,:) )438 CALL lbc_lnk(zmbk,'V',1.)439 mbkv(:,:) = MAX( INT( zmbk(:,:) ), 1 )440 441 # if defined key_trabbl_adv442 w_bbl(:,:,:) = 0.e0 ! initialisation of w_bbl to zero443 # endif444 595 ! 445 596 END SUBROUTINE tra_bbl_init … … 449 600 !! Dummy module : No bottom boundary layer scheme 450 601 !!---------------------------------------------------------------------- 451 LOGICAL, PUBLIC, PARAMETER :: lk_trabbl_dif = .FALSE. !: diff bbl flag 452 LOGICAL, PUBLIC, PARAMETER :: lk_trabbl_adv = .FALSE. !: adv bbl flag 602 LOGICAL, PUBLIC, PARAMETER :: lk_trabbl = .FALSE. !: bbl flag 453 603 CONTAINS 454 SUBROUTINE tra_bbl_dif( kt ) ! Empty routine 455 WRITE(*,*) 'tra_bbl_dif: You should not have seen this print! error?', kt 456 END SUBROUTINE tra_bbl_dif 457 SUBROUTINE tra_bbl_adv( kt ) ! Empty routine 458 WRITE(*,*) 'tra_bbl_adv: You should not have seen this print! error?', kt 459 END SUBROUTINE tra_bbl_adv 604 SUBROUTINE tra_bbl_init ! Dummy routine 605 END SUBROUTINE tra_bbl_init 606 SUBROUTINE tra_bbl( kt ) ! Dummy routine 607 WRITE(*,*) 'tra_bbl: You should not have seen this print! error?', kt 608 END SUBROUTINE tra_bbl 460 609 #endif 461 610
Note: See TracChangeset
for help on using the changeset viewer.
