Changeset 2528 for trunk/NEMOGCM/NEMO/TOP_SRC/TRP/trcbbl.F90
- Timestamp:
- 2010-12-27T18:33:53+01:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/TOP_SRC/TRP/trcbbl.F90
- Property svn:executable deleted
r1606 r2528 1 1 MODULE trcbbl 2 2 !!====================================================================== 3 3 !! *** MODULE trcbbl *** 4 4 !! Ocean passive tracers physics : advective and/or diffusive bottom boundary 5 5 !! layer scheme 6 6 !!====================================================================== 7 !! History : 8.0 ! 96-06 (L. Mortier) Original code 8 !! 8.0 ! 97-11 (G. Madec) Optimization 9 !! 8.5 ! 02-08 (G. Madec) free form + modules 10 !! 9.0 ! 04-03 (C. Ethe) Adaptation for passive tracers 11 !! ! 07-02 (C. Deltel) Diagnose ML trends for passive tracers 7 !!============================================================================== 8 !! History : OPA ! 1996-06 (L. Mortier) Original code 9 !! 8.0 ! 1997-11 (G. Madec) Optimization 10 !! NEMO 1.0 ! 2002-08 (G. Madec) free form + modules 11 !! - ! 2004-01 (A. de Miranda, G. Madec, J.M. Molines ) add advective bbl 12 !! 3.3 ! 2009-11 (G. Madec) merge trabbl and trabbl_adv + style + optimization 13 !! - ! 2010-04 (G. Madec) Campin & Goosse advective bbl 14 !! - ! 2010-06 (C. Ethe, G. Madec) merge TRA-TRC 12 15 !!---------------------------------------------------------------------- 13 #if defined key_top && ( defined key_trcbbl_dif || defined key_trcbbl_adv ) && ! defined key_c1d16 #if defined key_top && defined key_trabbl 14 17 !!---------------------------------------------------------------------- 15 !! 'key_trcbbl_dif' or diffusive bottom boundary layer 16 !! 'key_trcbbl_adv' advective bottom boundary layer 18 !! 'key_trabbl diffusive or/and adevective bottom boundary layer 17 19 !!---------------------------------------------------------------------- 18 !! trc_bbl_dif : update the passive tracer trends due to the bottom 19 !! boundary layer (diffusive only) 20 !! trc_bbl_adv : update the passive tracer trends due to the bottom 21 !! boundary layer (advective and/or diffusive) 20 !! trc_bbl : update the tracer trends due to the bottom boundary layer (advective and/or diffusive) 22 21 !!---------------------------------------------------------------------- 23 22 USE oce_trc ! ocean dynamics and active tracers variables 24 23 USE trc ! ocean passive tracers variables 25 USE trctrp_lec ! passive tracers transport 24 USE trcnam_trp ! passive tracers transport namelist variables 25 USE trabbl ! 26 26 USE prtctl_trc ! Print control for debbuging 27 USE eosbn2 28 USE lbclnk 29 USE trdmld_trc 30 USE trdmld_trc_oce 27 USE trdmod_oce 28 USE trdtra 31 29 32 IMPLICIT NONE 33 PRIVATE 30 PUBLIC trc_bbl ! routine called by step.F90 34 31 35 PUBLIC trc_bbl_dif ! routine called by step.F9036 PUBLIC trc_bbl_adv ! routine called by step.F9037 38 # if defined key_trcbbl_dif39 LOGICAL, PUBLIC, PARAMETER :: lk_trcbbl_dif = .TRUE. !: diffusive bottom boundary layer flag40 # else41 LOGICAL, PUBLIC, PARAMETER :: lk_trcbbl_dif = .FALSE. !: diffusive bottom boundary layer flag42 # endif43 44 # if defined key_trcbbl_adv45 LOGICAL, PUBLIC, PARAMETER :: lk_trcbbl_adv = .TRUE. !: advective bottom boundary layer flag46 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: u_trc_bbl !: veloc. involved in the advective BBL47 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: v_trc_bbl !: veloc. involved in the advective BBL48 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: w_trc_bbl !: vertic. increment of veloc. due to adv. BBL49 ! ! only affect tracer vertical advection50 # else51 LOGICAL, PUBLIC, PARAMETER :: lk_trcbbl_adv = .FALSE. !: advective bottom boundary layer flag52 # endif53 54 INTEGER, DIMENSION(jpi,jpj) :: mbkt, mbku, mbkv55 32 56 33 !! * Substitutions 57 34 # include "top_substitute.h90" 58 35 !!---------------------------------------------------------------------- 59 !! TOP 1.0 , LOCEAN-IPSL (2005)60 !! $ Header: /home/opalod/NEMOCVSROOT/NEMO/TOP_SRC/TRP/trcbbl.F90,v 1.12 2006/09/12 11:10:13 opalod Exp$61 !! Software governed by the CeCILL licence ( modipsl/doc/NEMO_CeCILL.txt)36 !! NEMO/TOP 3.3 , NEMO Consortium (2010) 37 !! $Id$ 38 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 62 39 !!---------------------------------------------------------------------- 63 40 64 41 CONTAINS 65 42 66 SUBROUTINE trc_bbl_dif( kt ) 43 44 SUBROUTINE trc_bbl( kt ) 67 45 !!---------------------------------------------------------------------- 68 !! *** ROUTINE trc_bbl_dif *** 46 !! *** ROUTINE bbl *** 47 !! 48 !! ** Purpose : Compute the before tracer (t & s) trend associated 49 !! with the bottom boundary layer and add it to the general trend 50 !! of tracer equations. 69 51 !! 70 !! ** Purpose : Compute the before tracer trend associated 71 !! with the bottom boundary layer and add it to the general trend 72 !! of tracer equations. The bottom boundary layer is supposed to be 73 !! a purely diffusive bottom boundary layer. 74 !! 75 !! ** Method : When the product grad( rho) * grad(h) < 0 (where grad 76 !! is an along bottom slope gradient) an additional lateral diffu- 77 !! sive trend along the bottom slope is added to the general tracer 78 !! trend, otherwise nothing is done. 79 !! Second order operator (laplacian type) with variable coefficient 80 !! computed as follow for temperature (idem on s): 81 !! difft = 1/(e1t*e2t*e3t) { di-1[ ahbt e2u*e3u/e1u di[ztb] ] 82 !! + dj-1[ ahbt e1v*e3v/e2v dj[ztb] ] } 83 !! where ztb is a 2D array: the bottom ocean temperature and ahtb 84 !! is a time and space varying diffusive coefficient defined by: 85 !! ahbt = zahbp if grad(rho).grad(h) < 0 86 !! = 0. otherwise. 87 !! Note that grad(.) is the along bottom slope gradient. grad(rho) 88 !! is evaluated using the local density (i.e. referenced at the 89 !! local depth). Typical value of ahbt is 2000 m2/s (equivalent to 90 !! a downslope velocity of 20 cm/s if the condition for slope 91 !! convection is satified) 92 !! Add this before trend to the general trend tra of the 93 !! botton ocean tracer point: 94 !! tra = tra + difft 95 !! 96 !! ** Action : - update tra at the bottom level with the bottom 97 !! boundary layer trend 98 !! 99 !! References : 100 !! Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591. 101 !!---------------------------------------------------------------------- 102 USE oce, ONLY : ztrtrd => ua ! use ua as 3D workspace 103 !! 104 INTEGER, INTENT( in ) :: kt ! ocean time-step 105 INTEGER :: ji, jj, jn ! dummy loop indices 106 INTEGER :: ik, iku, ikv 107 INTEGER :: ii0, ii1, ij0, ij1 ! temporary integers 108 INTEGER :: iku1, iku2, ikv1, ikv2 ! temporary intergers 109 REAL(wp) :: ze3u, ze3v ! temporary scalars 110 REAL(wp) :: zbtr, ztra 111 #if ! defined key_off_tra 112 REAL(wp) :: zgdrho, zalbet, zsign, zt, zs, zh 113 REAL(wp), DIMENSION(jpi,jpj) :: zki, zkj 114 #endif 115 REAL(wp), DIMENSION(jpi,jpj) :: zkx, zky ! temporary workspace arrays 116 REAL(wp), DIMENSION(jpi,jpj) :: ztnb, zsnb, zdep 117 REAL(wp), DIMENSION(jpi,jpj) :: ztrb, zahu, zahv 118 52 !!---------------------------------------------------------------------- 53 INTEGER, INTENT( in ) :: kt ! ocean time-step 119 54 CHARACTER (len=22) :: charout 120 REAL(wp) :: fsalbt, pft, pfs, pfh ! statement function 121 !!---------------------------------------------------------------------- 122 ! ratio alpha/beta 123 ! ================ 124 ! fsalbt: ratio of thermal over saline expension coefficients 125 ! pft : potential temperature in degrees celcius 126 ! pfs : salinity anomaly (s-35) in psu 127 ! pfh : depth in meters 128 129 fsalbt( pft, pfs, pfh ) = & 130 ( ( ( -0.255019e-07 * pft + 0.298357e-05 ) * pft & 131 - 0.203814e-03 ) * pft & 132 + 0.170907e-01 ) * pft & 133 + 0.665157e-01 & 134 +(-0.678662e-05 * pfs - 0.846960e-04 * pft + 0.378110e-02 ) * pfs & 135 + ( ( - 0.302285e-13 * pfh & 136 - 0.251520e-11 * pfs & 137 + 0.512857e-12 * pft * pft ) * pfh & 138 - 0.164759e-06 * pfs & 139 +( 0.791325e-08 * pft - 0.933746e-06 ) * pft & 140 + 0.380374e-04 ) * pfh 55 REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE :: ztrtrd 141 56 !!---------------------------------------------------------------------- 142 57 143 144 IF( kt == nittrc000 ) CALL trc_bbl_init 145 146 147 ! 0. 2D fields of bottom temperature and salinity, and bottom slope 148 ! ----------------------------------------------------------------- 149 ! mbathy= number of w-level, minimum value=1 (cf dommsk.F) 150 151 # if defined key_vectopt_loop 152 jj = 1 153 DO ji = 1, jpij ! vector opt. (forced unrolling) 154 # else 155 DO jj = 1, jpj 156 DO ji = 1, jpi 157 # endif 158 ik = mbkt(ji,jj) ! index of the bottom ocean T-level 159 ztnb(ji,jj) = tn(ji,jj,ik) * tmask(ji,jj,1) ! masked now T and S at ocean bottom 160 zsnb(ji,jj) = sn(ji,jj,ik) * tmask(ji,jj,1) 161 zdep(ji,jj) = fsdept(ji,jj,ik) ! depth of the ocean bottom T-level 162 # if ! defined key_vectopt_loop 163 END DO 164 # endif 165 END DO 166 167 IF( ln_zps ) THEN ! partial steps correction 168 169 # if defined key_vectopt_loop 170 jj = 1 171 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 172 # else 173 DO jj = 1, jpjm1 174 DO ji = 1, jpim1 175 # endif 176 iku1 = MAX( mbathy(ji+1,jj )-1, 1 ) 177 iku2 = MAX( mbathy(ji ,jj )-1, 1 ) 178 ikv1 = MAX( mbathy(ji ,jj+1)-1, 1 ) 179 ikv2 = MAX( mbathy(ji ,jj )-1, 1 ) 180 ze3u = MIN( fse3u(ji,jj,iku1), fse3u(ji,jj,iku2) ) 181 ze3v = MIN( fse3v(ji,jj,ikv1), fse3v(ji,jj,ikv2) ) 182 zahu(ji,jj) = atrcbbl * e2u(ji,jj) * ze3u / e1u(ji,jj) * umask(ji,jj,1) 183 zahv(ji,jj) = atrcbbl * e1v(ji,jj) * ze3v / e2v(ji,jj) * vmask(ji,jj,1) 184 # if ! defined key_vectopt_loop 185 END DO 186 # endif 187 END DO 188 ELSE ! z-coordinate - full steps or s-coordinate 189 # if defined key_vectopt_loop 190 jj = 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) = atrcbbl * e2u(ji,jj) * fse3u(ji,jj,iku) / e1u(ji,jj) * umask(ji,jj,1) 199 zahv(ji,jj) = atrcbbl * e1v(ji,jj) * fse3v(ji,jj,ikv) / e2v(ji,jj) * vmask(ji,jj,1) 200 # if ! defined key_vectopt_loop 201 END DO 202 # endif 203 END DO 58 IF( .NOT. lk_offline ) THEN 59 CALL bbl( kt, 'TRC' ) ! Online coupling with dynamics : Computation of bbl coef and bbl transport 60 l_bbl = .FALSE. ! Offline coupling with dynamics : Read bbl coef and bbl transport from input files 204 61 ENDIF 205 62 206 #if defined key_off_tra 207 !!===================================================================== 208 !! I. OFFLINE VERSION OF DIFFUSIVE BBL 209 !!===================================================================== 210 211 ! 1. Criteria of additional bottom diffusivity : grad(rho).grad(h) < 0 212 ! -------------------------------------------------------------------- 213 214 ! Only used in the online version of diffusive BBL (see below) 215 216 ! 2. Additional second order diffusive trends 217 ! ------------------------------------------- 218 ! ! =========== 219 DO jn = 1, jptra ! tracer loop 220 ! ! =========== 221 222 IF( l_trdtrc ) ztrtrd(:,:,:) = tra(:,:,:,jn) 223 224 ! first derivative (gradient) 225 # if defined key_vectopt_loop 226 jj = 1 227 DO ji = 1, jpij ! vector opt. (forced unrolling) 228 # else 229 DO jj = 1, jpj 230 DO ji = 1, jpi 231 # endif 232 ik = mbkt(ji,jj) 233 ztrb(ji,jj) = trb(ji,jj,ik,jn) * tmask(ji,jj,1) 234 # if ! defined key_vectopt_loop 235 END DO 236 # endif 237 END DO 238 239 # if defined key_vectopt_loop 240 jj = 1 241 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 242 # else 243 DO jj = 1, jpjm1 244 DO ji = 1, jpim1 245 # endif 246 zkx(ji,jj) = bblx(ji,jj) * zahu(ji,jj) * ( ztrb(ji+1,jj) - ztrb(ji,jj) ) 247 zky(ji,jj) = bbly(ji,jj) * zahv(ji,jj) * ( ztrb(ji,jj+1) - ztrb(ji,jj) ) 248 # if ! defined key_vectopt_loop 249 END DO 250 # endif 251 END DO 252 253 #else 254 !!===================================================================== 255 !! II. ONLINE VERSION OF DIFFUSIVE BBL 256 !!===================================================================== 257 258 ! 1. Criteria of additional bottom diffusivity : grad(rho).grad(h) < 0 259 ! -------------------------------------------------------------------- 260 ! Sign of the local density gradient along the i- and j-slopes 261 ! multiplied by the slope of the ocean bottom 262 SELECT CASE ( nn_eos ) 263 264 CASE ( 0 ) ! Jackett and McDougall (1994) formulation 265 266 # if defined key_vectopt_loop 267 jj = 1 268 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 269 # else 270 DO jj = 1, jpjm1 271 DO ji = 1, jpim1 272 # endif 273 ! temperature, salinity anomalie and depth 274 zt = 0.5 * ( ztnb(ji,jj) + ztnb(ji+1,jj) ) 275 zs = 0.5 * ( zsnb(ji,jj) + zsnb(ji+1,jj) ) - 35.0 276 zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) ) 277 ! masked ratio alpha/beta 278 zalbet = fsalbt( zt, zs, zh )*umask(ji,jj,1) 279 ! local density gradient along i-bathymetric slope 280 zgdrho = zalbet * ( ztnb(ji+1,jj) - ztnb(ji,jj) ) & 281 - ( zsnb(ji+1,jj) - zsnb(ji,jj) ) 282 ! sign of local i-gradient of density multiplied by the i-slope 283 zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) ) 284 zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj) 285 # if ! defined key_vectopt_loop 286 END DO 287 # endif 288 END DO 289 290 # if defined key_vectopt_loop 291 jj = 1 292 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 293 # else 294 DO jj = 1, jpjm1 295 DO ji = 1, jpim1 296 # endif 297 ! temperature, salinity anomalie and depth 298 zt = 0.5 * ( ztnb(ji,jj+1) + ztnb(ji,jj) ) 299 zs = 0.5 * ( zsnb(ji,jj+1) + zsnb(ji,jj) ) - 35.0 300 zh = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) ) 301 ! masked ratio alpha/beta 302 zalbet = fsalbt( zt, zs, zh )*vmask(ji,jj,1) 303 ! local density gradient along j-bathymetric slope 304 zgdrho = zalbet * ( ztnb(ji,jj+1) - ztnb(ji,jj) ) & 305 - ( zsnb(ji,jj+1) - zsnb(ji,jj) ) 306 ! sign of local j-gradient of density multiplied by the j-slope 307 zsign = SIGN( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) ) 308 zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj) 309 # if ! defined key_vectopt_loop 310 END DO 311 # endif 312 END DO 313 314 CASE ( 1 ) ! Linear formulation function of temperature only 315 316 # if defined key_vectopt_loop 317 jj = 1 318 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 319 # else 320 DO jj = 1, jpjm1 321 DO ji = 1, jpim1 322 # endif 323 ! local density gradient along i-bathymetric slope 324 zgdrho = ( ztnb(ji+1,jj) - ztnb(ji,jj) ) 325 ! sign of local i-gradient of density multiplied by the i-slope 326 zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) ) 327 zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj) 328 # if ! defined key_vectopt_loop 329 END DO 330 # endif 331 END DO 332 333 # if defined key_vectopt_loop 334 jj = 1 335 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 336 # else 337 DO jj = 1, jpjm1 338 DO ji = 1, jpim1 339 # endif 340 ! local density gradient along j-bathymetric slope 341 zgdrho = ( ztnb(ji,jj+1) - ztnb(ji,jj) ) 342 ! sign of local j-gradient of density multiplied by the j-slope 343 zsign = SIGN( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) ) 344 zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj) 345 346 # if ! defined key_vectopt_loop 347 END DO 348 # endif 349 END DO 350 351 CASE ( 2 ) ! Linear formulation function of temperature and salinity 352 353 DO jj = 1, jpjm1 354 DO ji = 1, fs_jpim1 ! vector opt. 355 ! local density gradient along i-bathymetric slope 356 zgdrho = - ( rn_beta*( zsnb(ji+1,jj) - zsnb(ji,jj) ) & 357 - rn_alpha*( ztnb(ji+1,jj) - ztnb(ji,jj) ) ) 358 ! sign of local i-gradient of density multiplied by the i-slope 359 zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) ) 360 zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj) 361 END DO 362 END DO 363 364 DO jj = 1, jpjm1 365 DO ji = 1, fs_jpim1 ! vector opt. 366 ! local density gradient along j-bathymetric slope 367 zgdrho = - ( rn_beta*( zsnb(ji,jj+1) - zsnb(ji,jj) ) & 368 - rn_alpha*( ztnb(ji,jj+1) - ztnb(ji,jj) ) ) 369 ! sign of local j-gradient of density multiplied by the j-slope 370 zsign = sign( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) ) 371 zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj) 372 END DO 373 END DO 374 375 CASE DEFAULT 376 377 WRITE(ctmp1,*) ' bad flag value for nn_eos = ', nn_eos 378 CALL ctl_stop( ctmp1 ) 379 380 END SELECT 381 382 ! 2. Additional second order diffusive trends 383 ! ------------------------------------------- 384 ! ! =========== 385 DO jn = 1, jptra ! tracer loop 386 ! ! =========== 387 IF( l_trdtrc ) ztrtrd(:,:,:) = tra(:,:,:,jn) 388 389 ! first derivative (gradient) 390 # if defined key_vectopt_loop 391 jj = 1 392 DO ji = 1, jpij ! vector opt. (forced unrolling) 393 # else 394 DO jj = 1, jpj 395 DO ji = 1, jpi 396 # endif 397 ik = mbkt(ji,jj) 398 ztrb(ji,jj) = trb(ji,jj,ik,jn) * tmask(ji,jj,1) 399 # if ! defined key_vectopt_loop 400 END DO 401 # endif 402 END DO 403 # if defined key_vectopt_loop 404 jj = 1 405 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 406 # else 407 DO jj = 1, jpjm1 408 DO ji = 1, jpim1 409 # endif 410 zkx(ji,jj) = zki(ji,jj) * ( ztrb(ji+1,jj) - ztrb(ji,jj) ) 411 zky(ji,jj) = zkj(ji,jj) * ( ztrb(ji,jj+1) - ztrb(ji,jj) ) 412 # if ! defined key_vectopt_loop 413 END DO 414 # endif 415 END DO 416 #endif 417 418 !!===================================================================== 419 !! III. COMMON CODE FOR OFFLINE/ONLINE VERSIONS OF DIFFUSIVE BBL 420 !!===================================================================== 421 422 IF( cp_cfg == "orca" ) THEN 423 424 SELECT CASE ( jp_cfg ) 425 ! ! ======================= 426 CASE ( 2 ) ! ORCA_R2 configuration 427 ! ! ======================= 428 ! Gibraltar enhancement of BBL 429 ij0 = 102 ; ij1 = 102 430 ii0 = 139 ; ii1 = 140 431 zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) 432 zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) 433 434 ! Red Sea enhancement of BBL 435 ij0 = 88 ; ij1 = 88 436 ii0 = 161 ; ii1 = 162 437 zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 10.e0 * zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) 438 zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 10.e0 * zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) 439 440 ! ! ======================= 441 CASE ( 4 ) ! ORCA_R4 configuration 442 ! ! ======================= 443 ! Gibraltar enhancement of BBL 444 ij0 = 52 ; ij1 = 52 445 ii0 = 70 ; ii1 = 71 446 zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zkx( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) 447 zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 * zky( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) 448 449 END SELECT 450 451 ENDIF 452 453 ! second derivative (divergence) and add to the general tracer trend 454 # if defined key_vectopt_loop 455 jj = 1 456 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) 457 # else 458 DO jj = 2, jpjm1 459 DO ji = 2, jpim1 460 # endif 461 ik = MAX( mbathy(ji,jj)-1, 1 ) 462 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,ik) ) 463 ztra = ( zkx(ji,jj) - zkx(ji-1,jj ) & 464 & + zky(ji,jj) - zky(ji ,jj-1) ) * zbtr 465 tra(ji,jj,ik,jn) = tra(ji,jj,ik,jn) + ztra 466 # if ! defined key_vectopt_loop 467 END DO 468 # endif 469 END DO 470 471 ! save the trends for diagnostic 472 IF( l_trdtrc ) THEN 473 ztrtrd(:,:,:) = tra(:,:,:,jn) - ztrtrd(:,:,:) 474 IF (luttrd(jn)) CALL trd_mod_trc( ztrtrd, jn, jptrc_trd_bbl, kt ) 475 END IF 476 ! ! =========== 477 END DO ! tracer loop 478 ! ! =========== 479 480 IF( ln_ctl ) THEN ! print mean trends (used for debugging) 481 WRITE(charout, FMT="('bbl - dif')") 482 CALL prt_ctl_trc_info(charout) 483 CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd') 63 IF( l_trdtrc ) THEN 64 ALLOCATE( ztrtrd(jpi,jpj,jpk,jptra) ) ! temporary save of trends 65 ztrtrd(:,:,:,:) = tra(:,:,:,:) 484 66 ENDIF 485 67 486 END SUBROUTINE trc_bbl_dif 68 !* Diffusive bbl : 69 IF( nn_bbl_ldf == 1 ) THEN 70 ! 71 CALL tra_bbl_dif( trb, tra, jptra ) 72 IF( ln_ctl ) THEN 73 WRITE(charout, FMT="(' bbl_dif')") ; CALL prt_ctl_trc_info(charout) 74 CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 75 ENDIF 76 ! 77 END IF 487 78 488 # if defined key_trcbbl_adv 489 !!---------------------------------------------------------------------- 490 !! 'key_trcbbl_adv' advective bottom boundary layer 491 !!---------------------------------------------------------------------- 492 # include "trcbbl_adv.h90" 493 # else 494 !!---------------------------------------------------------------------- 495 !! Default option : NO advective bottom boundary layer 496 !!---------------------------------------------------------------------- 497 SUBROUTINE trc_bbl_adv (kt ) ! Empty routine 498 INTEGER, INTENT(in) :: kt 499 WRITE(*,*) 'trc_bbl_adv: You should not have seen this print! error?', kt 500 END SUBROUTINE trc_bbl_adv 501 # endif 79 !* Advective bbl : bbl upstream advective trends added to the tracer trends 80 IF( nn_bbl_adv /= 0 ) THEN 81 ! 82 CALL tra_bbl_adv( trb, tra, jptra ) 83 IF( ln_ctl ) THEN 84 WRITE(charout, FMT="(' bbl_adv')") ; CALL prt_ctl_trc_info(charout) 85 CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 86 ENDIF 87 ! 88 END IF 502 89 503 SUBROUTINE trc_bbl_init 504 !!---------------------------------------------------------------------- 505 !! *** ROUTINE trc_bbl_init *** 506 !! 507 !! ** Purpose : Initialization for the bottom boundary layer scheme. 508 !!---------------------------------------------------------------------- 509 INTEGER :: ji, jj 510 REAL(wp), DIMENSION(jpi,jpj) :: zmbk 511 !!---------------------------------------------------------------------- 512 513 DO jj = 1, jpj 514 DO ji = 1, jpi 515 mbkt(ji,jj) = MAX( mbathy(ji,jj) - 1, 1 ) ! vertical index of the bottom ocean T-level 516 END DO 517 END DO 518 519 DO jj = 1, jpjm1 520 DO ji = 1, jpim1 521 mbku(ji,jj) = MAX( MIN( mbathy(ji+1,jj ), mbathy(ji,jj) ) - 1, 1 ) 522 mbkv(ji,jj) = MAX( MIN( mbathy(ji ,jj+1), mbathy(ji,jj) ) - 1, 1 ) 523 END DO 524 END DO 525 526 zmbk(:,:) = FLOAT( mbku (:,:) ) 527 CALL lbc_lnk(zmbk,'U',1.) 528 mbku (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 529 530 zmbk(:,:) = FLOAT( mbkv (:,:) ) 531 CALL lbc_lnk(zmbk,'V',1.) 532 mbkv (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 533 534 # if defined key_trcbbl_adv 535 w_trc_bbl(:,:,:) = 0.e0 ! initialisation of w_trc_bbl to zero 536 # endif 537 538 END SUBROUTINE trc_bbl_init 90 IF( l_trdtrc ) THEN ! save the horizontal diffusive trends for further diagnostics 91 DO jn = 1, jptra 92 ztrtrd(:,:,:,jn) = tra(:,:,:,jn) - ztrtrd(:,:,:,jn) 93 CALL trd_tra( kt, 'TRC', jn, jptra_trd_ldf, ztrtrd(:,:,:,jn) ) 94 END DO 95 DEALLOCATE( ztrtrd ) 96 ENDIF 97 ! 98 END SUBROUTINE trc_bbl 539 99 540 100 #else … … 542 102 !! Dummy module : No bottom boundary layer scheme 543 103 !!---------------------------------------------------------------------- 544 LOGICAL, PUBLIC, PARAMETER :: lk_trcbbl_dif = .FALSE. !: diff bbl flag545 LOGICAL, PUBLIC, PARAMETER :: lk_trcbbl_adv = .FALSE. !: adv bbl flag546 104 CONTAINS 547 SUBROUTINE trc_bbl_dif (kt ) ! Empty routine 548 INTEGER, INTENT(in) :: kt 549 WRITE(*,*) 'trc_bbl_dif: You should not have seen this print! error?', kt 550 END SUBROUTINE trc_bbl_dif 551 SUBROUTINE trc_bbl_adv (kt ) ! Empty routine 552 INTEGER, INTENT(in) :: kt 553 WRITE(*,*) 'trc_bbl_adv: You should not have seen this print! error?', kt 554 END SUBROUTINE trc_bbl_adv 105 SUBROUTINE trc_bbl( kt ) ! Empty routine 106 WRITE(*,*) 'tra_bbl: You should not have seen this print! error?', kt 107 END SUBROUTINE trc_bbl 555 108 #endif 556 109
Note: See TracChangeset
for help on using the changeset viewer.