- Timestamp:
- 2020-09-24T20:32:14+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r13383_HPC-02_Daley_Tiling/src/OCE/TRA/traldf_iso.F90
r13295 r13515 19 19 USE oce ! ocean dynamics and active tracers 20 20 USE dom_oce ! ocean space and time domain 21 USE domutl, ONLY : is_tile 21 22 USE trc_oce ! share passive tracers/Ocean variables 22 23 USE zdf_oce ! ocean vertical physics … … 49 50 CONTAINS 50 51 51 SUBROUTINE tra_ldf_iso( kt, Kmm, kit000, cdtype, pahu, pahv, & 52 & pgu , pgv , pgui, pgvi, & 53 & pt , pt2 , pt_rhs , kjpt , kpass ) 52 SUBROUTINE tra_ldf_iso( kt, Kmm, kit000, cdtype, pahu, pahv, & 53 & pgu , pgv , pgui, pgvi, & 54 & pt, pt2, pt_rhs, kjpt, kpass ) 55 !! 56 INTEGER , INTENT(in ) :: kt ! ocean time-step index 57 INTEGER , INTENT(in ) :: kit000 ! first time step index 58 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 59 INTEGER , INTENT(in ) :: kjpt ! number of tracers 60 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 61 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 62 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 63 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 64 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 65 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 66 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pt2 ! tracer (only used in kpass=2) 67 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pt_rhs ! tracer trend 68 !! 69 CALL tra_ldf_iso_t( kt, Kmm, kit000, cdtype, pahu, pahv, is_tile(pahu), & 70 & pgu , pgv , is_tile(pgu) , pgui, pgvi, is_tile(pgui), & 71 & pt, is_tile(pt), pt2, is_tile(pt2), pt_rhs, is_tile(pt_rhs), kjpt, kpass ) 72 END SUBROUTINE tra_ldf_iso 73 74 75 SUBROUTINE tra_ldf_iso_t( kt, Kmm, kit000, cdtype, pahu, pahv, ktah, & 76 & pgu , pgv , ktg , pgui, pgvi, ktgi, & 77 & pt, ktt, pt2, ktt2, pt_rhs, ktt_rhs, kjpt, kpass ) 54 78 !!---------------------------------------------------------------------- 55 79 !! *** ROUTINE tra_ldf_iso *** … … 92 116 !! ** Action : Update pt_rhs arrays with the before rotated diffusion 93 117 !!---------------------------------------------------------------------- 94 INTEGER , INTENT(in ) :: kt ! ocean time-step index 95 INTEGER , INTENT(in ) :: kit000 ! first time step index 96 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 97 INTEGER , INTENT(in ) :: kjpt ! number of tracers 98 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 99 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 100 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 101 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 102 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 103 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 104 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt2 ! tracer (only used in kpass=2) 105 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 118 INTEGER , INTENT(in ) :: kt ! ocean time-step index 119 INTEGER , INTENT(in ) :: kit000 ! first time step index 120 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 121 INTEGER , INTENT(in ) :: kjpt ! number of tracers 122 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 123 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 124 INTEGER , INTENT(in ) :: ktah, ktg, ktgi, ktt, ktt2, ktt_rhs 125 REAL(wp), DIMENSION(ST_2DT(ktah) ,jpk) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 126 REAL(wp), DIMENSION(ST_2DT(ktg) ,kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 127 REAL(wp), DIMENSION(ST_2DT(ktgi) ,kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 128 REAL(wp), DIMENSION(ST_2DT(ktt) ,jpk,kjpt), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 129 REAL(wp), DIMENSION(ST_2DT(ktt2) ,jpk,kjpt), INTENT(in ) :: pt2 ! tracer (only used in kpass=2) 130 REAL(wp), DIMENSION(ST_2DT(ktt_rhs),jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 106 131 ! 107 132 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 111 136 REAL(wp) :: zmskv, zahv_w, zabe2, zcof2, zcoef4 ! - - 112 137 REAL(wp) :: zcoef0, ze3w_2, zsign ! - - 113 REAL(wp), DIMENSION( jpi,jpj) :: zdkt, zdk1t, z2d114 REAL(wp), DIMENSION( jpi,jpj,jpk) :: zdit, zdjt, zftu, zftv, ztfw138 REAL(wp), DIMENSION(ST_2D(nn_hls)) :: zdkt, zdk1t, z2d 139 REAL(wp), DIMENSION(ST_2D(nn_hls),jpk) :: zdit, zdjt, zftu, zftv, ztfw 115 140 !!---------------------------------------------------------------------- 116 141 ! 117 142 IF( kpass == 1 .AND. kt == kit000 ) THEN 118 IF(lwp) WRITE(numout,*) 119 IF(lwp) WRITE(numout,*) 'tra_ldf_iso : rotated laplacian diffusion operator on ', cdtype 120 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 121 ! 122 akz (:,:,:) = 0._wp 123 ah_wslp2(:,:,:) = 0._wp 143 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 144 IF(lwp) WRITE(numout,*) 145 IF(lwp) WRITE(numout,*) 'tra_ldf_iso : rotated laplacian diffusion operator on ', cdtype 146 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 147 ENDIF 148 ! 149 DO_3D( 0, 0, 0, 0, 1, jpk ) 150 akz (ji,jj,jk) = 0._wp 151 ah_wslp2(ji,jj,jk) = 0._wp 152 END_3D 124 153 ENDIF 125 ! 126 l_hst = .FALSE. 127 l_ptr = .FALSE. 128 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) ) l_ptr = .TRUE. 129 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 130 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 154 ! 155 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 156 l_hst = .FALSE. 157 l_ptr = .FALSE. 158 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) ) l_ptr = .TRUE. 159 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 160 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 161 ENDIF 131 162 ! 132 163 ! … … 167 198 ! 168 199 IF( ln_traldf_blp ) THEN ! bilaplacian operator 169 DO_3D( 1, 0, 1, 0, 2, jpkm1 ) 170 akz(ji,jj,jk) = 16._wp & 171 & * ah_wslp2 (ji,jj,jk) & 172 & * ( akz (ji,jj,jk) & 173 & + ah_wslp2(ji,jj,jk) & 174 & / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) ) ) 200 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 201 akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk) & 202 & * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) ) ) 175 203 END_3D 176 204 ELSEIF( ln_traldf_lap ) THEN ! laplacian operator 177 DO_3D( 1, 0, 1, 0, 2, jpkm1 )205 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 178 206 ze3w_2 = e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) 179 207 zcoef0 = rDt * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2 ) … … 183 211 ! 184 212 ELSE ! 33 flux set to zero with akz=ah_wslp2 ==>> computed in full implicit 185 akz(:,:,:) = ah_wslp2(:,:,:) 213 DO_3D( 0, 0, 0, 0, 1, jpk ) 214 akz(ji,jj,jk) = ah_wslp2(ji,jj,jk) 215 END_3D 186 216 ENDIF 187 217 ENDIF … … 195 225 !!---------------------------------------------------------------------- 196 226 !!gm : bug.... why (x,:,:)? (1,jpj,:) and (jpi,1,:) should be sufficient.... 197 zdit ( 1,:,:) = 0._wp ; zdit (jpi,:,:) = 0._wp198 zdjt ( 1,:,:) = 0._wp ; zdjt (jpi,:,:) = 0._wp227 zdit (ntsi-nn_hls,:,:) = 0._wp ; zdit (ntei+nn_hls,:,:) = 0._wp 228 zdjt (ntsi-nn_hls,:,:) = 0._wp ; zdjt (ntei+nn_hls,:,:) = 0._wp 199 229 !!end 200 230 … … 204 234 zdjt(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 205 235 END_3D 236 ! TODO: NOT TESTED- requires zps 206 237 IF( ln_zps ) THEN ! botton and surface ocean correction of the horizontal gradient 207 238 DO_2D( 1, 0, 1, 0 ) … … 209 240 zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 210 241 END_2D 242 ! TODO: NOT TESTED- requires isf 211 243 IF( ln_isfcav ) THEN ! first wet level beneath a cavity 212 244 DO_2D( 1, 0, 1, 0 ) … … 223 255 DO jk = 1, jpkm1 ! Horizontal slab 224 256 ! 225 ! !== Vertical tracer gradient 226 zdk1t(:,:) = ( pt(:,:,jk,jn) - pt(:,:,jk+1,jn) ) * wmask(:,:,jk+1) ! level jk+1 227 ! 228 IF( jk == 1 ) THEN ; zdkt(:,:) = zdk1t(:,:) ! surface: zdkt(jk=1)=zdkt(jk=2) 229 ELSE ; zdkt(:,:) = ( pt(:,:,jk-1,jn) - pt(:,:,jk,jn) ) * wmask(:,:,jk) 230 ENDIF 257 DO_2D( 1, 1, 1, 1 ) 258 ! !== Vertical tracer gradient 259 zdk1t(ji,jj) = ( pt(ji,jj,jk,jn) - pt(ji,jj,jk+1,jn) ) * wmask(ji,jj,jk+1) ! level jk+1 260 ! 261 IF( jk == 1 ) THEN ; zdkt(ji,jj) = zdk1t(ji,jj) ! surface: zdkt(jk=1)=zdkt(jk=2) 262 ELSE ; zdkt(ji,jj) = ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) * wmask(ji,jj,jk) 263 ENDIF 264 END_2D 265 ! 231 266 DO_2D( 1, 0, 1, 0 ) 232 267 zabe1 = pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm) … … 330 365 END DO ! end tracer loop 331 366 ! 332 END SUBROUTINE tra_ldf_iso 367 END SUBROUTINE tra_ldf_iso_t 333 368 334 369 !!==============================================================================
Note: See TracChangeset
for help on using the changeset viewer.