Changeset 11082 for NEMO/branches/UKMO/NEMO_4.0_GO8_package/src/OCE/DYN
- Timestamp:
- 2019-06-06T16:21:52+02:00 (5 years ago)
- Location:
- NEMO/branches/UKMO/NEMO_4.0_GO8_package/src/OCE/DYN
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/UKMO/NEMO_4.0_GO8_package/src/OCE/DYN/dynkeg.F90
r10888 r11082 74 74 INTEGER, INTENT( in ) :: kscheme ! =0/1 type of KEG scheme 75 75 ! 76 INTEGER :: ji, jj, jk, jb ! dummy loop indices 77 INTEGER :: ii, ifu, ib_bdy ! local integers 78 INTEGER :: ij, ifv, igrd ! - - 79 REAL(wp) :: zu, zv ! local scalars 76 INTEGER :: ji, jj, jk, jb ! dummy loop indices 77 INTEGER :: ifu, ifv, igrd, ib_bdy ! local integers 78 REAL(wp) :: zu, zv ! local scalars 80 79 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zhke 81 80 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdu, ztrdv 81 REAL(wp) :: zweightu, zweightv 82 82 !!---------------------------------------------------------------------- 83 83 ! … … 97 97 98 98 zhke(:,:,jpk) = 0._wp 99 100 IF (ln_bdy) THEN101 ! Maria Luneva & Fred Wobus: July-2016102 ! compensate for lack of turbulent kinetic energy on liquid bdy points103 DO ib_bdy = 1, nb_bdy104 IF( cn_dyn3d(ib_bdy) /= 'none' ) THEN105 igrd = 2 ! Copying normal velocity into points outside bdy106 DO jb = 1, idx_bdy(ib_bdy)%nblenrim(igrd)107 DO jk = 1, jpkm1108 ii = idx_bdy(ib_bdy)%nbi(jb,igrd)109 ij = idx_bdy(ib_bdy)%nbj(jb,igrd)110 ifu = NINT( idx_bdy(ib_bdy)%flagu(jb,igrd) )111 un(ii-ifu,ij,jk) = un(ii,ij,jk) * umask(ii,ij,jk)112 END DO113 END DO114 !115 igrd = 3 ! Copying normal velocity into points outside bdy116 DO jb = 1, idx_bdy(ib_bdy)%nblenrim(igrd)117 DO jk = 1, jpkm1118 ii = idx_bdy(ib_bdy)%nbi(jb,igrd)119 ij = idx_bdy(ib_bdy)%nbj(jb,igrd)120 ifv = NINT( idx_bdy(ib_bdy)%flagv(jb,igrd) )121 vn(ii,ij-ifv,jk) = vn(ii,ij,jk) * vmask(ii,ij,jk)122 END DO123 END DO124 ENDIF125 ENDDO126 ENDIF127 99 128 100 SELECT CASE ( kscheme ) !== Horizontal kinetic energy at T-point ==! … … 140 112 END DO 141 113 END DO 114 ! 115 IF (ln_bdy) THEN 116 ! Maria Luneva & Fred Wobus: July-2016 117 ! compensate for lack of turbulent kinetic energy on liquid bdy points 118 DO ib_bdy = 1, nb_bdy 119 IF( cn_dyn3d(ib_bdy) /= 'none' ) THEN 120 igrd = 1 ! compensating null velocity on the bdy 121 DO jb = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 122 ji = idx_bdy(ib_bdy)%nbi(jb,igrd) ! maximum extent : from 2 to jpi-1 123 jj = idx_bdy(ib_bdy)%nbj(jb,igrd) ! maximum extent : from 2 to jpj-1 124 DO jk = 1, jpkm1 125 zhke(ji,jj,jk) = 0._wp 126 zweightu = umask(ji-1,jj ,jk) + umask(ji,jj,jk) 127 zweightv = vmask(ji ,jj-1,jk) + vmask(ji,jj,jk) 128 zu = un(ji-1,jj ,jk) * un(ji-1,jj ,jk) + un(ji ,jj ,jk) * un(ji ,jj ,jk) 129 zv = vn(ji ,jj-1,jk) * vn(ji ,jj-1,jk) + vn(ji ,jj ,jk) * vn(ji ,jj ,jk) 130 IF( zweightu > 0._wp ) zhke(ji,jj,jk) = zhke(ji,jj,jk) + zu / (2._wp * zweightu) 131 IF( zweightv > 0._wp ) zhke(ji,jj,jk) = zhke(ji,jj,jk) + zv / (2._wp * zweightv) 132 END DO 133 END DO 134 END IF 135 CALL lbc_bdy_lnk( 'dynkeg', zhke, 'T', 1., ib_bdy ) ! send 2 and recv jpi, jpj used in the computation of the speed tendencies 136 END DO 137 END IF 142 138 ! 143 139 CASE ( nkeg_HW ) !-- Hollingsworth scheme --! … … 158 154 END DO 159 155 END DO 156 IF (ln_bdy) THEN 157 ! Maria Luneva & Fred Wobus: July-2016 158 ! compensate for lack of turbulent kinetic energy on liquid bdy points 159 DO ib_bdy = 1, nb_bdy 160 IF( cn_dyn3d(ib_bdy) /= 'none' ) THEN 161 igrd = 1 ! compensation null velocity on land at the bdy 162 DO jb = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 163 ji = idx_bdy(ib_bdy)%nbi(jb,igrd) ! maximum extent : from 2 to jpi-1 164 jj = idx_bdy(ib_bdy)%nbj(jb,igrd) ! maximum extent : from 2 to jpj-1 165 DO jk = 1, jpkm1 166 zhke(ji,jj,jk) = 0._wp 167 zweightu = 8._wp * ( umask(ji-1,jj ,jk) + umask(ji ,jj ,jk) ) & 168 & + 2._wp * ( umask(ji-1,jj-1,jk) + umask(ji-1,jj+1,jk) + umask(ji ,jj-1,jk) + umask(ji ,jj+1,jk) ) 169 zweightv = 8._wp * ( vmask(ji ,jj-1,jk) + vmask(ji ,jj-1,jk) ) & 170 & + 2._wp * ( vmask(ji-1,jj-1,jk) + vmask(ji+1,jj-1,jk) + vmask(ji-1,jj ,jk) + vmask(ji+1,jj ,jk) ) 171 zu = 8._wp * ( un(ji-1,jj ,jk) * un(ji-1,jj ,jk) & 172 & + un(ji ,jj ,jk) * un(ji ,jj ,jk) ) & 173 & + ( un(ji-1,jj-1,jk) + un(ji-1,jj+1,jk) ) * ( un(ji-1,jj-1,jk) + un(ji-1,jj+1,jk) ) & 174 & + ( un(ji ,jj-1,jk) + un(ji ,jj+1,jk) ) * ( un(ji ,jj-1,jk) + un(ji ,jj+1,jk) ) 175 zv = 8._wp * ( vn(ji ,jj-1,jk) * vn(ji ,jj-1,jk) & 176 & + vn(ji ,jj ,jk) * vn(ji ,jj ,jk) ) & 177 & + ( vn(ji-1,jj-1,jk) + vn(ji+1,jj-1,jk) ) * ( vn(ji-1,jj-1,jk) + vn(ji+1,jj-1,jk) ) & 178 & + ( vn(ji-1,jj ,jk) + vn(ji+1,jj ,jk) ) * ( vn(ji-1,jj ,jk) + vn(ji+1,jj ,jk) ) 179 IF( zweightu > 0._wp ) zhke(ji,jj,jk) = zhke(ji,jj,jk) + zu / ( 2._wp * zweightu ) 180 IF( zweightv > 0._wp ) zhke(ji,jj,jk) = zhke(ji,jj,jk) + zv / ( 2._wp * zweightv ) 181 END DO 182 END DO 183 END IF 184 END DO 185 END IF 160 186 CALL lbc_lnk( 'dynkeg', zhke, 'T', 1. ) 161 187 ! 162 END SELECT 163 164 IF (ln_bdy) THEN 165 ! restore velocity masks at points outside boundary 166 un(:,:,:) = un(:,:,:) * umask(:,:,:) 167 vn(:,:,:) = vn(:,:,:) * vmask(:,:,:) 168 ENDIF 169 188 END SELECT 170 189 ! 171 190 DO jk = 1, jpkm1 !== grad( KE ) added to the general momentum trends ==! -
NEMO/branches/UKMO/NEMO_4.0_GO8_package/src/OCE/DYN/sshwzv.F90
r10888 r11082 297 297 IF(lwp) WRITE(numout,*) 'wAimp : Courant number-based partitioning of now vertical velocity ' 298 298 IF(lwp) WRITE(numout,*) '~~~~~ ' 299 ! 300 Cu_adv(:,:,jpk) = 0._wp ! bottom value : Cu_adv=0 (set once for all) 301 ENDIF 299 ENDIF 300 ! 302 301 ! 303 302 DO jk = 1, jpkm1 ! calculate Courant numbers … … 305 304 DO ji = 2, fs_jpim1 ! vector opt. 306 305 z1_e3w = 1._wp / e3w_n(ji,jj,jk) 307 Cu_adv(ji,jj,jk) = r2dt * ( ( MAX( wn(ji,jj,jk) , 0._wp ) - MIN( wn(ji,jj,jk+1) , 0._wp ) ) &308 & + ( MAX( e2u(ji ,jj)*e3uw_n(ji ,jj,jk)*un(ji ,jj,jk), 0._wp ) - &309 & MIN( e2u(ji-1,jj)*e3uw_n(ji-1,jj,jk)*un(ji-1,jj,jk), 0._wp ) ) &310 & * r1_e1e2t(ji,jj) &311 & + ( MAX( e1v(ji,jj )*e3vw_n(ji,jj ,jk)*vn(ji,jj ,jk), 0._wp ) - &312 & MIN( e1v(ji,jj-1)*e3vw_n(ji,jj-1,jk)*vn(ji,jj-1,jk), 0._wp ) ) &313 & * r1_e1e2t(ji,jj) &314 & ) * z1_e3w306 Cu_adv(ji,jj,jk) = 2._wp * rdt * ( ( MAX( wn(ji,jj,jk) , 0._wp ) - MIN( wn(ji,jj,jk+1) , 0._wp ) ) & ! 2*rdt and not r2dt (for restartability) 307 & + ( MAX( e2u(ji ,jj)*e3uw_n(ji ,jj,jk)*un(ji ,jj,jk), 0._wp ) - & 308 & MIN( e2u(ji-1,jj)*e3uw_n(ji-1,jj,jk)*un(ji-1,jj,jk), 0._wp ) ) & 309 & * r1_e1e2t(ji,jj) & 310 & + ( MAX( e1v(ji,jj )*e3vw_n(ji,jj ,jk)*vn(ji,jj ,jk), 0._wp ) - & 311 & MIN( e1v(ji,jj-1)*e3vw_n(ji,jj-1,jk)*vn(ji,jj-1,jk), 0._wp ) ) & 312 & * r1_e1e2t(ji,jj) & 313 & ) * z1_e3w 315 314 END DO 316 315 END DO 317 316 END DO 317 CALL lbc_lnk( 'sshwzv', Cu_adv, 'T', 1. ) 318 318 ! 319 319 CALL iom_put("Courant",Cu_adv) 320 320 ! 321 wi(:,:,:) = 0._wp ! Includes top and bottom values set to zero322 321 IF( MAXVAL( Cu_adv(:,:,:) ) > Cu_min ) THEN ! Quick check if any breaches anywhere 323 322 DO jk = 1, jpkm1 ! or scan Courant criterion and partition 324 DO jj = 2, jpjm1! w where necessary325 DO ji = 2, fs_jpim1 ! vector opt.323 DO jj = 1, jpj ! w where necessary 324 DO ji = 1, jpi 326 325 ! 327 326 zCu = MAX( Cu_adv(ji,jj,jk) , Cu_adv(ji,jj,jk+1) ) 328 327 ! 329 IF( zCu < Cu_min ) THEN!<-- Fully explicit328 IF( zCu <= Cu_min ) THEN !<-- Fully explicit 330 329 zcff = 0._wp 331 330 ELSEIF( zCu < Cu_cut ) THEN !<-- Mixed explicit … … 346 345 ELSE 347 346 ! Fully explicit everywhere 348 Cu_adv = 0.0_wp ! Reuse array to output coefficient 347 Cu_adv(:,:,:) = 0._wp ! Reuse array to output coefficient 348 wi (:,:,:) = 0._wp 349 349 ENDIF 350 350 CALL iom_put("wimp",wi)
Note: See TracChangeset
for help on using the changeset viewer.