- Timestamp:
- 2018-08-03T10:18:16+02:00 (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2018/dev_r9838_ENHANCE04_RK3/src/OCE/DYN/dynzdf.F90
r10001 r10030 36 36 37 37 PUBLIC dyn_zdf ! routine called by step.F90 38 39 REAL(wp) :: r_vvl ! non-linear free surface indicator: =0 if ln_linssh=T, =1 otherwise40 38 41 39 !! * Substitutions … … 70 68 INTEGER, INTENT(in) :: kt ! ocean time-step index 71 69 ! 72 INTEGER :: ji, jj, jk 73 INTEGER :: iku, ikv 74 REAL(wp) :: zzwi, z e3ua, z2dt_2 ! local scalars75 REAL(wp) :: zzws , ze3va! - -70 INTEGER :: ji, jj, jk ! dummy loop indices 71 INTEGER :: iku, ikv ! local integers 72 REAL(wp) :: zzwi, zDt_2 ! local scalars 73 REAL(wp) :: zzws ! - - 76 74 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwi, zwd, zws ! 3D workspace 77 75 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ztrdu, ztrdv ! - - 76 REAL(wp) :: ze3uw_A , ze3uw_Ap1 ! local real 77 REAL(wp) :: ze3vw_A , ze3vw_Ap1 ! local real 78 REAL(wp), DIMENSION(jpi,jpj) :: zsshu_hA, zsshv_hA ! 2D workspace 78 79 !!--------------------------------------------------------------------- 79 80 ! … … 84 85 IF(lwp) WRITE(numout,*) 'dyn_zdf_imp : vertical momentum diffusion implicit operator' 85 86 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ ' 86 ! 87 If( ln_linssh ) THEN ; r_vvl = 0._wp ! non-linear free surface indicator 88 ELSE ; r_vvl = 1._wp 89 ENDIF 90 ENDIF 91 ! 92 z2dt_2 = rDt * 0.5_wp !* =rn_Dt except in 1st Euler time step where it is equal to rn_Dt/2 87 ENDIF 88 ! 89 zDt_2 = rDt * 0.5_wp !* =rn_Dt except in 1st Euler time step where it is equal to rn_Dt/2 93 90 ! 94 91 ! … … 131 128 iku = mbku(ji,jj) ! ocean bottom level at u- and v-points 132 129 ikv = mbkv(ji,jj) ! (deepest ocean u- and v-points) 133 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,iku) + r_vvl * e3u_a(ji,jj,iku) 134 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,ikv) + r_vvl * e3v_a(ji,jj,ikv) 135 ua(ji,jj,iku) = ua(ji,jj,iku) + z2dt_2 * ( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * ua_b(ji,jj) / ze3ua 136 va(ji,jj,ikv) = va(ji,jj,ikv) + z2dt_2 * ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * va_b(ji,jj) / ze3va 130 ua(ji,jj,iku) = ua(ji,jj,iku) + zDt_2 * ( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * ua_b(ji,jj) / e3u_a(ji,jj,iku) 131 va(ji,jj,ikv) = va(ji,jj,ikv) + zDt_2 * ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * va_b(ji,jj) / e3v_a(ji,jj,ikv) 137 132 END DO 138 133 END DO … … 142 137 iku = miku(ji,jj) ! top ocean level at u- and v-points 143 138 ikv = mikv(ji,jj) ! (first wet ocean u- and v-points) 144 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,iku) + r_vvl * e3u_a(ji,jj,iku) 145 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,ikv) + r_vvl * e3v_a(ji,jj,ikv) 146 ua(ji,jj,iku) = ua(ji,jj,iku) + z2dt_2 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * ua_b(ji,jj) / ze3ua 147 va(ji,jj,ikv) = va(ji,jj,ikv) + z2dt_2 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * va_b(ji,jj) / ze3va 139 ua(ji,jj,iku) = ua(ji,jj,iku) + zDt_2 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * ua_b(ji,jj) / e3u_a(ji,jj,iku) 140 va(ji,jj,ikv) = va(ji,jj,ikv) + zDt_2 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * va_b(ji,jj) / e3v_a(ji,jj,ikv) 148 141 END DO 149 142 END DO … … 153 146 ! !== Vertical diffusion on u ==! 154 147 ! 148 ! !* multiplicative factors on e3uw(Naa) and e3vw(Naa) 149 ! 150 IF( ln_linssh ) THEN !-- linear ssh case 151 DO jj = 1, jpjm1 152 DO ji = 1, jpim1 153 zsshu_hA(ji,jj) = 0._wp ! no time variation in e3 154 zsshv_hA(ji,jj) = 0._wp 155 END DO 156 END DO 157 ELSE !-- Non linear ssh case 158 DO jj = 1, jpjm1 159 DO ji = 1, jpim1 160 zsshu_hA(ji,jj) = 0.5_wp * ( ssh(ji,jj,Naa) + ssh(ji+1,jj ,Naa) ) * r1_hu_0(ji,jj) * ssumask(ji,jj) 161 zsshv_hA(ji,jj) = 0.5_wp * ( ssh(ji,jj,Naa) + ssh(ji ,jj+1,Naa) ) * r1_hv_0(ji,jj) * ssvmask(ji,jj) 162 END DO 163 END DO 164 ENDIF 165 166 167 168 169 155 170 SELECT CASE( nldf_dyn ) !* Matrix construction 156 171 ! … … 159 174 DO jj = 2, jpjm1 160 175 DO ji = fs_2, fs_jpim1 ! vector opt. 161 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk) ! after scale factor at T-point 162 zzwi = - rDt * ( 0.5 * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) ) + akzu(ji,jj,jk ) ) & 163 & / ( ze3ua * e3uw_n(ji,jj,jk ) ) * wumask(ji,jj,jk ) 164 zzws = - rDt * ( 0.5 * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) + akzu(ji,jj,jk+1) ) & 165 & / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 176 !!gm Note that below, since ze3uw_A is used in a expression masked by wumask, 177 !! one can remove wumask from its expression (same for ze3uw_Ap1 178 ze3uw_A = e3uw_0(ji,jj,jk ) * ( 1._wp + zsshu_hA(ji,jj) * wumask(ji,jj,jk ) ) 179 ze3uw_Ap1 = e3uw_0(ji,jj,jk+1) * ( 1._wp + zsshu_hA(ji,jj) * wumask(ji,jj,jk+1) ) 180 ! 181 zzwi = - rDt * ( 0.5 * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) ) + akzu(ji,jj,jk ) ) & 182 & / ( e3u_a(ji ,jj,jk ) * ze3uw_A ) * wumask(ji,jj,jk ) 183 zzws = - rDt * ( 0.5 * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) + akzu(ji,jj,jk+1) ) & 184 & / ( e3u_a(ji ,jj,jk ) * ze3uw_Ap1 ) * wumask(ji,jj,jk+1) 166 185 zwi(ji,jj,jk) = zzwi 167 186 zws(ji,jj,jk) = zzws … … 174 193 DO jj = 2, jpjm1 175 194 DO ji = fs_2, fs_jpim1 ! vector opt. 176 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk) ! after scale factor at T-point 177 zzwi = - z2dt_2 * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) ) / ( ze3ua * e3uw_n(ji,jj,jk ) ) * wumask(ji,jj,jk ) 178 zzws = - z2dt_2 * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 195 !!gm Note that below, since ze3uw_A is used in a expression masked by wumask, 196 !! one can remove wumask from its expression (same for ze3uw_Ap1 197 ze3uw_A = e3uw_0(ji,jj,jk ) * ( 1._wp + zsshu_hA(ji,jj) * wumask(ji,jj,jk ) ) 198 ze3uw_Ap1 = e3uw_0(ji,jj,jk+1) * ( 1._wp + zsshu_hA(ji,jj) * wumask(ji,jj,jk+1) ) 199 ! 200 zzwi = - zDt_2 * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) ) / ( e3u_a(ji,jj,jk) * ze3uw_A ) * wumask(ji,jj,jk ) 201 zzws = - zDt_2 * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( e3u_a(ji,jj,jk) * ze3uw_Ap1 ) * wumask(ji,jj,jk+1) 179 202 zwi(ji,jj,jk) = zzwi 180 203 zws(ji,jj,jk) = zzws … … 202 225 DO ji = 2, jpim1 203 226 iku = mbku(ji,jj) ! ocean bottom level at u- and v-points 204 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,iku) + r_vvl * e3u_a(ji,jj,iku) ! after scale factor at T-point 205 zwd(ji,jj,iku) = zwd(ji,jj,iku) - z2dt_2 * ( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / ze3ua 227 zwd(ji,jj,iku) = zwd(ji,jj,iku) - zDt_2 * ( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / e3u_a(ji,jj,iku) 206 228 END DO 207 229 END DO … … 211 233 !!gm top Cd is masked (=0 outside cavities) no need of test on mik>=2 ==>> it has been suppressed 212 234 iku = miku(ji,jj) ! ocean top level at u- and v-points 213 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,iku) + r_vvl * e3u_a(ji,jj,iku) ! after scale factor at T-point 214 zwd(ji,jj,iku) = zwd(ji,jj,iku) - z2dt_2 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3ua 235 zwd(ji,jj,iku) = zwd(ji,jj,iku) - zDt_2 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / e3u_a(ji,jj,iku) 215 236 END DO 216 237 END DO … … 243 264 DO jj = 2, jpjm1 !== second recurrence: SOLk = RHSk - Lk / Dk-1 Lk-1 ==! 244 265 DO ji = fs_2, fs_jpim1 ! vector opt. 245 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,1) + r_vvl * e3u_a(ji,jj,1) 246 ua(ji,jj,1) = ua(ji,jj,1) + z2dt_2 * ( utau_b(ji,jj) + utau(ji,jj) ) / ( ze3ua * rho0 ) * umask(ji,jj,1) 266 ua(ji,jj,1) = ua(ji,jj,1) + zDt_2 * ( utau_b(ji,jj) + utau(ji,jj) ) / ( e3u_a(ji,jj,1) * rho0 ) * umask(ji,jj,1) 247 267 END DO 248 268 END DO … … 276 296 DO jj = 2, jpjm1 277 297 DO ji = fs_2, fs_jpim1 ! vector opt. 278 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk) ! after scale factor at T-point 279 zzwi = - rDt * ( 0.5 * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) ) + akzv(ji,jj,jk ) ) & 280 & / ( ze3va * e3vw_n(ji,jj,jk ) ) * wvmask(ji,jj,jk ) 281 zzws = - rDt * ( 0.5 * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) + akzv(ji,jj,jk+1) ) & 282 & / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 298 ze3vw_A = e3vw_0(ji,jj,jk ) * ( 1._wp + zsshv_hA(ji,jj) * wvmask(ji,jj,jk ) ) 299 ze3vw_Ap1 = e3vw_0(ji,jj,jk+1) * ( 1._wp + zsshv_hA(ji,jj) * wvmask(ji,jj,jk+1) ) 300 ! 301 zzwi = - rDt * ( 0.5 * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) ) + akzv(ji,jj,jk ) ) & 302 & / ( e3v_a(ji,jj ,jk ) * ze3vw_A ) * wvmask(ji,jj,jk ) 303 zzws = - rDt * ( 0.5 * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) + akzv(ji,jj,jk+1) ) & 304 & / ( e3v_a(ji,jj ,jk ) * ze3vw_Ap1 ) * wvmask(ji,jj,jk+1) 283 305 zwi(ji,jj,jk) = zzwi * wvmask(ji,jj,jk ) 284 306 zws(ji,jj,jk) = zzws * wvmask(ji,jj,jk+1) … … 291 313 DO jj = 2, jpjm1 292 314 DO ji = fs_2, fs_jpim1 ! vector opt. 293 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk) ! after scale factor at T-point 294 zzwi = - z2dt_2 * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) ) / ( ze3va * e3vw_n(ji,jj,jk ) ) * wvmask(ji,jj,jk ) 295 zzws = - z2dt_2 * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 315 ze3vw_A = e3vw_0(ji,jj,jk ) * ( 1._wp + zsshv_hA(ji,jj) * wvmask(ji,jj,jk ) ) 316 ze3vw_Ap1 = e3vw_0(ji,jj,jk+1) * ( 1._wp + zsshv_hA(ji,jj) * wvmask(ji,jj,jk+1) ) 317 ! 318 zzwi = - zDt_2 * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) ) / ( e3v_a(ji,jj,jk) * ze3vw_A ) * wvmask(ji,jj,jk ) 319 zzws = - zDt_2 * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( e3v_a(ji,jj,jk) * ze3vw_Ap1 ) * wvmask(ji,jj,jk+1) 296 320 zwi(ji,jj,jk) = zzwi * wvmask(ji,jj,jk ) 297 321 zws(ji,jj,jk) = zzws * wvmask(ji,jj,jk+1) … … 318 342 DO ji = 2, jpim1 319 343 ikv = mbkv(ji,jj) ! (deepest ocean u- and v-points) 320 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,ikv) + r_vvl * e3v_a(ji,jj,ikv) ! after scale factor at T-point 321 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - z2dt_2 * ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / ze3va 344 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - zDt_2 * ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / e3v_a(ji,jj,ikv) 322 345 END DO 323 346 END DO … … 326 349 DO ji = 2, jpim1 327 350 ikv = mikv(ji,jj) ! (first wet ocean u- and v-points) 328 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,ikv) + r_vvl * e3v_a(ji,jj,ikv) ! after scale factor at T-point 329 zwd(ji,jj,iku) = zwd(ji,jj,iku) - z2dt_2 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3va 351 zwd(ji,jj,iku) = zwd(ji,jj,iku) - zDt_2 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / e3v_a(ji,jj,ikv) 330 352 END DO 331 353 END DO … … 358 380 DO jj = 2, jpjm1 !== second recurrence: SOLk = RHSk - Lk / Dk-1 Lk-1 ==! 359 381 DO ji = fs_2, fs_jpim1 ! vector opt. 360 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,1) + r_vvl * e3v_a(ji,jj,1) 361 va(ji,jj,1) = va(ji,jj,1) + z2dt_2 * ( vtau_b(ji,jj) + vtau(ji,jj) ) / ( ze3va * rho0 ) * vmask(ji,jj,1) 382 va(ji,jj,1) = va(ji,jj,1) + zDt_2 * ( vtau_b(ji,jj) + vtau(ji,jj) ) / ( e3v_a(ji,jj,1) * rho0 ) * vmask(ji,jj,1) 362 383 END DO 363 384 END DO … … 402 423 END SUBROUTINE dyn_zdf 403 424 404 !!gm currently not used : just for memory to be able to add dissipation trend through vertical mixing405 406 SUBROUTINE zdf_trd( ptrdu, ptrdv ,kt )407 !!----------------------------------------------------------------------408 !! *** ROUTINE zdf_trd ***409 !!410 !! ** Purpose : compute the trend due to the vert. momentum diffusion411 !! together with the Leap-Frog time stepping using an412 !! implicit scheme.413 !!414 !! ** Method : - Leap-Frog time stepping on all trends but the vertical mixing415 !! ua = ub + 2*dt * ua vector form or linear free surf.416 !! ua = ( e3u_b*ub + 2*dt * e3u_n*ua ) / e3u_a otherwise417 !! - update the after velocity with the implicit vertical mixing.418 !! This requires to solver the following system:419 !! ua = ua + 1/e3u_a dk+1[ mi(avm) / e3uw_a dk[ua] ]420 !! with the following surface/top/bottom boundary condition:421 !! surface: wind stress input (averaged over kt-1/2 & kt+1/2)422 !! top & bottom : top stress (iceshelf-ocean) & bottom stress (cf zdfdrg.F90)423 !!424 !! ** Action : (ua,va) after velocity425 !!---------------------------------------------------------------------426 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: ptrdu, ptrdv ! 3D work arrays use for zdf trends diag427 INTEGER , INTENT(in ) :: kt ! ocean time-step index428 !429 INTEGER :: ji, jj, jk ! dummy loop indices430 REAL(wp) :: zzz ! local scalar431 REAL(wp) :: zavmu, zavmum1 ! - -432 REAL(wp) :: zavmv, zavmvm1 ! - -433 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: z2d ! - -434 !!---------------------------------------------------------------------435 !436 CALL lbc_lnk_multi( ua, 'U', -1., va, 'V', -1. ) ! apply lateral boundary condition on (ua,va)437 !438 !439 ! !== momentum trend due to vertical diffusion ==!440 !441 IF( ln_dynadv_vec .OR. ln_linssh ) THEN ! applied on velocity442 ptrdu(:,:,:) = ( ua(:,:,:) - ub(:,:,:) ) * r1_Dt - ptrdu(:,:,:)443 ptrdv(:,:,:) = ( va(:,:,:) - vb(:,:,:) ) * r1_Dt - ptrdv(:,:,:)444 ELSE ! applied on thickness weighted velocity445 ptrdu(:,:,:) = ( e3u_a(:,:,:)*ua(:,:,:) - e3u_b(:,:,:)*ub(:,:,:) ) / e3u_n(:,:,:) * r1_Dt - ptrdu(:,:,:)446 ptrdv(:,:,:) = ( e3v_a(:,:,:)*va(:,:,:) - e3v_b(:,:,:)*vb(:,:,:) ) / e3v_n(:,:,:) * r1_Dt - ptrdv(:,:,:)447 ENDIF448 CALL trd_dyn( ptrdu, ptrdv, jpdyn_zdf, kt )449 !450 !451 ! !== KE dissipation trend due to vertical diffusion ==!452 !453 IF( iom_use( 'dispkevfo' ) ) THEN ! ocean kinetic energy dissipation per unit area454 ! ! due to v friction (v=vertical)455 ! ! see NEMO_book appendix C, §C.8 (N.B. here averaged at t-points)456 ! ! Note that formally, in a Leap-Frog environment, the shear**2 should be the product of457 ! ! now by before shears, i.e. the source term of TKE (local positivity is not ensured).458 ! ! Note also that now e3 scale factors are used as after one are not computed !459 !460 ALLOCATE( z2d(jpi,jpj) )461 z2d(:,:) = 0._wp462 DO jk = 1, jpkm1463 DO jj = 2, jpjm1464 DO ji = 2, jpim1465 zavmu = 0.5 * ( avm(ji+1,jj,jk) + avm(ji ,jj,jk) )466 zavmum1 = 0.5 * ( avm(ji ,jj,jk) + avm(ji-1,jj,jk) )467 zavmv = 0.5 * ( avm(ji,jj+1,jk) + avm(ji,jj ,jk) )468 zavmvm1 = 0.5 * ( avm(ji,jj ,jk) + avm(ji,jj-1,jk) )469 470 z2d(ji,jj) = z2d(ji,jj) + ( &471 & zavmu * ( ua(ji ,jj,jk-1) - ua(ji ,jj,jk) )**2 / e3uw_n(ji ,jj,jk) * wumask(ji ,jj,jk) &472 & + zavmum1 * ( ua(ji-1,jj,jk-1) - ua(ji-1,jj,jk) )**2 / e3uw_n(ji-1,jj,jk) * wumask(ji-1,jj,jk) &473 & + zavmv * ( va(ji,jj ,jk-1) - va(ji,jj ,jk) )**2 / e3vw_n(ji,jj ,jk) * wvmask(ji,jj ,jk) &474 & + zavmvm1 * ( va(ji,jj-1,jk-1) - va(ji,jj-1,jk) )**2 / e3vw_n(ji,jj-1,jk) * wvmask(ji,jj-1,jk) &475 & )476 !!gm --- This trends is in fact properly computed in zdf_sh2 but with a backward shift of one time-step ===>>> use it ?477 !! No since in zdfshé only kz tke (or gls) is used478 !!479 !!gm --- formally, as done below, in a Leap-Frog environment, the shear**2 should be the product of480 !!gm now by before shears, i.e. the source term of TKE (local positivity is not ensured).481 !! CAUTION: requires to compute e3uw_a and e3vw_a !!!482 ! z2d(ji,jj) = z2d(ji,jj) + ( &483 ! & avmu(ji ,jj,jk) * ( un(ji ,jj,jk-1) - un(ji ,jj,jk) ) / e3uw_n(ji ,jj,jk) &484 ! & * ( ua(ji ,jj,jk-1) - ua(ji ,jj,jk) ) / e3uw_a(ji ,jj,jk) * wumask(ji ,jj,jk) &485 ! & + avmu(ji-1,jj,jk) * ( un(ji-1,jj,jk-1) - un(ji-1,jj,jk) ) / e3uw_n(ji-1,jj,jk) &486 ! & ( ua(ji-1,jj,jk-1) - ua(ji-1,jj,jk) ) / e3uw_a(ji-1,jj,jk) * wumask(ji-1,jj,jk) &487 ! & + avmv(ji,jj ,jk) * ( vn(ji,jj ,jk-1) - vn(ji,jj ,jk) ) / e3vw_n(ji,jj ,jk) &488 ! & ( va(ji,jj ,jk-1) - va(ji,jj ,jk) ) / e3vw_a(ji,jj ,jk) * wvmask(ji,jj ,jk) &489 ! & + avmv(ji,jj-1,jk) * ( vn(ji,jj-1,jk-1) - vn(ji,jj-1,jk) ) / e3vw_n(ji,jj-1,jk) &490 ! & ( va(ji,jj-1,jk-1) - va(ji,jj-1,jk) ) / e3vw_a(ji,jj-1,jk) * wvmask(ji,jj-1,jk) &491 ! & )492 !!gm end493 END DO494 END DO495 END DO496 zzz= - 0.5_wp* rho0 ! caution sign minus here497 z2d(:,:) = zzz * z2d(:,:)498 CALL lbc_lnk( z2d,'T', 1. )499 CALL iom_put( 'dispkevfo', z2d )500 DEALLOCATE( z2d )501 ENDIF502 !503 END SUBROUTINE zdf_trd504 505 !!gm end not used506 507 425 !!============================================================================== 508 426 END MODULE dynzdf
Note: See TracChangeset
for help on using the changeset viewer.
