Changeset 14789 for NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU/src/OCE/TRA/traadv_ubs.F90
- Timestamp:
- 2021-05-05T13:18:04+02:00 (3 years ago)
- Location:
- NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
src/OCE/TRA/traadv_ubs.F90 (modified) (20 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU
- Property svn:externals
-
old new 3 3 ^/utils/build/mk@HEAD mk 4 4 ^/utils/tools@HEAD tools 5 ^/vendors/AGRIF/dev _r12970_AGRIF_CMEMSext/AGRIF5 ^/vendors/AGRIF/dev@HEAD ext/AGRIF 6 6 ^/vendors/FCM@HEAD ext/FCM 7 7 ^/vendors/IOIPSL@HEAD ext/IOIPSL 8 ^/vendors/PPR@HEAD ext/PPR 8 9 9 10 # SETTE 10 ^/utils/CI/sette@1 3559sette11 ^/utils/CI/sette@14244 sette
-
- Property svn:externals
-
NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU/src/OCE/TRA/traadv_ubs.F90
r13497 r14789 10 10 !!---------------------------------------------------------------------- 11 11 !! tra_adv_ubs : update the tracer trend with the horizontal 12 !! advection trends using a third order biaised scheme 12 !! advection trends using a third order biaised scheme 13 13 !!---------------------------------------------------------------------- 14 14 USE oce ! ocean dynamics and active tracers … … 16 16 USE trc_oce ! share passive tracers/Ocean variables 17 17 USE trd_oce ! trends: ocean variables 18 USE traadv_fct ! acces to routine interp_4th_cpt 19 USE trdtra ! trends manager: tracers 18 USE traadv_fct ! acces to routine interp_4th_cpt 19 USE trdtra ! trends manager: tracers 20 20 USE diaptr ! poleward transport diagnostics 21 21 USE diaar5 ! AR5 diagnostics … … 25 25 USE lib_mpp ! massively parallel library 26 26 USE lbclnk ! ocean lateral boundary condition (or mpp link) 27 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 27 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 28 28 29 29 IMPLICIT NONE … … 51 51 !!---------------------------------------------------------------------- 52 52 !! *** ROUTINE tra_adv_ubs *** 53 !! 53 !! 54 54 !! ** Purpose : Compute the now trend due to the advection of tracers 55 55 !! and add it to the general trend of passive tracer equations. … … 60 60 !! For example the i-component of the advective fluxes are given by : 61 61 !! ! e2u e3u un ( mi(Tn) - zltu(i ) ) if un(i) >= 0 62 !! ztu = ! or 62 !! ztu = ! or 63 63 !! ! e2u e3u un ( mi(Tn) - zltu(i+1) ) if un(i) < 0 64 64 !! where zltu is the second derivative of the before temperature field: 65 65 !! zltu = 1/e3t di[ e2u e3u / e1u di[Tb] ] 66 !! This results in a dissipatively dominant (i.e. hyper-diffusive) 67 !! truncation error. The overall performance of the advection scheme 68 !! is similar to that reported in (Farrow and Stevens, 1995). 66 !! This results in a dissipatively dominant (i.e. hyper-diffusive) 67 !! truncation error. The overall performance of the advection scheme 68 !! is similar to that reported in (Farrow and Stevens, 1995). 69 69 !! For stability reasons, the first term of the fluxes which corresponds 70 !! to a second order centered scheme is evaluated using the now velocity 71 !! (centered in time) while the second term which is the diffusive part 72 !! of the scheme, is evaluated using the before velocity (forward in time). 70 !! to a second order centered scheme is evaluated using the now velocity 71 !! (centered in time) while the second term which is the diffusive part 72 !! of the scheme, is evaluated using the before velocity (forward in time). 73 73 !! Note that UBS is not positive. Do not use it on passive tracers. 74 74 !! On the vertical, the advection is evaluated using a FCT scheme, 75 !! as the UBS have been found to be too diffusive. 76 !! kn_ubs_v argument controles whether the FCT is based on 77 !! a 2nd order centrered scheme (kn_ubs_v=2) or on a 4th order compact 75 !! as the UBS have been found to be too diffusive. 76 !! kn_ubs_v argument controles whether the FCT is based on 77 !! a 2nd order centrered scheme (kn_ubs_v=2) or on a 4th order compact 78 78 !! scheme (kn_ubs_v=4). 79 79 !! … … 82 82 !! - poleward advective heat and salt transport (ln_diaptr=T) 83 83 !! 84 !! Reference : Shchepetkin, A. F., J. C. McWilliams, 2005, Ocean Modelling, 9, 347-404. 84 !! Reference : Shchepetkin, A. F., J. C. McWilliams, 2005, Ocean Modelling, 9, 347-404. 85 85 !! Farrow, D.E., Stevens, D.P., 1995, J. Phys. Ocean. 25, 1731�1741. 86 86 !!---------------------------------------------------------------------- … … 92 92 INTEGER , INTENT(in ) :: kn_ubs_v ! number of tracers 93 93 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 94 ! TEMP: [tiling] This can be A2D(nn_hls) if using XIOS (subdomain support) 94 95 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pU, pV, pW ! 3 ocean volume transport components 95 96 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! tracers and RHS of tracer equation … … 99 100 REAL(wp) :: zfp_ui, zfm_ui, zcenut, ztak, zfp_wk, zfm_wk ! - - 100 101 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,*) '~~~~~~~~~~~~' 102 REAL(wp), DIMENSION(A2D(nn_hls),jpk) :: ztu, ztv, zltu, zltv, zti, ztw ! 3D workspace 103 !!---------------------------------------------------------------------- 104 ! 105 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 106 IF( kt == kit000 ) THEN 107 IF(lwp) WRITE(numout,*) 108 IF(lwp) WRITE(numout,*) 'tra_adv_ubs : horizontal UBS advection scheme on ', cdtype 109 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 110 ENDIF 111 ! 112 l_trd = .FALSE. 113 l_hst = .FALSE. 114 l_ptr = .FALSE. 115 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 116 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE. 117 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 118 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 108 119 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 120 ! 118 121 ztw (:,:, 1 ) = 0._wp ! surface & bottom value : set to zero for all tracers 119 122 zltu(:,:,jpk) = 0._wp ; zltv(:,:,jpk) = 0._wp 120 123 ztw (:,:,jpk) = 0._wp ; zti (:,:,jpk) = 0._wp 121 !122 124 ! ! =========== 123 125 DO jn = 1, kjpt ! tracer loop 124 126 ! ! =========== 125 ! 127 ! 126 128 DO jk = 1, jpkm1 !== horizontal laplacian of before tracer ==! 127 DO_2D( 1, 0, 1, 0) ! First derivative (masked gradient)129 DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) ! First derivative (masked gradient) 128 130 zeeu = e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm) * umask(ji,jj,jk) 129 131 zeev = e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm) * vmask(ji,jj,jk) … … 131 133 ztv(ji,jj,jk) = zeev * ( pt(ji ,jj+1,jk,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) ) 132 134 END_2D 133 DO_2D( 0, 0, 0, 0) ! Second derivative (divergence)135 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) ! Second derivative (divergence) 134 136 zcoef = 1._wp / ( 6._wp * e3t(ji,jj,jk,Kmm) ) 135 137 zltu(ji,jj,jk) = ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) ) * zcoef 136 138 zltv(ji,jj,jk) = ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) * zcoef 137 139 END_2D 138 ! 139 END DO 140 CALL lbc_lnk( 'traadv_ubs', zltu, 'T', 1.0_wp ) ; CALL lbc_lnk( 'traadv_ubs', zltv, 'T', 1.0_wp ) ! Lateral boundary cond. (unchanged sgn)141 ! 140 ! 141 END DO 142 IF (nn_hls.EQ.1) CALL lbc_lnk( 'traadv_ubs', zltu, 'T', 1.0_wp, zltv, 'T', 1.0_wp ) ! Lateral boundary cond. (unchanged sgn) 143 ! 142 144 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) !== Horizontal advective fluxes ==! (UBS) 143 145 zfp_ui = pU(ji,jj,jk) + ABS( pU(ji,jj,jk) ) ! upstream transport (x2) … … 153 155 END_3D 154 156 ! 155 zltu(:,:,:) = pt(:,:,:,jn,Krhs) ! store the initial trends before its update 157 DO_3D( 1, 1, 1, 1, 1, jpk ) 158 zltu(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) ! store the initial trends before its update 159 END_3D 156 160 ! 157 161 DO jk = 1, jpkm1 !== add the horizontal advective trend ==! … … 162 166 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 163 167 END_2D 164 ! 168 ! 165 169 END DO 166 170 ! 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 ! 171 DO_3D( 1, 1, 1, 1, 1, jpk ) 172 zltu(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) - zltu(ji,jj,jk) ! Horizontal advective trend used in vertical 2nd order FCT case 173 END_3D ! and/or in trend diagnostic (l_trd=T) 174 ! 170 175 IF( l_trd ) THEN ! trend diagnostics 171 176 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztu, pU, pt(:,:,:,jn,Kmm) ) 172 177 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztv, pV, pt(:,:,:,jn,Kmm) ) 173 178 END IF 174 ! 179 ! 175 180 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 176 181 IF( l_ptr ) CALL dia_ptr_hst( jn, 'adv', ztv(:,:,:) ) … … 183 188 SELECT CASE( kn_ubs_v ) ! select the vertical advection scheme 184 189 ! 185 CASE( 2 ) ! 2nd order FCT 186 ! 187 IF( l_trd ) zltv(:,:,:) = pt(:,:,:,jn,Krhs) ! store pt(:,:,:,:,Krhs) if trend diag. 190 CASE( 2 ) ! 2nd order FCT 191 ! 192 IF( l_trd ) THEN 193 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 194 zltv(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) ! store pt(:,:,:,:,Krhs) if trend diag. 195 END_3D 196 ENDIF 188 197 ! 189 198 ! !* upstream advection with initial mass fluxes & intermediate update ==! … … 196 205 IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface 197 206 DO_2D( 1, 1, 1, 1 ) 198 ztw(ji,jj, mikt(ji,jj) ) = pW(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kbb) ! linear free surface 207 ztw(ji,jj, mikt(ji,jj) ) = pW(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kbb) ! linear free surface 199 208 END_2D 200 209 ELSE ! no cavities: only at the ocean surface 201 ztw(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kbb) 210 DO_2D( 1, 1, 1, 1 ) 211 ztw(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kbb) 212 END_2D 202 213 ENDIF 203 214 ENDIF … … 206 217 ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) & 207 218 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 208 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) + ztak 219 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) + ztak 209 220 zti(ji,jj,jk) = ( pt(ji,jj,jk,jn,Kbb) + p2dt * ( ztak + zltu(ji,jj,jk) ) ) * tmask(ji,jj,jk) 210 221 END_3D 211 CALL lbc_lnk( 'traadv_ubs', zti, 'T', 1.0_wp ) ! Lateral boundary conditions on zti, zsi (unchanged sign)212 222 ! 213 223 ! !* anti-diffusive flux : high order minus low order … … 226 236 ztw(ji,jj,jk) = pW(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 227 237 END_3D 228 IF( ln_linssh ) ztw(:,:, 1 ) = pW(:,:,1) * pt(:,:,1,jn,Kmm) !!gm ISF & 4th COMPACT doesn't work 238 IF( ln_linssh ) THEN 239 DO_2D( 1, 1, 1, 1 ) 240 ztw(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kmm) !!gm ISF & 4th COMPACT doesn't work 241 END_2D 242 ENDIF 229 243 ! 230 244 END SELECT … … 252 266 !!--------------------------------------------------------------------- 253 267 !! *** ROUTINE nonosc_z *** 254 !! 255 !! ** Purpose : compute monotonic tracer fluxes from the upstream 256 !! scheme and the before field by a nonoscillatory algorithm 268 !! 269 !! ** Purpose : compute monotonic tracer fluxes from the upstream 270 !! scheme and the before field by a nonoscillatory algorithm 257 271 !! 258 272 !! ** Method : ... ??? … … 262 276 !! in-space based differencing for fluid 263 277 !!---------------------------------------------------------------------- 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 direction278 INTEGER , INTENT(in ) :: Kmm ! time level index 279 REAL(wp), INTENT(in ) :: p2dt ! tracer time-step 280 REAL(wp), DIMENSION(jpi,jpj,jpk) :: pbef ! before field 281 REAL(wp), INTENT(inout), DIMENSION(A2D(nn_hls) ,jpk) :: paft ! after field 282 REAL(wp), INTENT(inout), DIMENSION(A2D(nn_hls) ,jpk) :: pcc ! monotonic flux in the k direction 269 283 ! 270 284 INTEGER :: ji, jj, jk ! dummy loop indices 271 285 INTEGER :: ikm1 ! local integer 272 286 REAL(wp) :: zpos, zneg, zbt, za, zb, zc, zbig, zrtrn ! local scalars 273 REAL(wp), DIMENSION( jpi,jpj,jpk) :: zbetup, zbetdo! 3D workspace287 REAL(wp), DIMENSION(A2D(nn_hls),jpk) :: zbetup, zbetdo ! 3D workspace 274 288 !!---------------------------------------------------------------------- 275 289 ! … … 281 295 ! -------------------- 282 296 ! ! large negative value (-zbig) inside land 283 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) ) 284 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) ) 297 DO_3D( 0, 0, 0, 0, 1, jpk ) 298 pbef(ji,jj,jk) = pbef(ji,jj,jk) * tmask(ji,jj,jk) - zbig * ( 1.e0 - tmask(ji,jj,jk) ) 299 paft(ji,jj,jk) = paft(ji,jj,jk) * tmask(ji,jj,jk) - zbig * ( 1.e0 - tmask(ji,jj,jk) ) 300 END_3D 285 301 ! 286 302 DO jk = 1, jpkm1 ! search maximum in neighbourhood … … 293 309 END DO 294 310 ! ! large positive value (+zbig) inside land 295 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) ) 296 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) ) 311 DO_3D( 0, 0, 0, 0, 1, jpk ) 312 pbef(ji,jj,jk) = pbef(ji,jj,jk) * tmask(ji,jj,jk) + zbig * ( 1.e0 - tmask(ji,jj,jk) ) 313 paft(ji,jj,jk) = paft(ji,jj,jk) * tmask(ji,jj,jk) + zbig * ( 1.e0 - tmask(ji,jj,jk) ) 314 END_3D 297 315 ! 298 316 DO jk = 1, jpkm1 ! search minimum in neighbourhood … … 305 323 END DO 306 324 ! ! restore masked values to zero 307 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) 308 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) 325 DO_3D( 0, 0, 0, 0, 1, jpk ) 326 pbef(ji,jj,jk) = pbef(ji,jj,jk) * tmask(ji,jj,jk) 327 paft(ji,jj,jk) = paft(ji,jj,jk) * tmask(ji,jj,jk) 328 END_3D 309 329 ! 310 330 ! Positive and negative part of fluxes and beta terms
Note: See TracChangeset
for help on using the changeset viewer.
