Changeset 5965 for branches/2014/dev_r4650_UKMO14.5_SST_BIAS_CORRECTION/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfbfr.F90
- Timestamp:
- 2015-12-01T16:35:30+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4650_UKMO14.5_SST_BIAS_CORRECTION/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfbfr.F90
r4624 r5965 41 41 REAL(wp), PUBLIC :: rn_bfrien ! local factor to enhance coefficient bfri (PUBLIC for TAM) 42 42 LOGICAL , PUBLIC :: ln_bfr2d ! logical switch for 2D enhancement (PUBLIC for TAM) 43 REAL(wp), PUBLIC :: rn_tfri1 ! top drag coefficient (linear case) (PUBLIC for TAM) 44 REAL(wp), PUBLIC :: rn_tfri2 ! top drag coefficient (non linear case) (PUBLIC for TAM) 45 REAL(wp), PUBLIC :: rn_tfri2_max ! Maximum top drag coefficient (non linear case and ln_loglayer=T) (PUBLIC for TAM) 46 REAL(wp), PUBLIC :: rn_tfeb2 ! background top turbulent kinetic energy [m2/s2] (PUBLIC for TAM) 47 REAL(wp), PUBLIC :: rn_tfrien ! local factor to enhance coefficient tfri (PUBLIC for TAM) 48 LOGICAL , PUBLIC :: ln_tfr2d ! logical switch for 2D enhancement (PUBLIC for TAM) 49 43 50 LOGICAL , PUBLIC :: ln_loglayer ! switch for log layer bfr coeff. (PUBLIC for TAM) 44 51 REAL(wp), PUBLIC :: rn_bfrz0 ! bottom roughness for loglayer bfr coeff (PUBLIC for TAM) 52 REAL(wp), PUBLIC :: rn_tfrz0 ! bottom roughness for loglayer bfr coeff (PUBLIC for TAM) 45 53 LOGICAL , PUBLIC :: ln_bfrimp ! logical switch for implicit bottom friction 46 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: bfrcoef2d ! 2D bottomdrag coefficient (PUBLIC for TAM)54 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: bfrcoef2d, tfrcoef2d ! 2D bottom/top drag coefficient (PUBLIC for TAM) 47 55 48 56 !! * Substitutions … … 60 68 !! *** FUNCTION zdf_bfr_alloc *** 61 69 !!---------------------------------------------------------------------- 62 ALLOCATE( bfrcoef2d(jpi,jpj), STAT=zdf_bfr_alloc )70 ALLOCATE( bfrcoef2d(jpi,jpj), tfrcoef2d(jpi,jpj), STAT=zdf_bfr_alloc ) 63 71 ! 64 72 IF( lk_mpp ) CALL mpp_sum ( zdf_bfr_alloc ) … … 88 96 INTEGER :: ikbt, ikbu, ikbv ! local integers 89 97 REAL(wp) :: zvu, zuv, zecu, zecv, ztmp ! temporary scalars 90 REAL(wp), POINTER, DIMENSION(:,:) :: zbfrt 98 REAL(wp), POINTER, DIMENSION(:,:) :: zbfrt, ztfrt 91 99 !!---------------------------------------------------------------------- 92 100 ! … … 101 109 IF( nn_bfr == 2 ) THEN ! quadratic bottom friction only 102 110 ! 103 CALL wrk_alloc( jpi, jpj, zbfrt )111 CALL wrk_alloc( jpi, jpj, zbfrt, ztfrt ) 104 112 105 113 IF ( ln_loglayer.AND.lk_vvl ) THEN ! "log layer" bottom friction coefficient 106 114 107 # if defined key_vectopt_loop108 DO jj = 1, 1109 !CDIR NOVERRCHK110 DO ji = 1, jpij ! vector opt. (forced unrolling)111 # else112 !CDIR NOVERRCHK113 115 DO jj = 1, jpj 114 !CDIR NOVERRCHK115 116 DO ji = 1, jpi 116 # endif117 117 ikbt = mbkt(ji,jj) 118 ! JC: possible WAD implementation should modify line below if layers vanish118 !! JC: possible WAD implementation should modify line below if layers vanish 119 119 ztmp = tmask(ji,jj,ikbt) * ( vkarmn / LOG( 0.5_wp * fse3t_n(ji,jj,ikbt) / rn_bfrz0 ))**2._wp 120 120 zbfrt(ji,jj) = MAX(bfrcoef2d(ji,jj), ztmp) … … 122 122 END DO 123 123 END DO 124 ! (ISF) 125 IF ( ln_isfcav ) THEN 126 DO jj = 1, jpj 127 DO ji = 1, jpi 128 ikbt = mikt(ji,jj) 129 ! JC: possible WAD implementation should modify line below if layers vanish 130 ztmp = (1-tmask(ji,jj,1)) * ( vkarmn / LOG( 0.5_wp * fse3t_n(ji,jj,ikbt) / rn_bfrz0 ))**2._wp 131 ztfrt(ji,jj) = MAX(tfrcoef2d(ji,jj), ztmp) 132 ztfrt(ji,jj) = MIN(ztfrt(ji,jj), rn_tfri2_max) 133 END DO 134 END DO 135 END IF 124 136 ! 125 137 ELSE 126 138 zbfrt(:,:) = bfrcoef2d(:,:) 127 ENDIF 128 129 # if defined key_vectopt_loop 130 DO jj = 1, 1 131 !CDIR NOVERRCHK 132 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) 133 # else 134 !CDIR NOVERRCHK 139 ztfrt(:,:) = tfrcoef2d(:,:) 140 ENDIF 141 135 142 DO jj = 2, jpjm1 136 !CDIR NOVERRCHK137 143 DO ji = 2, jpim1 138 # endif139 144 ikbu = mbku(ji,jj) ! ocean bottom level at u- and v-points 140 145 ikbv = mbkv(ji,jj) ! (deepest ocean u- and v-points) … … 150 155 bfrua(ji,jj) = - 0.5_wp * ( zbfrt(ji,jj) + zbfrt(ji+1,jj ) ) * zecu 151 156 bfrva(ji,jj) = - 0.5_wp * ( zbfrt(ji,jj) + zbfrt(ji ,jj+1) ) * zecv 157 ! 158 ! in case of 2 cell water column, we assume each cell feels the top and bottom friction 159 IF ( ln_isfcav ) THEN 160 IF ( miku(ji,jj) + 1 .GE. mbku(ji,jj) ) THEN 161 bfrua(ji,jj) = - 0.5_wp * ( ( zbfrt(ji,jj) + zbfrt(ji+1,jj ) ) & 162 & + ( ztfrt(ji,jj) + ztfrt(ji+1,jj ) ) ) & 163 & * zecu * (1._wp - umask(ji,jj,1)) 164 END IF 165 IF ( mikv(ji,jj) + 1 .GE. mbkv(ji,jj) ) THEN 166 bfrva(ji,jj) = - 0.5_wp * ( ( zbfrt(ji,jj) + zbfrt(ji ,jj+1) ) & 167 & + ( ztfrt(ji,jj) + ztfrt(ji ,jj+1) ) ) & 168 & * zecv * (1._wp - vmask(ji,jj,1)) 169 END IF 170 END IF 152 171 END DO 153 172 END DO 154 155 !156 173 CALL lbc_lnk( bfrua, 'U', 1. ) ; CALL lbc_lnk( bfrva, 'V', 1. ) ! Lateral boundary condition 174 175 IF ( ln_isfcav ) THEN 176 DO jj = 2, jpjm1 177 DO ji = 2, jpim1 178 ! (ISF) ======================================================================== 179 ikbu = miku(ji,jj) ! ocean top level at u- and v-points 180 ikbv = mikv(ji,jj) ! (1st wet ocean u- and v-points) 181 ! 182 zvu = 0.25 * ( vn(ji,jj ,ikbu) + vn(ji+1,jj ,ikbu) & 183 & + vn(ji,jj-1,ikbu) + vn(ji+1,jj-1,ikbu) ) 184 zuv = 0.25 * ( un(ji,jj ,ikbv) + un(ji-1,jj ,ikbv) & 185 & + un(ji,jj+1,ikbv) + un(ji-1,jj+1,ikbv) ) 186 ! 187 zecu = SQRT( un(ji,jj,ikbu) * un(ji,jj,ikbu) + zvu*zvu + rn_tfeb2 ) 188 zecv = SQRT( vn(ji,jj,ikbv) * vn(ji,jj,ikbv) + zuv*zuv + rn_tfeb2 ) 189 ! 190 tfrua(ji,jj) = - 0.5_wp * ( ztfrt(ji,jj) + ztfrt(ji+1,jj ) ) * zecu * (1._wp - umask(ji,jj,1)) 191 tfrva(ji,jj) = - 0.5_wp * ( ztfrt(ji,jj) + ztfrt(ji ,jj+1) ) * zecv * (1._wp - vmask(ji,jj,1)) 192 ! (ISF) END ==================================================================== 193 ! in case of 2 cell water column, we assume each cell feels the top and bottom friction 194 IF ( miku(ji,jj) + 1 .GE. mbku(ji,jj) ) THEN 195 tfrua(ji,jj) = - 0.5_wp * ( ( ztfrt(ji,jj) + ztfrt(ji+1,jj ) ) & 196 & + ( zbfrt(ji,jj) + zbfrt(ji+1,jj ) ) ) & 197 & * zecu * (1._wp - umask(ji,jj,1)) 198 END IF 199 IF ( mikv(ji,jj) + 1 .GE. mbkv(ji,jj) ) THEN 200 tfrva(ji,jj) = - 0.5_wp * ( ( ztfrt(ji,jj) + ztfrt(ji ,jj+1) ) & 201 & + ( zbfrt(ji,jj) + zbfrt(ji ,jj+1) ) ) & 202 & * zecv * (1._wp - vmask(ji,jj,1)) 203 END IF 204 END DO 205 END DO 206 CALL lbc_lnk( tfrua, 'U', 1. ) ; CALL lbc_lnk( tfrva, 'V', 1. ) ! Lateral boundary condition 207 END IF 208 ! 157 209 ! 158 210 IF(ln_ctl) CALL prt_ctl( tab2d_1=bfrua, clinfo1=' bfr - u: ', mask1=umask, & 159 211 & tab2d_2=bfrva, clinfo2= ' v: ', mask2=vmask,ovlap=1 ) 160 CALL wrk_dealloc( jpi,jpj,zbfrt )212 CALL wrk_dealloc( jpi,jpj,zbfrt, ztfrt ) 161 213 ENDIF 162 214 ! … … 183 235 INTEGER :: ios 184 236 REAL(wp) :: zminbfr, zmaxbfr ! temporary scalars 237 REAL(wp) :: zmintfr, zmaxtfr ! temporary scalars 185 238 REAL(wp) :: ztmp, zfru, zfrv ! - - 186 239 !! 187 240 NAMELIST/nambfr/ nn_bfr, rn_bfri1, rn_bfri2, rn_bfri2_max, rn_bfeb2, rn_bfrz0, ln_bfr2d, & 188 & rn_bfrien, ln_bfrimp, ln_loglayer 241 & rn_tfri1, rn_tfri2, rn_tfri2_max, rn_tfeb2, rn_tfrz0, ln_tfr2d, & 242 & rn_bfrien, rn_tfrien, ln_bfrimp, ln_loglayer 189 243 !!---------------------------------------------------------------------- 190 244 ! … … 215 269 bfrua(:,:) = 0.e0 216 270 bfrva(:,:) = 0.e0 271 tfrua(:,:) = 0.e0 272 tfrva(:,:) = 0.e0 217 273 ! 218 274 CASE( 1 ) 219 275 IF(lwp) WRITE(numout,*) ' linear botton friction' 220 IF(lwp) WRITE(numout,*) ' friction coef. rn_bfri1 = ', rn_bfri1276 IF(lwp) WRITE(numout,*) ' bottom friction coef. rn_bfri1 = ', rn_bfri1 221 277 IF( ln_bfr2d ) THEN 222 278 IF(lwp) WRITE(numout,*) ' read coef. enhancement distribution from file ln_bfr2d = ', ln_bfr2d 223 279 IF(lwp) WRITE(numout,*) ' coef rn_bfri2 enhancement factor rn_bfrien = ',rn_bfrien 224 280 ENDIF 281 IF ( ln_isfcav ) THEN 282 IF(lwp) WRITE(numout,*) ' top friction coef. rn_bfri1 = ', rn_tfri1 283 IF( ln_tfr2d ) THEN 284 IF(lwp) WRITE(numout,*) ' read coef. enhancement distribution from file ln_tfr2d = ', ln_tfr2d 285 IF(lwp) WRITE(numout,*) ' coef rn_tfri2 enhancement factor rn_tfrien = ',rn_tfrien 286 ENDIF 287 END IF 225 288 ! 226 289 IF(ln_bfr2d) THEN … … 236 299 bfrua(:,:) = - bfrcoef2d(:,:) 237 300 bfrva(:,:) = - bfrcoef2d(:,:) 301 ! 302 IF ( ln_isfcav ) THEN 303 IF(ln_tfr2d) THEN 304 ! tfr_coef is a coefficient in [0,1] giving the mask where to apply the bfr enhancement 305 CALL iom_open('tfr_coef.nc',inum) 306 CALL iom_get (inum, jpdom_data, 'tfr_coef',tfrcoef2d,1) ! tfrcoef2d is used as tmp array 307 CALL iom_close(inum) 308 tfrcoef2d(:,:) = rn_tfri1 * ( 1 + rn_tfrien * tfrcoef2d(:,:) ) 309 ELSE 310 tfrcoef2d(:,:) = rn_tfri1 ! initialize tfrcoef2d to the namelist variable 311 ENDIF 312 ! 313 tfrua(:,:) = - tfrcoef2d(:,:) 314 tfrva(:,:) = - tfrcoef2d(:,:) 315 END IF 238 316 ! 239 317 CASE( 2 ) … … 252 330 IF(lwp) WRITE(numout,*) ' coef rn_bfri2 enhancement factor rn_bfrien = ',rn_bfrien 253 331 ENDIF 332 IF ( ln_isfcav ) THEN 333 IF(lwp) WRITE(numout,*) ' quadratic top friction' 334 IF(lwp) WRITE(numout,*) ' friction coef. rn_tfri2 = ', rn_tfri2 335 IF(lwp) WRITE(numout,*) ' Max. coef. (log case) rn_tfri2_max = ', rn_tfri2_max 336 IF(lwp) WRITE(numout,*) ' background tke rn_tfeb2 = ', rn_tfeb2 337 IF(lwp) WRITE(numout,*) ' log formulation ln_tfr2d = ', ln_loglayer 338 IF(lwp) WRITE(numout,*) ' top roughness rn_tfrz0 [m] = ', rn_tfrz0 339 IF( rn_tfrz0<=0.e0 ) THEN 340 WRITE(ctmp1,*) ' top roughness must be strictly positive' 341 CALL ctl_stop( ctmp1 ) 342 ENDIF 343 IF( ln_tfr2d ) THEN 344 IF(lwp) WRITE(numout,*) ' read coef. enhancement distribution from file ln_tfr2d = ', ln_tfr2d 345 IF(lwp) WRITE(numout,*) ' coef rn_tfri2 enhancement factor rn_tfrien = ',rn_tfrien 346 ENDIF 347 END IF 254 348 ! 255 349 IF(ln_bfr2d) THEN … … 263 357 bfrcoef2d(:,:) = rn_bfri2 ! initialize bfrcoef2d to the namelist variable 264 358 ENDIF 359 360 IF ( ln_isfcav ) THEN 361 IF(ln_tfr2d) THEN 362 ! tfr_coef is a coefficient in [0,1] giving the mask where to apply the bfr enhancement 363 CALL iom_open('tfr_coef.nc',inum) 364 CALL iom_get (inum, jpdom_data, 'tfr_coef',tfrcoef2d,1) ! tfrcoef2d is used as tmp array 365 CALL iom_close(inum) 366 ! 367 tfrcoef2d(:,:) = rn_tfri2 * ( 1 + rn_tfrien * tfrcoef2d(:,:) ) 368 ELSE 369 tfrcoef2d(:,:) = rn_tfri2 ! initialize tfrcoef2d to the namelist variable 370 ENDIF 371 END IF 265 372 ! 266 373 IF ( ln_loglayer.AND.(.NOT.lk_vvl) ) THEN ! set "log layer" bottom friction once for all 267 # if defined key_vectopt_loop268 DO jj = 1, 1269 !CDIR NOVERRCHK270 DO ji = 1, jpij ! vector opt. (forced unrolling)271 # else272 !CDIR NOVERRCHK273 374 DO jj = 1, jpj 274 !CDIR NOVERRCHK275 375 DO ji = 1, jpi 276 # endif277 376 ikbt = mbkt(ji,jj) 278 377 ztmp = tmask(ji,jj,ikbt) * ( vkarmn / LOG( 0.5_wp * fse3t_n(ji,jj,ikbt) / rn_bfrz0 ))**2._wp … … 281 380 END DO 282 381 END DO 382 IF ( ln_isfcav ) THEN 383 DO jj = 1, jpj 384 DO ji = 1, jpi 385 ikbt = mikt(ji,jj) 386 ztmp = tmask(ji,jj,ikbt) * ( vkarmn / LOG( 0.5_wp * fse3t_n(ji,jj,ikbt) / rn_tfrz0 ))**2._wp 387 tfrcoef2d(ji,jj) = MAX(tfrcoef2d(ji,jj), ztmp) 388 tfrcoef2d(ji,jj) = MIN(tfrcoef2d(ji,jj), rn_tfri2_max) 389 END DO 390 END DO 391 END IF 283 392 ENDIF 284 393 ! … … 308 417 zminbfr = 1.e10_wp ! initialise tracker for minimum of bottom friction coefficient 309 418 zmaxbfr = -1.e10_wp ! initialise tracker for maximum of bottom friction coefficient 310 ! 311 # if defined key_vectopt_loop 312 DO jj = 1, 1 313 !CDIR NOVERRCHK 314 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) 315 # else 316 !CDIR NOVERRCHK 419 zmintfr = 1.e10_wp ! initialise tracker for minimum of bottom friction coefficient 420 zmaxtfr = -1.e10_wp ! initialise tracker for maximum of bottom friction coefficient 421 ! 317 422 DO jj = 2, jpjm1 318 !CDIR NOVERRCHK319 423 DO ji = 2, jpim1 320 # endif321 424 ikbu = mbku(ji,jj) ! deepest ocean level at u- and v-points 322 425 ikbv = mbkv(ji,jj) … … 339 442 zminbfr = MIN( zminbfr, MIN( zfru, ABS( bfrcoef2d(ji,jj) ) ) ) 340 443 zmaxbfr = MAX( zmaxbfr, MIN( zfrv, ABS( bfrcoef2d(ji,jj) ) ) ) 444 ! (ISF) 445 IF ( ln_isfcav ) THEN 446 ikbu = miku(ji,jj) ! 1st wet ocean level at u- and v-points 447 ikbv = mikv(ji,jj) 448 zfru = 0.5 * fse3u(ji,jj,ikbu) / rdt 449 zfrv = 0.5 * fse3v(ji,jj,ikbv) / rdt 450 IF( ABS( tfrcoef2d(ji,jj) ) > zfru ) THEN 451 IF( ln_ctl ) THEN 452 WRITE(numout,*) 'TFR ', narea, nimpp+ji, njmpp+jj, ikbu 453 WRITE(numout,*) 'TFR ', ABS( tfrcoef2d(ji,jj) ), zfru 454 ENDIF 455 ictu = ictu + 1 456 ENDIF 457 IF( ABS( tfrcoef2d(ji,jj) ) > zfrv ) THEN 458 IF( ln_ctl ) THEN 459 WRITE(numout,*) 'TFR ', narea, nimpp+ji, njmpp+jj, ikbv 460 WRITE(numout,*) 'TFR ', tfrcoef2d(ji,jj), zfrv 461 ENDIF 462 ictv = ictv + 1 463 ENDIF 464 zmintfr = MIN( zmintfr, MIN( zfru, ABS( tfrcoef2d(ji,jj) ) ) ) 465 zmaxtfr = MAX( zmaxtfr, MIN( zfrv, ABS( tfrcoef2d(ji,jj) ) ) ) 466 END IF 467 ! END ISF 341 468 END DO 342 469 END DO … … 346 473 CALL mpp_min( zminbfr ) 347 474 CALL mpp_max( zmaxbfr ) 475 IF ( ln_isfcav) CALL mpp_min( zmintfr ) 476 IF ( ln_isfcav) CALL mpp_max( zmaxtfr ) 348 477 ENDIF 349 478 IF( .NOT.ln_bfrimp) THEN 350 479 IF( lwp .AND. ictu + ictv > 0 ) THEN 351 WRITE(numout,*) ' Bottom friction stability check failed at ', ictu, ' U-points '352 WRITE(numout,*) ' Bottom friction stability check failed at ', ictv, ' V-points '480 WRITE(numout,*) ' Bottom/Top friction stability check failed at ', ictu, ' U-points ' 481 WRITE(numout,*) ' Bottom/Top friction stability check failed at ', ictv, ' V-points ' 353 482 WRITE(numout,*) ' Bottom friction coefficient now ranges from: ', zminbfr, ' to ', zmaxbfr 354 WRITE(numout,*) ' Bottom friction coefficient will be reduced where necessary' 483 IF ( ln_isfcav ) WRITE(numout,*) ' Top friction coefficient now ranges from: ', zmintfr, ' to ', zmaxtfr 484 WRITE(numout,*) ' Bottom/Top friction coefficient will be reduced where necessary' 355 485 ENDIF 356 486 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.
