- 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_cen.F90
r13295 r13516 12 12 !!---------------------------------------------------------------------- 13 13 USE dom_oce ! ocean space and time domain 14 ! TEMP: This change not necessary after trd_tra is tiled 15 USE domain, ONLY : dom_tile 14 16 USE eosbn2 ! equation of state 15 17 USE traadv_fct ! acces to routine interp_4th_cpt … … 71 73 INTEGER , INTENT(in ) :: kn_cen_h ! =2/4 (2nd or 4th order scheme) 72 74 INTEGER , INTENT(in ) :: kn_cen_v ! =2/4 (2nd or 4th order scheme) 75 ! TEMP: This can be ST_2D(nn_hls) after trd_tra is tiled 73 76 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pU, pV, pW ! 3 ocean volume flux components 74 77 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! tracers and RHS of tracer equation … … 76 79 INTEGER :: ji, jj, jk, jn ! dummy loop indices 77 80 INTEGER :: ierr ! local integer 78 REAL(wp) :: zC2t_u, zC4t_u ! local scalars 79 REAL(wp) :: zC2t_v, zC4t_v ! - - 80 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwx, zwy, zwz, ztu, ztv, ztw 81 ! TEMP: This change not necessary after trd_tra is tiled 82 INTEGER :: itile 83 REAL(wp) :: zC2t_u, zC2t_v ! local scalars 84 REAL(wp), DIMENSION(ST_2D(nn_hls),jpk) :: zwx, zwy, zwz, ztu, ztv, ztw, zltu, zltv 85 ! TEMP: This change not necessary after trd_tra is tiled 86 REAL(wp), DIMENSION(:,:,:,:), SAVE, ALLOCATABLE :: ztrdx, ztrdy, ztrdz 81 87 !!---------------------------------------------------------------------- 82 ! 83 IF( kt == kit000 ) THEN 84 IF(lwp) WRITE(numout,*) 85 IF(lwp) WRITE(numout,*) 'tra_adv_cen : centered advection scheme on ', cdtype, ' order h/v =', kn_cen_h,'/', kn_cen_v 86 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~ ' 88 ! TEMP: This change not necessary after trd_tra is tiled 89 itile = ntile 90 ! 91 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 92 IF( kt == kit000 ) THEN 93 IF(lwp) WRITE(numout,*) 94 IF(lwp) WRITE(numout,*) 'tra_adv_cen : centered advection scheme on ', cdtype, ' order h/v =', kn_cen_h,'/', kn_cen_v 95 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~ ' 96 ENDIF 97 ! ! set local switches 98 l_trd = .FALSE. 99 l_hst = .FALSE. 100 l_ptr = .FALSE. 101 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 102 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE. 103 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 104 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 105 106 ! TEMP: This can be ST_2D(nn_hls) after trd_tra is tiled 107 IF( kt == kit000 .AND. l_trd ) THEN 108 ALLOCATE( ztrdx(jpi,jpj,jpk,jpts), ztrdy(jpi,jpj,jpk,jpts), ztrdz(jpi,jpj,jpk,jpts) ) 109 ENDIF 87 110 ENDIF 88 ! ! set local switches89 l_trd = .FALSE.90 l_hst = .FALSE.91 l_ptr = .FALSE.92 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.93 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE.94 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &95 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE.96 111 ! 97 112 ! … … 110 125 ! 111 126 CASE( 4 ) !* 4th order centered 112 ztu(:,:,jpk) = 0._wp ! Bottom value : flux set to zero 113 ztv(:,:,jpk) = 0._wp 114 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 115 ztu(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk) 116 ztv(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk) 117 END_3D 118 CALL lbc_lnk_multi( 'traadv_cen', ztu, 'U', -1.0_wp , ztv, 'V', -1.0_wp ) ! Lateral boundary cond. 119 ! 120 DO_3D( 0, 0, 1, 0, 1, jpkm1 ) 127 zltu(:,:,jpk) = 0._wp ! Bottom value : flux set to zero 128 zltv(:,:,jpk) = 0._wp 129 DO jk = 1, jpkm1 ! Laplacian 130 DO_2D( 1, 0, 1, 0 ) 131 ztu(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk) 132 ztv(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk) 133 END_2D 134 DO_2D( 0, 0, 0, 0 ) 135 zltu(ji,jj,jk) = ztu(ji,jj,jk) + ztu(ji-1,jj,jk) 136 zltv(ji,jj,jk) = ztv(ji,jj,jk) + ztv(ji,jj-1,jk) 137 END_2D 138 END DO 139 CALL lbc_lnk_multi( 'traadv_cen', zltu, 'T', 1. , zltv, 'T', 1. ) ! Lateral boundary cond. (unchanged sgn) 140 ! 141 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 121 142 zC2t_u = pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj ,jk,jn,Kmm) ! C2 interpolation of T at u- & v-points (x2) 122 143 zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji ,jj+1,jk,jn,Kmm) 123 ! ! C4 interpolation of T at u- & v-points (x2)124 zC4t_u = zC2t_u + r1_6 * ( ztu(ji-1,jj,jk) - ztu(ji+1,jj,jk) )125 zC4t_v = zC2t_v + r1_6 * ( ztv(ji,jj-1,jk) - ztv(ji,jj+1,jk) )126 144 ! ! C4 fluxes 127 zwx(ji,jj,jk) = 0.5_wp * pU(ji,jj,jk) * zC4t_u128 zwy(ji,jj,jk) = 0.5_wp * pV(ji,jj,jk) * zC4t_v145 zwx(ji,jj,jk) = 0.5_wp * pU(ji,jj,jk) * ( zC2t_u + r1_6 * (zltu(ji,jj,jk) - zltu(ji+1,jj,jk)) ) 146 zwy(ji,jj,jk) = 0.5_wp * pV(ji,jj,jk) * ( zC2t_v + r1_6 * (zltv(ji,jj,jk) - zltv(ji,jj+1,jk)) ) 129 147 END_3D 130 148 ! … … 149 167 ! 150 168 IF( ln_linssh ) THEN !* top value (linear free surf. only as zwz is multiplied by wmask) 169 ! TODO: NOT TESTED- requires isf 151 170 IF( ln_isfcav ) THEN ! ice-shelf cavities (top of the ocean) 152 171 DO_2D( 1, 1, 1, 1 ) … … 154 173 END_2D 155 174 ELSE ! no ice-shelf cavities (only ocean surface) 156 zwz(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kmm) 175 DO_2D( 1, 1, 1, 1 ) 176 zwz(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kmm) 177 END_2D 157 178 ENDIF 158 179 ENDIF … … 166 187 END_3D 167 188 ! ! trend diagnostics 189 ! TEMP: These changes not necessary after trd_tra is tiled 168 190 IF( l_trd ) THEN 169 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kmm) ) 170 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kmm) ) 171 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zwz, pW, pt(:,:,:,jn,Kmm) ) 172 END IF 173 ! ! "Poleward" heat and salt transports 191 DO_3D( 1, 0, 1, 0, 1, jpk ) 192 ztrdx(ji,jj,jk,jn) = zwx(ji,jj,jk) 193 ztrdy(ji,jj,jk,jn) = zwy(ji,jj,jk) 194 ztrdz(ji,jj,jk,jn) = zwz(ji,jj,jk) 195 END_3D 196 197 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 198 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 199 200 ! TODO: TO BE TILED- trd_tra 201 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztrdx(:,:,:,jn), pU, pt(:,:,:,jn,Kmm) ) 202 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztrdy(:,:,:,jn), pV, pt(:,:,:,jn,Kmm) ) 203 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, ztrdz(:,:,:,jn), pW, pt(:,:,:,jn,Kmm) ) 204 205 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 206 ENDIF 207 ENDIF 208 ! ! "Poleward" heat and salt transports 174 209 IF( l_ptr ) CALL dia_ptr_hst( jn, 'adv', zwy(:,:,:) ) 175 210 ! ! heat and salt transport
Note: See TracChangeset
for help on using the changeset viewer.