- Timestamp:
- 2020-09-24T20:38:10+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r13383_HPC-02_Daley_Tiling/src/OCE/TRA/traadv_ubs.F90
r13295 r13516 14 14 USE oce ! ocean dynamics and active tracers 15 15 USE dom_oce ! ocean space and time domain 16 ! TEMP: This change not necessary after trd_tra is tiled 17 USE domain, ONLY : dom_tile 16 18 USE trc_oce ! share passive tracers/Ocean variables 17 19 USE trd_oce ! trends: ocean variables … … 92 94 INTEGER , INTENT(in ) :: kn_ubs_v ! number of tracers 93 95 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 96 ! TEMP: This can be ST_2D(nn_hls) after trd_tra is tiled 94 97 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pU, pV, pW ! 3 ocean volume transport components 95 98 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! tracers and RHS of tracer equation 96 99 ! 100 ! TEMP: This change not necessary after trd_tra is tiled 101 INTEGER :: itile 97 102 INTEGER :: ji, jj, jk, jn ! dummy loop indices 98 103 REAL(wp) :: ztra, zbtr, zcoef ! local scalars 99 104 REAL(wp) :: zfp_ui, zfm_ui, zcenut, ztak, zfp_wk, zfm_wk ! - - 100 105 REAL(wp) :: zfp_vj, zfm_vj, zcenvt, zeeu, zeev, z_hdivn ! - - 101 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztu, ztv, zltu, zltv, zti, ztw ! 3D workspace 102 !!---------------------------------------------------------------------- 103 ! 104 IF( kt == kit000 ) THEN 105 IF(lwp) WRITE(numout,*) 106 IF(lwp) WRITE(numout,*) 'tra_adv_ubs : horizontal UBS advection scheme on ', cdtype 107 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 106 REAL(wp), DIMENSION(ST_2D(nn_hls),jpk) :: ztu, ztv, zltu, zltv, zti, ztw ! 3D workspace 107 ! TEMP: This change not necessary after trd_tra is tiled 108 REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE :: ztrdx, ztrdy, ztrdz 109 !!---------------------------------------------------------------------- 110 ! TEMP: This change not necessary after trd_tra is tiled 111 itile = ntile 112 ! 113 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 114 IF( kt == kit000 ) THEN 115 IF(lwp) WRITE(numout,*) 116 IF(lwp) WRITE(numout,*) 'tra_adv_ubs : horizontal UBS advection scheme on ', cdtype 117 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 118 ENDIF 119 ! 120 l_trd = .FALSE. 121 l_hst = .FALSE. 122 l_ptr = .FALSE. 123 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 124 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE. 125 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 126 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 127 128 ! TEMP: This can be ST_2D(nn_hls) after trd_tra is tiled 129 IF( kt == kit000 .AND. l_trd ) THEN 130 ALLOCATE( ztrdx(jpi,jpj,jpk), ztrdy(jpi,jpj,jpk), ztrdz(jpi,jpj,jpk) ) 131 ENDIF 108 132 ENDIF 109 !110 l_trd = .FALSE.111 l_hst = .FALSE.112 l_ptr = .FALSE.113 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.114 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE.115 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &116 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE.117 133 ! 118 134 ztw (:,:, 1 ) = 0._wp ! surface & bottom value : set to zero for all tracers … … 153 169 END_3D 154 170 ! 155 zltu(:,:,:) = pt(:,:,:,jn,Krhs) ! store the initial trends before its update 171 DO_3D( 1, 1, 1, 1, 1, jpk ) 172 zltu(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) ! store the initial trends before its update 173 END_3D 156 174 ! 157 175 DO jk = 1, jpkm1 !== add the horizontal advective trend ==! … … 165 183 END DO 166 184 ! 167 zltu(:,:,:) = pt(:,:,:,jn,Krhs) - zltu(:,:,:) ! Horizontal advective trend used in vertical 2nd order FCT case 168 ! ! and/or in trend diagnostic (l_trd=T) 169 ! 185 DO_3D( 1, 1, 1, 1, 1, jpk ) 186 zltu(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) - zltu(ji,jj,jk) ! Horizontal advective trend used in vertical 2nd order FCT case 187 END_3D ! and/or in trend diagnostic (l_trd=T) 188 ! 189 ! TEMP: These changes not necessary after trd_tra is tiled 170 190 IF( l_trd ) THEN ! trend diagnostics 171 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztu, pU, pt(:,:,:,jn,Kmm) ) 172 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztv, pV, pt(:,:,:,jn,Kmm) ) 191 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 192 ztrdx(ji,jj,jk) = ztu(ji,jj,jk) 193 ztrdy(ji,jj,jk) = ztv(ji,jj,jk) 194 END_3D 195 196 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 197 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 198 199 ! TODO: TO BE TILED- trd_tra 200 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztrdx, pU, pt(:,:,:,jn,Kmm) ) 201 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztrdy, pV, pt(:,:,:,jn,Kmm) ) 202 203 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 204 ENDIF 173 205 END IF 174 206 ! … … 185 217 CASE( 2 ) ! 2nd order FCT 186 218 ! 187 IF( l_trd ) zltv(:,:,:) = pt(:,:,:,jn,Krhs) ! store pt(:,:,:,:,Krhs) if trend diag. 219 IF( l_trd ) THEN 220 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 221 zltv(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) ! store pt(:,:,:,:,Krhs) if trend diag. 222 END_3D 223 ENDIF 188 224 ! 189 225 ! !* upstream advection with initial mass fluxes & intermediate update ==! … … 194 230 END_3D 195 231 IF( ln_linssh ) THEN ! top ocean value (only in linear free surface as ztw has been w-masked) 232 ! TODO: NOT TESTED- requires isf 196 233 IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface 197 234 DO_2D( 1, 1, 1, 1 ) … … 199 236 END_2D 200 237 ELSE ! no cavities: only at the ocean surface 201 ztw(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kbb) 238 DO_2D( 1, 1, 1, 1 ) 239 ztw(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kbb) 240 END_2D 202 241 ENDIF 203 242 ENDIF … … 226 265 ztw(ji,jj,jk) = pW(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 227 266 END_3D 228 IF( ln_linssh ) ztw(:,:, 1 ) = pW(:,:,1) * pt(:,:,1,jn,Kmm) !!gm ISF & 4th COMPACT doesn't work 267 IF( ln_linssh ) THEN 268 DO_2D( 1, 1, 1, 1 ) 269 ztw(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kmm) !!gm ISF & 4th COMPACT doesn't work 270 END_2D 271 ENDIF 229 272 ! 230 273 END SELECT … … 235 278 END_3D 236 279 ! 280 ! TEMP: These changes not necessary after trd_tra is tiled 237 281 IF( l_trd ) THEN ! vertical advective trend diagnostics 238 282 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 239 z ltv(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) - zltv(ji,jj,jk) &283 ztrdz(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) - zltv(ji,jj,jk) & 240 284 & + pt(ji,jj,jk,jn,Kmm) * ( pW(ji,jj,jk) - pW(ji,jj,jk+1) ) & 241 285 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 242 286 END_3D 243 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zltv ) 287 288 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 289 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 290 291 ! TODO: TO BE TILED- trd_tra 292 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, ztrdz ) 293 294 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 295 ENDIF 244 296 ENDIF 245 297 ! … … 262 314 !! in-space based differencing for fluid 263 315 !!---------------------------------------------------------------------- 264 INTEGER , INTENT(in ) 265 REAL(wp), INTENT(in ) 266 REAL(wp), DIMENSION 267 REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) :: paft ! after field268 REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) :: pcc ! monotonic flux in the k direction316 INTEGER , INTENT(in ) :: Kmm ! time level index 317 REAL(wp), INTENT(in ) :: p2dt ! tracer time-step 318 REAL(wp), DIMENSION(jpi,jpj,jpk) :: pbef ! before field 319 REAL(wp), INTENT(inout), DIMENSION(ST_2D(nn_hls) ,jpk) :: paft ! after field 320 REAL(wp), INTENT(inout), DIMENSION(ST_2D(nn_hls) ,jpk) :: pcc ! monotonic flux in the k direction 269 321 ! 270 322 INTEGER :: ji, jj, jk ! dummy loop indices 271 323 INTEGER :: ikm1 ! local integer 272 324 REAL(wp) :: zpos, zneg, zbt, za, zb, zc, zbig, zrtrn ! local scalars 273 REAL(wp), DIMENSION( jpi,jpj,jpk) :: zbetup, zbetdo! 3D workspace325 REAL(wp), DIMENSION(ST_2D(nn_hls),jpk) :: zbetup, zbetdo ! 3D workspace 274 326 !!---------------------------------------------------------------------- 275 327 ! … … 281 333 ! -------------------- 282 334 ! ! large negative value (-zbig) inside land 283 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) ) 284 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) ) 335 DO_3D( 0, 0, 0, 0, 1, jpk ) 336 pbef(ji,jj,jk) = pbef(ji,jj,jk) * tmask(ji,jj,jk) - zbig * ( 1.e0 - tmask(ji,jj,jk) ) 337 paft(ji,jj,jk) = paft(ji,jj,jk) * tmask(ji,jj,jk) - zbig * ( 1.e0 - tmask(ji,jj,jk) ) 338 END_3D 285 339 ! 286 340 DO jk = 1, jpkm1 ! search maximum in neighbourhood … … 293 347 END DO 294 348 ! ! large positive value (+zbig) inside land 295 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) ) 296 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) ) 349 DO_3D( 0, 0, 0, 0, 1, jpk ) 350 pbef(ji,jj,jk) = pbef(ji,jj,jk) * tmask(ji,jj,jk) + zbig * ( 1.e0 - tmask(ji,jj,jk) ) 351 paft(ji,jj,jk) = paft(ji,jj,jk) * tmask(ji,jj,jk) + zbig * ( 1.e0 - tmask(ji,jj,jk) ) 352 END_3D 297 353 ! 298 354 DO jk = 1, jpkm1 ! search minimum in neighbourhood … … 305 361 END DO 306 362 ! ! restore masked values to zero 307 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) 308 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) 363 DO_3D( 0, 0, 0, 0, 1, jpk ) 364 pbef(ji,jj,jk) = pbef(ji,jj,jk) * tmask(ji,jj,jk) 365 paft(ji,jj,jk) = paft(ji,jj,jk) * tmask(ji,jj,jk) 366 END_3D 309 367 ! 310 368 ! Positive and negative part of fluxes and beta terms
Note: See TracChangeset
for help on using the changeset viewer.