Changeset 6140 for trunk/NEMOGCM/NEMO/OPA_SRC/TRA/zpshde.F90
- Timestamp:
- 2015-12-21T12:35:23+01:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/TRA/zpshde.F90
r5836 r6140 32 32 33 33 !! * Substitutions 34 # include "domzgr_substitute.h90"35 34 # include "vectopt_loop_substitute.h90" 36 35 !!---------------------------------------------------------------------- … … 111 110 iku = mbku(ji,jj) ; ikum1 = MAX( iku - 1 , 1 ) ! last and before last ocean level at u- & v-points 112 111 ikv = mbkv(ji,jj) ; ikvm1 = MAX( ikv - 1 , 1 ) ! if level first is a p-step, ik.m1=1 113 ze3wu = fse3w(ji+1,jj ,iku) - fse3w(ji,jj,iku) 114 ze3wv = fse3w(ji ,jj+1,ikv) - fse3w(ji,jj,ikv) 112 !!gm BUG ? when applied to before fields, e3w_b should be used.... 113 ze3wu = e3w_n(ji+1,jj ,iku) - e3w_n(ji,jj,iku) 114 ze3wv = e3w_n(ji ,jj+1,ikv) - e3w_n(ji,jj,ikv) 115 115 ! 116 116 ! i- direction 117 117 IF( ze3wu >= 0._wp ) THEN ! case 1 118 zmaxu = ze3wu / fse3w(ji+1,jj,iku)118 zmaxu = ze3wu / e3w_n(ji+1,jj,iku) 119 119 ! interpolated values of tracers 120 120 zti (ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikum1,jn) - pta(ji+1,jj,iku,jn) ) … … 122 122 pgtu(ji,jj,jn) = umask(ji,jj,1) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) ) 123 123 ELSE ! case 2 124 zmaxu = -ze3wu / fse3w(ji,jj,iku)124 zmaxu = -ze3wu / e3w_n(ji,jj,iku) 125 125 ! interpolated values of tracers 126 126 zti (ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikum1,jn) - pta(ji,jj,iku,jn) ) … … 131 131 ! j- direction 132 132 IF( ze3wv >= 0._wp ) THEN ! case 1 133 zmaxv = ze3wv / fse3w(ji,jj+1,ikv)133 zmaxv = ze3wv / e3w_n(ji,jj+1,ikv) 134 134 ! interpolated values of tracers 135 135 ztj (ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvm1,jn) - pta(ji,jj+1,ikv,jn) ) … … 137 137 pgtv(ji,jj,jn) = vmask(ji,jj,1) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) ) 138 138 ELSE ! case 2 139 zmaxv = -ze3wv / fse3w(ji,jj,ikv)139 zmaxv = -ze3wv / e3w_n(ji,jj,ikv) 140 140 ! interpolated values of tracers 141 141 ztj (ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvm1,jn) - pta(ji,jj,ikv,jn) ) … … 156 156 iku = mbku(ji,jj) 157 157 ikv = mbkv(ji,jj) 158 ze3wu = fse3w(ji+1,jj ,iku) - fse3w(ji,jj,iku)159 ze3wv = fse3w(ji ,jj+1,ikv) - fse3w(ji,jj,ikv)160 IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = fsdept(ji ,jj,iku) ! i-direction: case 1161 ELSE ; zhi(ji,jj) = fsdept(ji+1,jj,iku) ! - - case 2162 ENDIF 163 IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = fsdept(ji,jj ,ikv) ! j-direction: case 1164 ELSE ; zhj(ji,jj) = fsdept(ji,jj+1,ikv) ! - - case 2158 ze3wu = e3w_n(ji+1,jj ,iku) - e3w_n(ji,jj,iku) 159 ze3wv = e3w_n(ji ,jj+1,ikv) - e3w_n(ji,jj,ikv) 160 IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = gdept_n(ji ,jj,iku) ! i-direction: case 1 161 ELSE ; zhi(ji,jj) = gdept_n(ji+1,jj,iku) ! - - case 2 162 ENDIF 163 IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = gdept_n(ji,jj ,ikv) ! j-direction: case 1 164 ELSE ; zhj(ji,jj) = gdept_n(ji,jj+1,ikv) ! - - case 2 165 165 ENDIF 166 166 END DO … … 174 174 iku = mbku(ji,jj) 175 175 ikv = mbkv(ji,jj) 176 ze3wu = fse3w(ji+1,jj ,iku) - fse3w(ji,jj,iku)177 ze3wv = fse3w(ji ,jj+1,ikv) - fse3w(ji,jj,ikv)176 ze3wu = e3w_n(ji+1,jj ,iku) - e3w_n(ji,jj,iku) 177 ze3wv = e3w_n(ji ,jj+1,ikv) - e3w_n(ji,jj,ikv) 178 178 IF( ze3wu >= 0._wp ) THEN ; pgru(ji,jj) = umask(ji,jj,1) * ( zri(ji ,jj ) - prd(ji,jj,iku) ) ! i: 1 179 179 ELSE ; pgru(ji,jj) = umask(ji,jj,1) * ( prd(ji+1,jj,iku) - zri(ji,jj ) ) ! i: 2 … … 191 191 ! 192 192 END SUBROUTINE zps_hde 193 194 195 SUBROUTINE zps_hde_isf( kt, kjpt, pta, pgtu , pgtv , pgtui, pgtvi, & 196 & prd, pgru , pgrv , pmru , pmrv , pgzu , pgzv , pge3ru , pge3rv , & 197 & pgrui, pgrvi, pmrui, pmrvi, pgzui, pgzvi, pge3rui, pge3rvi ) 198 !!---------------------------------------------------------------------- 199 !! *** ROUTINE zps_hde *** 193 ! 194 SUBROUTINE zps_hde_isf( kt, kjpt, pta, pgtu, pgtv, pgtui, pgtvi, & 195 & prd, pgru, pgrv, pgrui, pgrvi ) 196 !!---------------------------------------------------------------------- 197 !! *** ROUTINE zps_hde_isf *** 200 198 !! 201 199 !! ** Purpose : Compute the horizontal derivative of T, S and rho 202 200 !! at u- and v-points with a linear interpolation for z-coordinate 203 !! with partial steps .201 !! with partial steps for top (ice shelf) and bottom. 204 202 !! 205 203 !! ** Method : In z-coord with partial steps, scale factors on last 206 204 !! levels are different for each grid point, so that T, S and rd 207 205 !! points are not at the same depth as in z-coord. To have horizontal 208 !! gradients again, we interpolate T and S at the good depth : 206 !! gradients again, we interpolate T and S at the good depth : 207 !! For the bottom case: 209 208 !! Linear interpolation of T, S 210 209 !! Computation of di(tb) and dj(tb) by vertical interpolation: … … 235 234 !! di(rho) = rd~ - rd(i,j,k) or rd(i+1,j,k) - rd~ 236 235 !! 236 !! For the top case (ice shelf): As for the bottom case but upside down 237 !! 237 238 !! ** Action : compute for top and bottom interfaces 238 239 !! - pgtu, pgtv, pgtui, pgtvi: horizontal gradient of tracer at u- & v-points 239 240 !! - pgru, pgrv, pgrui, pgtvi: horizontal gradient of rho (if present) at u- & v-points 240 !! - pmru, pmrv, pmrui, pmrvi: horizontal sum of rho at u- & v- point (used in dynhpg with vvl) 241 !! - pgzu, pgzv, pgzui, pgzvi: horizontal gradient of z at u- and v- point (used in dynhpg with vvl) 242 !! - pge3ru, pge3rv, pge3rui, pge3rvi: horizontal gradient of rho weighted by local e3w at u- & v-points 243 !!---------------------------------------------------------------------- 244 INTEGER , INTENT(in ) :: kt ! ocean time-step index 245 INTEGER , INTENT(in ) :: kjpt ! number of tracers 246 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pta ! 4D tracers fields 247 ! !! u-point ! v-point ! 248 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT( out) :: pgtu , pgtv ! bottom GRADh( ptra ) 249 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT( out) :: pgtui , pgtvi ! top GRADh( ptra ) 250 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ), OPTIONAL :: prd ! 3D density anomaly fields 251 ! !! u-point ! v-point ! 252 REAL(wp), DIMENSION(jpi,jpj ), INTENT( out), OPTIONAL :: pgru , pgrv ! bottom GRADh( prd ) 253 REAL(wp), DIMENSION(jpi,jpj ), INTENT( out), OPTIONAL :: pmru , pmrv ! bottom SUM ( prd ) 254 REAL(wp), DIMENSION(jpi,jpj ), INTENT( out), OPTIONAL :: pgzu , pgzv ! bottom GRADh( z ) 255 REAL(wp), DIMENSION(jpi,jpj ), INTENT( out), OPTIONAL :: pge3ru , pge3rv ! bottom GRADh( prd ) weighted by e3w 256 REAL(wp), DIMENSION(jpi,jpj ), INTENT( out), OPTIONAL :: pgrui , pgrvi ! top GRADh( prd ) 257 REAL(wp), DIMENSION(jpi,jpj ), INTENT( out), OPTIONAL :: pmrui , pmrvi ! top SUM ( prd ) 258 REAL(wp), DIMENSION(jpi,jpj ), INTENT( out), OPTIONAL :: pgzui , pgzvi ! top GRADh( z ) 259 REAL(wp), DIMENSION(jpi,jpj ), INTENT( out), OPTIONAL :: pge3rui , pge3rvi ! top GRADh( prd ) weighted by e3w 241 !!---------------------------------------------------------------------- 242 INTEGER , INTENT(in ) :: kt ! ocean time-step index 243 INTEGER , INTENT(in ) :: kjpt ! number of tracers 244 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pta ! 4D tracers fields 245 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT( out) :: pgtu, pgtv ! hor. grad. of ptra at u- & v-pts 246 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT( out) :: pgtui, pgtvi ! hor. grad. of stra at u- & v-pts (ISF) 247 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ), OPTIONAL :: prd ! 3D density anomaly fields 248 REAL(wp), DIMENSION(jpi,jpj ), INTENT( out), OPTIONAL :: pgru, pgrv ! hor. grad of prd at u- & v-pts (bottom) 249 REAL(wp), DIMENSION(jpi,jpj ), INTENT( out), OPTIONAL :: pgrui, pgrvi ! hor. grad of prd at u- & v-pts (top) 260 250 ! 261 251 INTEGER :: ji, jj, jn ! Dummy loop indices 262 252 INTEGER :: iku, ikv, ikum1, ikvm1,ikup1, ikvp1 ! partial step level (ocean bottom level) at u- and v-points 263 REAL(wp) :: ze3wu, ze3wv, zmaxu, zmaxv , zdzwu, zdzwv, zdzwuip1, zdzwvjp1! temporary scalars253 REAL(wp) :: ze3wu, ze3wv, zmaxu, zmaxv ! temporary scalars 264 254 REAL(wp), DIMENSION(jpi,jpj) :: zri, zrj, zhi, zhj ! NB: 3rd dim=1 to use eos 265 255 REAL(wp), DIMENSION(jpi,jpj,kjpt) :: zti, ztj ! … … 277 267 DO jj = 1, jpjm1 278 268 DO ji = 1, jpim1 279 iku = mbku(ji,jj) ; ikum1 = MAX( iku - 1 , 1 ) ! last and before last ocean level at u- & v-points 280 ikv = mbkv(ji,jj) ; ikvm1 = MAX( ikv - 1 , 1 ) ! if level first is a p-step, ik.m1=1 269 270 iku = mbku(ji,jj); ikum1 = MAX( iku - 1 , 1 ) ! last and before last ocean level at u- & v-points 271 ikv = mbkv(ji,jj); ikvm1 = MAX( ikv - 1 , 1 ) ! if level first is a p-step, ik.m1=1 272 ze3wu = gdept_n(ji+1,jj,iku) - gdept_n(ji,jj,iku) 273 ze3wv = gdept_n(ji,jj+1,ikv) - gdept_n(ji,jj,ikv) 274 ! 275 ! i- direction 276 IF( ze3wu >= 0._wp ) THEN ! case 1 277 zmaxu = ze3wu / e3w_n(ji+1,jj,iku) 278 ! interpolated values of tracers 279 zti (ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikum1,jn) - pta(ji+1,jj,iku,jn) ) 280 ! gradient of tracers 281 pgtu(ji,jj,jn) = ssumask(ji,jj) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) ) 282 ELSE ! case 2 283 zmaxu = -ze3wu / e3w_n(ji,jj,iku) 284 ! interpolated values of tracers 285 zti (ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikum1,jn) - pta(ji,jj,iku,jn) ) 286 ! gradient of tracers 287 pgtu(ji,jj,jn) = ssumask(ji,jj) * ( pta(ji+1,jj,iku,jn) - zti(ji,jj,jn) ) 288 ENDIF 289 ! 290 ! j- direction 291 IF( ze3wv >= 0._wp ) THEN ! case 1 292 zmaxv = ze3wv / e3w_n(ji,jj+1,ikv) 293 ! interpolated values of tracers 294 ztj (ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvm1,jn) - pta(ji,jj+1,ikv,jn) ) 295 ! gradient of tracers 296 pgtv(ji,jj,jn) = ssvmask(ji,jj) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) ) 297 ELSE ! case 2 298 zmaxv = -ze3wv / e3w_n(ji,jj,ikv) 299 ! interpolated values of tracers 300 ztj (ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvm1,jn) - pta(ji,jj,ikv,jn) ) 301 ! gradient of tracers 302 pgtv(ji,jj,jn) = ssvmask(ji,jj) * ( pta(ji,jj+1,ikv,jn) - ztj(ji,jj,jn) ) 303 ENDIF 304 305 END DO 306 END DO 307 CALL lbc_lnk( pgtu(:,:,jn), 'U', -1. ) ; CALL lbc_lnk( pgtv(:,:,jn), 'V', -1. ) ! Lateral boundary cond. 308 ! 309 END DO 310 311 ! horizontal derivative of density anomalies (rd) 312 IF( PRESENT( prd ) ) THEN ! depth of the partial step level 313 pgru(:,:)=0.0_wp ; pgrv(:,:)=0.0_wp ; 314 ! 315 DO jj = 1, jpjm1 316 DO ji = 1, jpim1 317 318 iku = mbku(ji,jj) 319 ikv = mbkv(ji,jj) 320 ze3wu = gdept_n(ji+1,jj,iku) - gdept_n(ji,jj,iku) 321 ze3wv = gdept_n(ji,jj+1,ikv) - gdept_n(ji,jj,ikv) 322 ! 323 IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = gdept_n(ji ,jj,iku) ! i-direction: case 1 324 ELSE ; zhi(ji,jj) = gdept_n(ji+1,jj,iku) ! - - case 2 325 ENDIF 326 IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = gdept_n(ji,jj ,ikv) ! j-direction: case 1 327 ELSE ; zhj(ji,jj) = gdept_n(ji,jj+1,ikv) ! - - case 2 328 ENDIF 329 330 END DO 331 END DO 332 333 ! Compute interpolated rd from zti, ztj for the 2 cases at the depth of the partial 334 ! step and store it in zri, zrj for each case 335 CALL eos( zti, zhi, zri ) 336 CALL eos( ztj, zhj, zrj ) 337 338 DO jj = 1, jpjm1 ! Gradient of density at the last level 339 DO ji = 1, jpim1 340 iku = mbku(ji,jj) 341 ikv = mbkv(ji,jj) 342 ze3wu = gdept_n(ji+1,jj,iku) - gdept_n(ji,jj,iku) 343 ze3wv = gdept_n(ji,jj+1,ikv) - gdept_n(ji,jj,ikv) 344 345 IF( ze3wu >= 0._wp ) THEN ; pgru(ji,jj) = ssumask(ji,jj) * ( zri(ji ,jj ) - prd(ji,jj,iku) ) ! i: 1 346 ELSE ; pgru(ji,jj) = ssumask(ji,jj) * ( prd(ji+1,jj,iku) - zri(ji,jj ) ) ! i: 2 347 ENDIF 348 IF( ze3wv >= 0._wp ) THEN ; pgrv(ji,jj) = ssvmask(ji,jj) * ( zrj(ji,jj ) - prd(ji,jj,ikv) ) ! j: 1 349 ELSE ; pgrv(ji,jj) = ssvmask(ji,jj) * ( prd(ji,jj+1,ikv) - zrj(ji,jj ) ) ! j: 2 350 ENDIF 351 352 END DO 353 END DO 354 355 CALL lbc_lnk( pgru , 'U', -1. ) ; CALL lbc_lnk( pgrv , 'V', -1. ) ! Lateral boundary conditions 356 ! 357 END IF 358 ! 359 ! !== (ISH) compute grui and gruvi ==! 360 ! 361 DO jn = 1, kjpt !== Interpolation of tracers at the last ocean level ==! ! 362 DO jj = 1, jpjm1 363 DO ji = 1, jpim1 364 iku = miku(ji,jj); ikup1 = miku(ji,jj) + 1 365 ikv = mikv(ji,jj); ikvp1 = mikv(ji,jj) + 1 366 ! 281 367 ! (ISF) case partial step top and bottom in adjacent cell in vertical 282 368 ! cannot used e3w because if 2 cell water column, we have ps at top and bottom 283 369 ! in this case e3w(i,j) - e3w(i,j+1) is not the distance between Tj~ and Tj 284 370 ! the only common depth between cells (i,j) and (i,j+1) is gdepw_0 285 ze3wu = (gdept_0(ji+1,jj,iku) - gdepw_0(ji+1,jj,iku)) - (gdept_0(ji,jj,iku) - gdepw_0(ji,jj,iku))286 ze3wv = (gdept_0(ji,jj+1,ikv) - gdepw_0(ji,jj+1,ikv)) - (gdept_0(ji,jj,ikv) - gdepw_0(ji,jj,ikv))287 ! 371 ze3wu = gdept_n(ji,jj,iku) - gdept_n(ji+1,jj,iku) 372 ze3wv = gdept_n(ji,jj,ikv) - gdept_n(ji,jj+1,ikv) 373 288 374 ! i- direction 289 375 IF( ze3wu >= 0._wp ) THEN ! case 1 290 zmaxu = ze3wu / fse3w(ji+1,jj,iku) 291 ! interpolated values of tracers 292 zti (ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikum1,jn) - pta(ji+1,jj,iku,jn) ) 376 zmaxu = ze3wu / e3w_n(ji+1,jj,ikup1) 377 ! interpolated values of tracers 378 zti(ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikup1,jn) - pta(ji+1,jj,iku,jn) ) 379 ! gradient of tracers 380 pgtui(ji,jj,jn) = ssumask(ji,jj) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) ) 381 ELSE ! case 2 382 zmaxu = - ze3wu / e3w_n(ji,jj,ikup1) 383 ! interpolated values of tracers 384 zti(ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikup1,jn) - pta(ji,jj,iku,jn) ) 293 385 ! gradient of tracers 294 pgtu(ji,jj,jn) = umask(ji,jj,iku) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) ) 295 ELSE ! case 2 296 zmaxu = -ze3wu / fse3w(ji,jj,iku) 297 ! interpolated values of tracers 298 zti (ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikum1,jn) - pta(ji,jj,iku,jn) ) 299 ! gradient of tracers 300 pgtu(ji,jj,jn) = umask(ji,jj,iku) * ( pta(ji+1,jj,iku,jn) - zti(ji,jj,jn) ) 386 pgtui(ji,jj,jn) = ssumask(ji,jj) * ( pta(ji+1,jj,iku,jn) - zti(ji,jj,jn) ) 301 387 ENDIF 302 388 ! 303 389 ! j- direction 304 390 IF( ze3wv >= 0._wp ) THEN ! case 1 305 zmaxv = ze3wv / fse3w(ji,jj+1,ikv) 306 ! interpolated values of tracers 307 ztj (ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvm1,jn) - pta(ji,jj+1,ikv,jn) ) 308 ! gradient of tracers 309 pgtv(ji,jj,jn) = vmask(ji,jj,ikv) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) ) 310 ELSE ! case 2 311 zmaxv = -ze3wv / fse3w(ji,jj,ikv) 312 ! interpolated values of tracers 313 ztj (ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvm1,jn) - pta(ji,jj,ikv,jn) ) 314 ! gradient of tracers 315 pgtv(ji,jj,jn) = vmask(ji,jj,ikv) * ( pta(ji,jj+1,ikv,jn) - ztj(ji,jj,jn) ) 316 ENDIF 317 END DO 318 END DO 319 CALL lbc_lnk( pgtu(:,:,jn), 'U', -1. ) ; CALL lbc_lnk( pgtv(:,:,jn), 'V', -1. ) ! Lateral boundary cond. 391 zmaxv = ze3wv / e3w_n(ji,jj+1,ikvp1) 392 ! interpolated values of tracers 393 ztj(ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvp1,jn) - pta(ji,jj+1,ikv,jn) ) 394 ! gradient of tracers 395 pgtvi(ji,jj,jn) = ssvmask(ji,jj) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) ) 396 ELSE ! case 2 397 zmaxv = - ze3wv / e3w_n(ji,jj,ikvp1) 398 ! interpolated values of tracers 399 ztj(ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvp1,jn) - pta(ji,jj,ikv,jn) ) 400 ! gradient of tracers 401 pgtvi(ji,jj,jn) = ssvmask(ji,jj) * ( pta(ji,jj+1,ikv,jn) - ztj(ji,jj,jn) ) 402 ENDIF 403 404 END DO 405 END DO 406 CALL lbc_lnk( pgtui(:,:,jn), 'U', -1. ); CALL lbc_lnk( pgtvi(:,:,jn), 'V', -1. ) ! Lateral boundary cond. 320 407 ! 321 408 END DO … … 323 410 IF( PRESENT( prd ) ) THEN !== horizontal derivative of density anomalies (rd) ==! (optional part) 324 411 ! 325 pgru (:,:)=0._wp ; pgrv (:,:) = 0._wp 326 pgzu (:,:)=0._wp ; pgzv (:,:) = 0._wp 327 pmru (:,:)=0._wp ; pmru (:,:) = 0._wp 328 pge3ru(:,:)=0._wp ; pge3rv(:,:) = 0._wp 329 ! 330 DO jj = 1, jpjm1 ! depth of the partial step level 331 DO ji = 1, jpim1 332 iku = mbku(ji,jj) 333 ikv = mbkv(ji,jj) 334 ze3wu = (gdept_0(ji+1,jj,iku) - gdepw_0(ji+1,jj,iku)) - (gdept_0(ji,jj,iku) - gdepw_0(ji,jj,iku)) 335 ze3wv = (gdept_0(ji,jj+1,ikv) - gdepw_0(ji,jj+1,ikv)) - (gdept_0(ji,jj,ikv) - gdepw_0(ji,jj,ikv)) 336 ! 337 IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = fsdept(ji+1,jj,iku) - ze3wu ! i-direction: case 1 338 ELSE ; zhi(ji,jj) = fsdept(ji ,jj,iku) + ze3wu ! - - case 2 339 ENDIF 340 IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = fsdept(ji,jj+1,ikv) - ze3wv ! j-direction: case 1 341 ELSE ; zhj(ji,jj) = fsdept(ji,jj ,ikv) + ze3wv ! - - case 2 342 ENDIF 412 pgrui(:,:) =0.0_wp; pgrvi(:,:) =0.0_wp; 413 DO jj = 1, jpjm1 414 DO ji = 1, jpim1 415 416 iku = miku(ji,jj) 417 ikv = mikv(ji,jj) 418 ze3wu = gdept_n(ji,jj,iku) - gdept_n(ji+1,jj,iku) 419 ze3wv = gdept_n(ji,jj,ikv) - gdept_n(ji,jj+1,ikv) 420 ! 421 IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = gdept_n(ji ,jj,iku) ! i-direction: case 1 422 ELSE ; zhi(ji,jj) = gdept_n(ji+1,jj,iku) ! - - case 2 423 ENDIF 424 425 IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = gdept_n(ji,jj ,ikv) ! j-direction: case 1 426 ELSE ; zhj(ji,jj) = gdept_n(ji,jj+1,ikv) ! - - case 2 427 ENDIF 428 343 429 END DO 344 430 END DO … … 346 432 CALL eos( zti, zhi, zri ) ! interpolated density from zti, ztj 347 433 CALL eos( ztj, zhj, zrj ) ! at the partial step depth output in zri, zrj 348 434 ! 349 435 DO jj = 1, jpjm1 ! Gradient of density at the last level 350 436 DO ji = 1, jpim1 351 iku = mbku(ji,jj) ; ikum1 = MAX( iku - 1 , 1 ) ! last and before last ocean level at u- & v-points 352 ikv = mbkv(ji,jj) ; ikvm1 = MAX( ikv - 1 , 1 ) ! last and before last ocean level at u- & v-points 353 ze3wu = (gdept_0(ji+1,jj,iku) - gdepw_0(ji+1,jj,iku)) - (gdept_0(ji,jj,iku) - gdepw_0(ji,jj,iku)) 354 ze3wv = (gdept_0(ji,jj+1,ikv) - gdepw_0(ji,jj+1,ikv)) - (gdept_0(ji,jj,ikv) - gdepw_0(ji,jj,ikv)) 355 IF( ze3wu >= 0._wp ) THEN 356 pgzu(ji,jj) = (fsde3w(ji+1,jj,iku) - ze3wu) - fsde3w(ji,jj,iku) 357 pgru(ji,jj) = umask(ji,jj,iku) * ( zri(ji ,jj) - prd(ji,jj,iku) ) ! i: 1 358 pmru(ji,jj) = umask(ji,jj,iku) * ( zri(ji ,jj) + prd(ji,jj,iku) ) ! i: 1 359 pge3ru(ji,jj) = umask(ji,jj,iku) & 360 * ( (fse3w(ji+1,jj,iku) - ze3wu )* ( zri(ji ,jj ) + prd(ji+1,jj,ikum1) + 2._wp) & 361 - fse3w(ji ,jj,iku) * ( prd(ji ,jj,iku) + prd(ji ,jj,ikum1) + 2._wp) ) ! j: 2 362 ELSE 363 pgzu(ji,jj) = fsde3w(ji+1,jj,iku) - (fsde3w(ji,jj,iku) + ze3wu) 364 pgru(ji,jj) = umask(ji,jj,iku) * ( prd(ji+1,jj,iku) - zri(ji,jj) ) ! i: 2 365 pmru(ji,jj) = umask(ji,jj,iku) * ( prd(ji+1,jj,iku) + zri(ji,jj) ) ! i: 2 366 pge3ru(ji,jj) = umask(ji,jj,iku) & 367 * ( fse3w(ji+1,jj,iku) * ( prd(ji+1,jj,iku) + prd(ji+1,jj,ikum1) + 2._wp) & 368 -(fse3w(ji ,jj,iku) + ze3wu) * ( zri(ji ,jj ) + prd(ji ,jj,ikum1) + 2._wp) ) ! j: 2 369 ENDIF 370 IF( ze3wv >= 0._wp ) THEN 371 pgzv(ji,jj) = (fsde3w(ji,jj+1,ikv) - ze3wv) - fsde3w(ji,jj,ikv) 372 pgrv(ji,jj) = vmask(ji,jj,ikv) * ( zrj(ji,jj ) - prd(ji,jj,ikv) ) ! j: 1 373 pmrv(ji,jj) = vmask(ji,jj,ikv) * ( zrj(ji,jj ) + prd(ji,jj,ikv) ) ! j: 1 374 pge3rv(ji,jj) = vmask(ji,jj,ikv) & 375 * ( (fse3w(ji,jj+1,ikv) - ze3wv )* ( zrj(ji,jj ) + prd(ji,jj+1,ikvm1) + 2._wp) & 376 - fse3w(ji,jj ,ikv) * ( prd(ji,jj ,ikv) + prd(ji,jj ,ikvm1) + 2._wp) ) ! j: 2 377 ELSE 378 pgzv(ji,jj) = fsde3w(ji,jj+1,ikv) - (fsde3w(ji,jj,ikv) + ze3wv) 379 pgrv(ji,jj) = vmask(ji,jj,ikv) * ( prd(ji,jj+1,ikv) - zrj(ji,jj) ) ! j: 2 380 pmrv(ji,jj) = vmask(ji,jj,ikv) * ( prd(ji,jj+1,ikv) + zrj(ji,jj) ) ! j: 2 381 pge3rv(ji,jj) = vmask(ji,jj,ikv) & 382 * ( fse3w(ji,jj+1,ikv) * ( prd(ji,jj+1,ikv) + prd(ji,jj+1,ikvm1) + 2._wp) & 383 -(fse3w(ji,jj ,ikv) + ze3wv) * ( zrj(ji,jj ) + prd(ji,jj ,ikvm1) + 2._wp) ) ! j: 2 384 ENDIF 385 END DO 386 END DO 387 CALL lbc_lnk( pgru , 'U', -1. ) ; CALL lbc_lnk( pgrv , 'V', -1. ) ! Lateral boundary conditions 388 CALL lbc_lnk( pmru , 'U', 1. ) ; CALL lbc_lnk( pmrv , 'V', 1. ) ! Lateral boundary conditions 389 CALL lbc_lnk( pgzu , 'U', -1. ) ; CALL lbc_lnk( pgzv , 'V', -1. ) ! Lateral boundary conditions 390 CALL lbc_lnk( pge3ru , 'U', -1. ) ; CALL lbc_lnk( pge3rv , 'V', -1. ) ! Lateral boundary conditions 391 ! 392 END IF 393 ! 394 ! !== (ISH) compute grui and gruvi ==! 395 ! 396 DO jn = 1, kjpt !== Interpolation of tracers at the last ocean level ==! ! 397 DO jj = 1, jpjm1 398 DO ji = 1, jpim1 399 iku = miku(ji,jj) ; ikup1 = miku(ji,jj) + 1 400 ikv = mikv(ji,jj) ; ikvp1 = mikv(ji,jj) + 1 401 ! 402 ! (ISF) case partial step top and bottom in adjacent cell in vertical 403 ! cannot used e3w because if 2 cell water column, we have ps at top and bottom 404 ! in this case e3w(i,j) - e3w(i,j+1) is not the distance between Tj~ and Tj 405 ! the only common depth between cells (i,j) and (i,j+1) is gdepw_0 406 ze3wu = (gdepw_0(ji+1,jj,iku+1) - gdept_0(ji+1,jj,iku)) - (gdepw_0(ji,jj,iku+1) - gdept_0(ji,jj,iku)) 407 ze3wv = (gdepw_0(ji,jj+1,ikv+1) - gdept_0(ji,jj+1,ikv)) - (gdepw_0(ji,jj,ikv+1) - gdept_0(ji,jj,ikv)) 408 ! i- direction 409 IF( ze3wu >= 0._wp ) THEN ! case 1 410 zmaxu = ze3wu / fse3w(ji+1,jj,iku+1) 411 ! interpolated values of tracers 412 zti(ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,iku+1,jn) - pta(ji+1,jj,iku,jn) ) 413 ! gradient of tracers 414 pgtui(ji,jj,jn) = umask(ji,jj,iku) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) ) 415 ELSE ! case 2 416 zmaxu = - ze3wu / fse3w(ji,jj,iku+1) 417 ! interpolated values of tracers 418 zti(ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,iku+1,jn) - pta(ji,jj,iku,jn) ) 419 ! gradient of tracers 420 pgtui(ji,jj,jn) = umask(ji,jj,iku) * ( pta(ji+1,jj,iku,jn) - zti(ji,jj,jn) ) 421 ENDIF 422 ! 423 ! j- direction 424 IF( ze3wv >= 0._wp ) THEN ! case 1 425 zmaxv = ze3wv / fse3w(ji,jj+1,ikv+1) 426 ! interpolated values of tracers 427 ztj(ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikv+1,jn) - pta(ji,jj+1,ikv,jn) ) 428 ! gradient of tracers 429 pgtvi(ji,jj,jn) = vmask(ji,jj,ikv) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) ) 430 ELSE ! case 2 431 zmaxv = - ze3wv / fse3w(ji,jj,ikv+1) 432 ! interpolated values of tracers 433 ztj(ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikv+1,jn) - pta(ji,jj,ikv,jn) ) 434 ! gradient of tracers 435 pgtvi(ji,jj,jn) = vmask(ji,jj,ikv) * ( pta(ji,jj+1,ikv,jn) - ztj(ji,jj,jn) ) 436 ENDIF 437 END DO!! 438 END DO!! 439 CALL lbc_lnk( pgtui(:,:,jn), 'U', -1. ) ; CALL lbc_lnk( pgtvi(:,:,jn), 'V', -1. ) ! Lateral boundary cond. 440 ! 441 END DO 442 443 IF( PRESENT( prd ) ) THEN !== horizontal derivative of density anomalies (rd) ==! (optional part) 444 ! 445 pgrui(:,:) =0.0_wp ; pgrvi(:,:) =0.0_wp ; 446 pgzui(:,:) =0.0_wp ; pgzvi(:,:) =0.0_wp ; 447 pmrui(:,:) =0.0_wp ; pmrui(:,:) =0.0_wp ; 448 pge3rui(:,:)=0.0_wp ; pge3rvi(:,:)=0.0_wp ; 449 ! 450 DO jj = 1, jpjm1 ! depth of the partial step level 451 DO ji = 1, jpim1 452 iku = miku(ji,jj) 453 ikv = mikv(ji,jj) 454 ze3wu = (gdepw_0(ji+1,jj,iku+1) - gdept_0(ji+1,jj,iku)) - (gdepw_0(ji,jj,iku+1) - gdept_0(ji,jj,iku)) 455 ze3wv = (gdepw_0(ji,jj+1,ikv+1) - gdept_0(ji,jj+1,ikv)) - (gdepw_0(ji,jj,ikv+1) - gdept_0(ji,jj,ikv)) 456 ! 457 IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = fsdept(ji+1,jj,iku) + ze3wu ! i-direction: case 1 458 ELSE ; zhi(ji,jj) = fsdept(ji ,jj,iku) - ze3wu ! - - case 2 459 ENDIF 460 IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = fsdept(ji,jj+1,ikv) + ze3wv ! j-direction: case 1 461 ELSE ; zhj(ji,jj) = fsdept(ji,jj ,ikv) - ze3wv ! - - case 2 462 ENDIF 463 END DO 464 END DO 465 ! 466 CALL eos( zti, zhi, zri ) ! interpolated density from zti, ztj 467 CALL eos( ztj, zhj, zrj ) ! at the partial step depth output in zri, zrj 468 ! 469 DO jj = 1, jpjm1 ! Gradient of density at the last level 470 DO ji = 1, jpim1 471 iku = miku(ji,jj) ; ikup1 = miku(ji,jj) + 1 472 ikv = mikv(ji,jj) ; ikvp1 = mikv(ji,jj) + 1 473 ze3wu = (gdepw_0(ji+1,jj,iku+1) - gdept_0(ji+1,jj,iku)) - (gdepw_0(ji,jj,iku+1) - gdept_0(ji,jj,iku)) 474 ze3wv = (gdepw_0(ji,jj+1,ikv+1) - gdept_0(ji,jj+1,ikv)) - (gdepw_0(ji,jj,ikv+1) - gdept_0(ji,jj,ikv)) 475 IF( ze3wu >= 0._wp ) THEN 476 pgzui (ji,jj) = (fsde3w(ji+1,jj,iku) + ze3wu) - fsde3w(ji,jj,iku) 477 pgrui (ji,jj) = umask(ji,jj,iku) * ( zri(ji,jj) - prd(ji,jj,iku) ) ! i: 1 478 pmrui (ji,jj) = umask(ji,jj,iku) * ( zri(ji,jj) + prd(ji,jj,iku) ) ! i: 1 479 pge3rui(ji,jj) = umask(ji,jj,iku+1) & 480 & * ( (fse3w(ji+1,jj,iku+1) - ze3wu) * (zri(ji,jj ) + prd(ji+1,jj,iku+1) + 2._wp) & 481 & - fse3w(ji ,jj,iku+1) * (prd(ji,jj,iku) + prd(ji ,jj,iku+1) + 2._wp) ) ! i: 1 482 ELSE 483 pgzui (ji,jj) = fsde3w(ji+1,jj,iku) - (fsde3w(ji,jj,iku) - ze3wu) 484 pgrui (ji,jj) = umask(ji,jj,iku) * ( prd(ji+1,jj,iku) - zri(ji,jj) ) ! i: 2 485 pmrui (ji,jj) = umask(ji,jj,iku) * ( prd(ji+1,jj,iku) + zri(ji,jj) ) ! i: 2 486 pge3rui(ji,jj) = umask(ji,jj,iku+1) & 487 & * ( fse3w(ji+1,jj,iku+1) * (prd(ji+1,jj,iku) + prd(ji+1,jj,iku+1) + 2._wp) & 488 & -(fse3w(ji ,jj,iku+1) + ze3wu) * (zri(ji,jj ) + prd(ji ,jj,iku+1) + 2._wp) ) ! i: 2 489 ENDIF 490 IF( ze3wv >= 0._wp ) THEN 491 pgzvi (ji,jj) = (fsde3w(ji,jj+1,ikv) + ze3wv) - fsde3w(ji,jj,ikv) 492 pgrvi (ji,jj) = vmask(ji,jj,ikv) * ( zrj(ji,jj ) - prd(ji,jj,ikv) ) ! j: 1 493 pmrvi (ji,jj) = vmask(ji,jj,ikv) * ( zrj(ji,jj ) + prd(ji,jj,ikv) ) ! j: 1 494 pge3rvi(ji,jj) = vmask(ji,jj,ikv+1) & 495 & * ( (fse3w(ji,jj+1,ikv+1) - ze3wv) * ( zrj(ji,jj ) + prd(ji,jj+1,ikv+1) + 2._wp) & 496 & - fse3w(ji,jj ,ikv+1) * ( prd(ji,jj,ikv) + prd(ji,jj ,ikv+1) + 2._wp) ) ! j: 1 497 ! + 2 due to the formulation in density and not in anomalie in hpg sco 498 ELSE 499 pgzvi (ji,jj) = fsde3w(ji,jj+1,ikv) - (fsde3w(ji,jj,ikv) - ze3wv) 500 pgrvi (ji,jj) = vmask(ji,jj,ikv) * ( prd(ji,jj+1,ikv) - zrj(ji,jj) ) ! j: 2 501 pmrvi (ji,jj) = vmask(ji,jj,ikv) * ( prd(ji,jj+1,ikv) + zrj(ji,jj) ) ! j: 2 502 pge3rvi(ji,jj) = vmask(ji,jj,ikv+1) & 503 & * ( fse3w(ji,jj+1,ikv+1) * ( prd(ji,jj+1,ikv) + prd(ji,jj+1,ikv+1) + 2._wp) & 504 & -(fse3w(ji,jj ,ikv+1) + ze3wv) * ( zrj(ji,jj ) + prd(ji,jj ,ikv+1) + 2._wp) ) ! j: 2 505 ENDIF 506 END DO 507 END DO 508 CALL lbc_lnk( pgrui , 'U', -1. ) ; CALL lbc_lnk( pgrvi , 'V', -1. ) ! Lateral boundary conditions 509 CALL lbc_lnk( pmrui , 'U', 1. ) ; CALL lbc_lnk( pmrvi , 'V', 1. ) ! Lateral boundary conditions 510 CALL lbc_lnk( pgzui , 'U', -1. ) ; CALL lbc_lnk( pgzvi , 'V', -1. ) ! Lateral boundary conditions 511 CALL lbc_lnk( pge3rui , 'U', -1. ) ; CALL lbc_lnk( pge3rvi , 'V', -1. ) ! Lateral boundary conditions 437 iku = miku(ji,jj) 438 ikv = mikv(ji,jj) 439 ze3wu = gdept_n(ji,jj,iku) - gdept_n(ji+1,jj,iku) 440 ze3wv = gdept_n(ji,jj,ikv) - gdept_n(ji,jj+1,ikv) 441 442 IF( ze3wu >= 0._wp ) THEN ; pgrui(ji,jj) = ssumask(ji,jj) * ( zri(ji ,jj ) - prd(ji,jj,iku) ) ! i: 1 443 ELSE ; pgrui(ji,jj) = ssumask(ji,jj) * ( prd(ji+1,jj ,iku) - zri(ji,jj ) ) ! i: 2 444 ENDIF 445 IF( ze3wv >= 0._wp ) THEN ; pgrvi(ji,jj) = ssvmask(ji,jj) * ( zrj(ji ,jj ) - prd(ji,jj,ikv) ) ! j: 1 446 ELSE ; pgrvi(ji,jj) = ssvmask(ji,jj) * ( prd(ji ,jj+1,ikv) - zrj(ji,jj ) ) ! j: 2 447 ENDIF 448 449 END DO 450 END DO 451 CALL lbc_lnk( pgrui , 'U', -1. ); CALL lbc_lnk( pgrvi , 'V', -1. ) ! Lateral boundary conditions 512 452 ! 513 453 END IF
Note: See TracChangeset
for help on using the changeset viewer.