- Timestamp:
- 2010-12-06T06:19:30+01:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/TRA/traldf_iso_grif.F90
r2450 r2454 11 11 !! 'key_ldfslp' slope of the lateral diffusive direction 12 12 !!---------------------------------------------------------------------- 13 !! tra_ldf_iso_grif : update the tracer trend with the horizontal 14 !! component of a iso-neutral laplacian operator 15 !! and with the vertical part of 16 !! the isopycnal or geopotential s-coord. operator 13 !! tra_ldf_iso_grif : update the tracer trend with the horizontal component 14 !! of the Griffies iso-neutral laplacian operator 17 15 !!---------------------------------------------------------------------- 18 16 USE oce ! ocean dynamics and active tracers 19 17 USE dom_oce ! ocean space and time domain 18 USE trc_oce ! share passive tracers/Ocean variables 19 USE zdf_oce ! ocean vertical physics 20 20 USE ldftra_oce ! ocean active tracers: lateral physics 21 USE zdf_oce ! ocean vertical physics22 USE in_out_manager ! I/O manager23 USE iom !24 21 USE ldfslp ! iso-neutral slopes 25 22 USE diaptr ! poleward transport diagnostics 26 USE trc_oce ! share passive tracers/Ocean variables 23 USE in_out_manager ! I/O manager 24 USE iom ! I/O library 27 25 #if defined key_diaar5 28 26 USE phycst ! physical constants … … 112 110 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdit, zdjt, ztfw ! 3D workspace 113 111 ! 114 REAL(wp) :: zslope_skew, zslope_iso,zslope2, zbu, zbv115 REAL(wp) :: ze1ur, zdxt,ze2vr,ze3wr,zdyt,zdzt116 REAL(wp) :: ah,zah_slp,zaei_slp112 REAL(wp) :: zslope_skew, zslope_iso, zslope2, zbu, zbv 113 REAL(wp) :: ze1ur, zdxt, ze2vr, ze3wr, zdyt, zdzt 114 REAL(wp) :: zah, zah_slp, zaei_slp 117 115 #if defined key_diaar5 118 116 REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace … … 144 142 !!gm Future development: consider using Ah defined at T-points and attached to the 4 t-point triads 145 143 146 ah_wslp2(:,:,:) = 0. 0144 ah_wslp2(:,:,:) = 0._wp 147 145 IF( ln_traldf_gdia ) THEN 148 psix_eiv(:,:,:) = 0. 0149 psiy_eiv(:,:,:) = 0. 0146 psix_eiv(:,:,:) = 0._wp 147 psiy_eiv(:,:,:) = 0._wp 150 148 ENDIF 151 149 152 DO ip =0,1153 DO kp =0,1154 DO jk =1,jpkm1155 DO jj =1,jpjm1156 DO ji =1,fs_jpim1157 ze3wr =1.0_wp/fse3w(ji+ip,jj,jk+kp)158 zbu = 0.25_wp*e1u(ji,jj)*e2u(ji,jj)*fse3u(ji,jj,jk)159 ah= fsahtu(ji,jj,jk) ! fsaht(ji+ip,jj,jk)150 DO ip = 0, 1 151 DO kp = 0, 1 152 DO jk = 1, jpkm1 153 DO jj = 1, jpjm1 154 DO ji = 1, fs_jpim1 155 ze3wr = 1._wp / fse3w(ji+ip,jj,jk+kp) 156 zbu = 0.25_wp * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) 157 zah = fsahtu(ji,jj,jk) ! fsaht(ji+ip,jj,jk) 160 158 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 161 zslope2 = ( zslope_skew &162 & - umask(ji,jj,jk+kp)*(fsdept(ji+1,jj,jk) - fsdept(ji ,jj ,jk))*ze1ur )**2163 ah_wslp2(ji+ip,jj,jk+kp) =ah_wslp2(ji+ip,jj,jk+kp) +&164 & ah*((zbu*ze3wr)/(e1t(ji+ip,jj)*e2t(ji+ip,jj)))*zslope2159 zslope2 = zslope_skew - ( fsdept(ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * ze1ur * umask(ji,jj,jk+kp) 160 zslope2 = zslope2 *zslope2 161 ah_wslp2(ji+ip,jj,jk+kp) = ah_wslp2(ji+ip,jj,jk+kp) & 162 & + ah * ( zbu * ze3wr / ( e1t(ji+ip,jj) * e2t(ji+ip,jj) ) ) * zslope2 165 163 IF( ln_traldf_gdia ) THEN 166 zaei_slp = fsaeiw(ji+ip,jj,jk) *zslope_skew!fsaeit(ji+ip,jj,jk)*zslope_skew167 psix_eiv(ji,jj,jk+kp) = psix_eiv(ji,jj,jk+kp) + 0.25_wp *zaei_slp164 zaei_slp = fsaeiw(ji+ip,jj,jk) * zslope_skew !fsaeit(ji+ip,jj,jk)*zslope_skew 165 psix_eiv(ji,jj,jk+kp) = psix_eiv(ji,jj,jk+kp) + 0.25_wp * zaei_slp 168 166 ENDIF 169 167 END DO … … 173 171 END DO 174 172 ! 175 DO jp =0,1176 DO kp =0,1177 DO jk =1,jpkm1178 DO jj =1,jpjm1173 DO jp = 0, 1 174 DO kp = 0, 1 175 DO jk = 1, jpkm1 176 DO jj = 1, jpjm1 179 177 DO ji=1,fs_jpim1 180 ze3wr =1.0_wp/fse3w(ji,jj+jp,jk+kp)181 zbv = 0.25_wp*e1v(ji,jj)*e2v(ji,jj)*fse3v(ji,jj,jk)182 ah = fsahtu(ji,jj,jk)!fsaht(ji,jj+jp,jk)178 ze3wr = 1.0_wp / fse3w(ji,jj+jp,jk+kp) 179 zbv = 0.25_wp * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) 180 zah = fsahtu(ji,jj,jk) !fsaht(ji,jj+jp,jk) 183 181 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 184 zslope2 =(zslope_skew - vmask(ji,jj,jk+kp)*(fsdept(ji,jj+1,jk) - fsdept(ji ,jj ,jk))*ze2vr &185 & )**2186 ah_wslp2(ji,jj+jp,jk+kp) =ah_wslp2(ji,jj+jp,jk+kp) +&187 & ah*((zbv*ze3wr)/(e1t(ji,jj+jp)*e2t(ji,jj+jp)))*zslope2182 zslope2 = zslope_skew - ( fsdept(ji,jj+1,jk) - fsdept(ji,jj,jk) ) * ze2vr * vmask(ji,jj,jk+kp) 183 zslope2 = zslope2 * zslope2 184 ah_wslp2(ji,jj+jp,jk+kp) = ah_wslp2(ji,jj+jp,jk+kp) & 185 & + zah * ( zbv * ze3wr / ( e1t(ji,jj+jp) * e2t(ji,jj+jp) ) ) * zslope2 188 186 IF( ln_traldf_gdia ) THEN 189 187 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew !fsaeit(ji,jj+jp,jk)*zslope_skew 190 psiy_eiv(ji,jj,jk+kp) = psiy_eiv(ji,jj,jk+kp) + 0.25_wp *zaei_slp188 psiy_eiv(ji,jj,jk+kp) = psiy_eiv(ji,jj,jk+kp) + 0.25_wp * zaei_slp 191 189 ENDIF 192 190 END DO … … 245 243 ENDIF 246 244 247 DO ip =0,1 !== Horizontal & vertical fluxes248 DO kp =0,1249 DO jj =1,jpjm1250 DO ji =1,fs_jpim1245 DO ip = 0, 1 !== Horizontal & vertical fluxes 246 DO kp = 0, 1 247 DO jj = 1, jpjm1 248 DO ji = 1, fs_jpim1 251 249 ze1ur = 1._wp / e1u(ji,jj) 252 250 zdxt = zdit(ji,jj,jk) * ze1ur 253 ze3wr = 1._wp /fse3w(ji+ip,jj,jk+kp)251 ze3wr = 1._wp / fse3w(ji+ip,jj,jk+kp) 254 252 zdzt = zdkt(ji+ip,jj,kp) * ze3wr 255 253 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) … … 257 255 258 256 zbu = 0.25_wp * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) 259 ah = fsahtu(ji,jj,jk)!*umask(ji,jj,jk+kp) !fsaht(ji+ip,jj,jk)===>> ????260 zah_slp = ah*zslope_iso257 zah = fsahtu(ji,jj,jk) !*umask(ji,jj,jk+kp) !fsaht(ji+ip,jj,jk) ===>> ???? 258 zah_slp = zah * zslope_iso 261 259 zaei_slp = fsaeiw(ji+ip,jj,jk) * zslope_skew !fsaeit(ji+ip,jj,jk)*zslope_skew 262 zftu(ji,jj,jk) = zftu(ji,jj,jk) - (ah*zdxt + (zah_slp - zaei_slp)*zdzt)*zbu*ze1ur263 ztfw(ji+ip,jj,jk+kp) = ztfw(ji+ip,jj,jk+kp) - (zah_slp + zaei_slp) *zdxt *zbu*ze3wr260 zftu(ji,jj,jk) = zftu(ji,jj,jk) - ( zah * zdxt + (zah_slp - zaei_slp) * zdzt ) * zbu * ze1ur 261 ztfw(ji+ip,jj,jk+kp) = ztfw(ji+ip,jj,jk+kp) - (zah_slp + zaei_slp) * zdxt * zbu * ze3wr 264 262 END DO 265 263 END DO … … 267 265 END DO 268 266 269 DO jp =0,1270 DO kp =0,1271 DO jj =1,jpjm1272 DO ji =1,fs_jpim1273 ze2vr = 1._wp /e2v(ji,jj)274 zdyt = zdjt(ji,jj,jk) *ze2vr275 ze3wr = 1._wp /fse3w(ji,jj+jp,jk+kp)267 DO jp = 0, 1 268 DO kp = 0, 1 269 DO jj = 1, jpjm1 270 DO ji = 1, fs_jpim1 271 ze2vr = 1._wp / e2v(ji,jj) 272 zdyt = zdjt(ji,jj,jk) * ze2vr 273 ze3wr = 1._wp / fse3w(ji,jj+jp,jk+kp) 276 274 zdzt = zdkt(ji,jj+jp,kp) * ze3wr 277 275 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 278 zslope_iso = triadj(ji,jj+jp,jk,1-jp,kp)276 zslope_iso = triadj(ji,jj+jp,jk,1-jp,kp) 279 277 zbv = 0.25_wp * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) 280 ah = fsahtv(ji,jj,jk)!*vmask(ji,jj,jk+kp) !fsaht(ji,jj+jp,jk)281 zah_slp = ah * zslope_iso282 zaei_slp = fsaeiw(ji,jj+jp,jk) *zslope_skew!fsaeit(ji,jj+jp,jk)*zslope_skew283 zftv(ji,jj,jk) = zftv(ji,jj,jk) - ( ah*zdyt + (zah_slp - zaei_slp)*zdzt)*zbv*ze2vr284 ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - (zah_slp + zaei_slp) *zdyt*zbv*ze3wr278 zah = fsahtv(ji,jj,jk) !*vmask(ji,jj,jk+kp) !fsaht(ji,jj+jp,jk) 279 zah_slp = zah * zslope_iso 280 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew !fsaeit(ji,jj+jp,jk)*zslope_skew 281 zftv(ji,jj,jk) = zftv(ji,jj,jk) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr 282 ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - (zah_slp + zaei_slp) * zdyt * zbv * ze3wr 285 283 END DO 286 284 END DO … … 350 348 !!---------------------------------------------------------------------- 351 349 CONTAINS 352 SUBROUTINE tra_ldf_iso_grif( kt ) ! Empty routine 353 WRITE(*,*) 'tra_ldf_iso: You should not have seen this print! error?', kt 350 SUBROUTINE tra_ldf_iso_grif( kt, cdtype, pgu, pgv, ptb, pta, kjpt, pahtb0 ) ! Empty routine 351 CHARACTER(len=3) :: cdtype 352 REAL, DIMENSION(:,:,:) :: pgu, pgv ! tracer gradient at pstep levels 353 REAL, DIMENSION(:,:,:,:) :: ptb, pta 354 WRITE(*,*) 'tra_ldf_iso_grif: You should not have seen this print! error?', kt, cdtype, & 355 & pgu(1,1,1), pgv(1,1,1), ptb(1,1,1,1), pta(1,1,1,1), kjpt, pahtb0 354 356 END SUBROUTINE tra_ldf_iso_grif 355 357 #endif
Note: See TracChangeset
for help on using the changeset viewer.
