- Timestamp:
- 2016-11-28T17:04:10+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_mus.F90
r5930 r7351 41 41 42 42 !! * Substitutions 43 # include "domzgr_substitute.h90"44 43 # include "vectopt_loop_substitute.h90" 45 44 !!---------------------------------------------------------------------- … … 62 61 !! ld_msc_ups=T : 63 62 !! 64 !! ** Action : - update (ta,sa) with the now advective tracer trends 65 !! - send trends to trdtra module for further diagnostcs 63 !! ** Action : - update pta with the now advective tracer trends 64 !! - send trends to trdtra module for further diagnostcs (l_trdtra=T) 65 !! - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 66 66 !! 67 67 !! References : Estubier, A., and M. Levy, Notes Techn. Pole de Modelisation … … 73 73 INTEGER , INTENT(in ) :: kjpt ! number of tracers 74 74 LOGICAL , INTENT(in ) :: ld_msc_ups ! use upstream scheme within muscl 75 REAL(wp) , DIMENSION( jpk ), INTENT(in ) :: p2dt ! vertical profile oftracer time-step75 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 76 76 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pun, pvn, pwn ! 3 ocean velocity components 77 77 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before tracer field … … 82 82 REAL(wp) :: zu, z0u, zzwx, zw ! local scalars 83 83 REAL(wp) :: zv, z0v, zzwy, z0w ! - - 84 REAL(wp) :: z dt, zalpha! - -84 REAL(wp) :: zalpha ! - - 85 85 REAL(wp), POINTER, DIMENSION(:,:,:) :: zslpx, zslpy ! 3D workspace 86 86 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwx , zwy ! - - … … 116 116 ENDIF 117 117 ! 118 ! ! =========== 119 DO jn = 1, kjpt ! tracer loop 120 ! ! =========== 121 ! I. Horizontal advective fluxes 122 ! ------------------------------ 123 ! first guess of the slopes 124 zwx(:,:,jpk) = 0.e0 ; zwy(:,:,jpk) = 0.e0 ! bottom values 125 ! interior values 126 DO jk = 1, jpkm1 118 DO jn = 1, kjpt !== loop over the tracers ==! 119 ! 120 ! !* Horizontal advective fluxes 121 ! 122 ! !-- first guess of the slopes 123 zwx(:,:,jpk) = 0._wp ! bottom values 124 zwy(:,:,jpk) = 0._wp 125 DO jk = 1, jpkm1 ! interior values 127 126 DO jj = 1, jpjm1 128 127 DO ji = 1, fs_jpim1 ! vector opt. … … 132 131 END DO 133 132 END DO 134 ! 135 CALL lbc_lnk( zwx, 'U', -1. ) ! lateral boundary conditions on zwx, zwy (changed sign) 133 CALL lbc_lnk( zwx, 'U', -1. ) ! lateral boundary conditions (changed sign) 136 134 CALL lbc_lnk( zwy, 'V', -1. ) 137 ! !-- Slopes of tracer 138 zslpx(:,:,jpk) = 0.e0 ; zslpy(:,:,jpk) = 0.e0 ! bottom values 139 DO jk = 1, jpkm1 ! interior values 135 ! !-- Slopes of tracer 136 zslpx(:,:,jpk) = 0._wp ! bottom values 137 zslpy(:,:,jpk) = 0._wp 138 DO jk = 1, jpkm1 ! interior values 140 139 DO jj = 2, jpj 141 140 DO ji = fs_2, jpi ! vector opt. … … 148 147 END DO 149 148 ! 150 DO jk = 1, jpkm1 !Slopes limitation149 DO jk = 1, jpkm1 !-- Slopes limitation 151 150 DO jj = 2, jpj 152 151 DO ji = fs_2, jpi ! vector opt. … … 161 160 END DO 162 161 ! 163 ! !-- MUSCL horizontal advective fluxes 164 DO jk = 1, jpkm1 ! interior values 165 zdt = p2dt(jk) 162 DO jk = 1, jpkm1 !-- MUSCL horizontal advective fluxes 166 163 DO jj = 2, jpjm1 167 164 DO ji = fs_2, fs_jpim1 ! vector opt. … … 169 166 z0u = SIGN( 0.5, pun(ji,jj,jk) ) 170 167 zalpha = 0.5 - z0u 171 zu = z0u - 0.5 * pun(ji,jj,jk) * zdt / ( e1e2u(ji,jj) * fse3u(ji,jj,jk))168 zu = z0u - 0.5 * pun(ji,jj,jk) * p2dt * r1_e1e2u(ji,jj) / e3u_n(ji,jj,jk) 172 169 zzwx = ptb(ji+1,jj,jk,jn) + xind(ji,jj,jk) * zu * zslpx(ji+1,jj,jk) 173 170 zzwy = ptb(ji ,jj,jk,jn) + xind(ji,jj,jk) * zu * zslpx(ji ,jj,jk) … … 176 173 z0v = SIGN( 0.5, pvn(ji,jj,jk) ) 177 174 zalpha = 0.5 - z0v 178 zv = z0v - 0.5 * pvn(ji,jj,jk) * zdt / ( e1e2v(ji,jj) * fse3v(ji,jj,jk))175 zv = z0v - 0.5 * pvn(ji,jj,jk) * p2dt * r1_e1e2v(ji,jj) / e3v_n(ji,jj,jk) 179 176 zzwx = ptb(ji,jj+1,jk,jn) + xind(ji,jj,jk) * zv * zslpy(ji,jj+1,jk) 180 177 zzwy = ptb(ji,jj ,jk,jn) + xind(ji,jj,jk) * zv * zslpy(ji,jj ,jk) … … 185 182 CALL lbc_lnk( zwx, 'U', -1. ) ; CALL lbc_lnk( zwy, 'V', -1. ) ! lateral boundary conditions (changed sign) 186 183 ! 187 DO jk = 1, jpkm1 ! Tracer flux divergence at t-point added to the generaltrend184 DO jk = 1, jpkm1 !-- Tracer advective trend 188 185 DO jj = 2, jpjm1 189 186 DO ji = fs_2, fs_jpim1 ! vector opt. 190 187 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 191 188 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) ) & 192 & / ( e1e2t(ji,jj) * fse3t(ji,jj,jk))189 & * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 193 190 END DO 194 191 END DO 195 192 END DO 196 ! ! trend diagnostics (contribution of upstream fluxes)193 ! ! trend diagnostics 197 194 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. & 198 195 &( cdtype == 'TRC' .AND. l_trdtrc ) ) THEN … … 200 197 CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pvn, ptb(:,:,:,jn) ) 201 198 END IF 202 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes)199 ! ! "Poleward" heat and salt transports 203 200 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 204 201 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) 205 202 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) 206 203 ENDIF 207 208 ! II.Vertical advective fluxes209 ! -----------------------------204 ! 205 ! !* Vertical advective fluxes 206 ! 210 207 ! !-- first guess of the slopes 211 208 zwx(:,:, 1 ) = 0._wp ! surface & bottom boundary conditions 212 zwx(:,:,jpk) = 0._wp ! surface & bottom boundary conditions213 DO jk = 2, jpkm1 209 zwx(:,:,jpk) = 0._wp 210 DO jk = 2, jpkm1 ! interior values 214 211 zwx(:,:,jk) = tmask(:,:,jk) * ( ptb(:,:,jk-1,jn) - ptb(:,:,jk,jn) ) 215 212 END DO 216 217 213 ! !-- Slopes of tracer 218 214 zslpx(:,:,1) = 0._wp ! surface values … … 220 216 DO jj = 1, jpj 221 217 DO ji = 1, jpi 222 zslpx(ji,jj,jk) = ( zwx(ji,jj,jk) + zwx(ji,jj,jk+1) ) & 223 & * ( 0.25 + SIGN( 0.25, zwx(ji,jj,jk) * zwx(ji,jj,jk+1) ) ) 224 END DO 225 END DO 226 END DO 227 ! !-- Slopes limitation 228 DO jk = 2, jpkm1 ! interior values 229 DO jj = 1, jpj 218 zslpx(ji,jj,jk) = ( zwx(ji,jj,jk) + zwx(ji,jj,jk+1) ) & 219 & * ( 0.25 + SIGN( 0.25, zwx(ji,jj,jk) * zwx(ji,jj,jk+1) ) ) 220 END DO 221 END DO 222 END DO 223 DO jk = 2, jpkm1 !-- Slopes limitation 224 DO jj = 1, jpj ! interior values 230 225 DO ji = 1, jpi 231 226 zslpx(ji,jj,jk) = SIGN( 1., zslpx(ji,jj,jk) ) * MIN( ABS( zslpx(ji,jj,jk ) ), & … … 235 230 END DO 236 231 END DO 237 ! !-- vertical advective flux 238 DO jk = 1, jpkm1 ! interior values 239 zdt = p2dt(jk) 232 DO jk = 1, jpk-2 !-- vertical advective flux 240 233 DO jj = 2, jpjm1 241 234 DO ji = fs_2, fs_jpim1 ! vector opt. 242 235 z0w = SIGN( 0.5, pwn(ji,jj,jk+1) ) 243 236 zalpha = 0.5 + z0w 244 zw = z0w - 0.5 * pwn(ji,jj,jk+1) * zdt / ( e1e2t(ji,jj) * fse3w(ji,jj,jk+1))237 zw = z0w - 0.5 * pwn(ji,jj,jk+1) * p2dt * r1_e1e2t(ji,jj) / e3w_n(ji,jj,jk+1) 245 238 zzwx = ptb(ji,jj,jk+1,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk+1) 246 239 zzwy = ptb(ji,jj,jk ,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk ) … … 249 242 END DO 250 243 END DO 251 ! ! top values (bottom already set to zero) 252 IF( lk_vvl ) THEN ! variable volume 253 zwx(:,:, 1 ) = 0._wp ! k=1 only as zwx has been multiplied by wmask 254 ELSE ! linear free surface 255 IF( ln_isfcav ) THEN ! ice-shelf cavities (top of the ocean) 244 IF( ln_linssh ) THEN ! top values, linear free surface only 245 IF( ln_isfcav ) THEN ! ice-shelf cavities (top of the ocean) 256 246 DO jj = 1, jpj 257 247 DO ji = 1, jpi … … 259 249 END DO 260 250 END DO 261 ELSE 251 ELSE ! no cavities: only at the ocean surface 262 252 zwx(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn) 263 253 ENDIF 264 254 ENDIF 265 255 ! 266 DO jk = 1, jpkm1 ! Compute & add thevertical advective trend256 DO jk = 1, jpkm1 !-- vertical advective trend 267 257 DO jj = 2, jpjm1 268 258 DO ji = fs_2, fs_jpim1 ! vector opt. 269 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) ) / ( e1e2t(ji,jj) * fse3t(ji,jj,jk))270 END DO 271 END DO 272 END DO 273 ! ! Save the vertical advectivetrends for diagnostic259 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 260 END DO 261 END DO 262 END DO 263 ! ! send trends for diagnostic 274 264 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. & 275 265 &( cdtype == 'TRC' .AND. l_trdtrc ) ) & 276 266 CALL trd_tra( kt, cdtype, jn, jptra_zad, zwx, pwn, ptb(:,:,:,jn) ) 277 267 ! 278 END DO 268 END DO ! end of tracer loop 279 269 ! 280 270 CALL wrk_dealloc( jpi,jpj,jpk, zslpx, zslpy, zwx, zwy )
Note: See TracChangeset
for help on using the changeset viewer.