Changeset 4488
- Timestamp:
- 2014-02-06T11:43:09+01:00 (10 years ago)
- Location:
- branches/2013/dev_MERGE_2013/NEMOGCM
- Files:
-
- 9 edited
-
CONFIG/AMM12/EXP00/namelist_cfg (modified) (1 diff)
-
NEMO/OPA_SRC/DOM/dom_oce.F90 (modified) (3 diffs)
-
NEMO/OPA_SRC/DOM/domvvl.F90 (modified) (4 diffs)
-
NEMO/OPA_SRC/DOM/domzgr_substitute.h90 (modified) (2 diffs)
-
NEMO/OPA_SRC/DYN/dynldf.F90 (modified) (2 diffs)
-
NEMO/OPA_SRC/DYN/dynldf_bilapg.F90 (modified) (3 diffs)
-
NEMO/OPA_SRC/DYN/dynldf_iso.F90 (modified) (1 diff)
-
NEMO/OPA_SRC/LDF/ldfslp.F90 (modified) (3 diffs)
-
NEMO/OPA_SRC/TRA/traldf.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
branches/2013/dev_MERGE_2013/NEMOGCM/CONFIG/AMM12/EXP00/namelist_cfg
r4383 r4488 323 323 ! ! Type of the operator : 324 324 ln_dynldf_bilap = .true. ! bilaplacian operator 325 ln_dynldf_lap = .false. ! bilaplacian operator 326 ! ! Direction of action : 327 ln_dynldf_level = .true. ! iso-level 328 ln_dynldf_hor = .false. ! horizontal (geopotential) (require "key_ldfslp" in s-coord.) 325 329 ! Coefficient 326 330 rn_ahm_0_lap = 60.0 ! horizontal laplacian eddy viscosity [m2/s] 327 rn_ahmb_0 = 0.0 ! background eddy viscosity for ldf_iso [m2/s]328 331 rn_ahm_0_blp = -1.0e+10 ! horizontal bilaplacian eddy viscosity [m4/s] 329 332 / -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90
r4370 r4488 191 191 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdep3w_n !: now depth of T-points (sum of e3w) (m) 192 192 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdept_n, gdepw_n !: now depth at T-W points (m) 193 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdept_b, gdepw_b !: before depth at T-W points (m) 193 194 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: e3t_n !: now vertical scale factors at t point (m) 194 195 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: e3u_n , e3v_n !: - - - - u --v points (m) … … 196 197 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: e3uw_n , e3vw_n !: - - - - uw--vw points (m) 197 198 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: e3t_b !: before - - - - t points (m) 199 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: e3w_b !: before - - - - t points (m) 198 200 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: e3u_b , e3v_b !: - - - - - u --v points (m) 199 201 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: e3uw_b , e3vw_b !: - - - - - uw--vw points (m) … … 350 352 & e3t_b (jpi,jpj,jpk) , e3u_b (jpi,jpj,jpk) , e3v_b (jpi,jpj,jpk) , & 351 353 & e3uw_b (jpi,jpj,jpk) , e3vw_b(jpi,jpj,jpk) , & 352 & e3t_a (jpi,jpj,jpk) , e3u_a (jpi,jpj,jpk) , e3v_a (jpi,jpj,jpk), & 354 & gdept_b (jpi,jpj,jpk) ,gdepw_b(jpi,jpj,jpk) , e3w_b (jpi,jpj,jpk) , & 355 & e3t_a (jpi,jpj,jpk) , e3u_a (jpi,jpj,jpk) , e3v_a (jpi,jpj,jpk) , & 353 356 & ehu_a (jpi,jpj) , ehv_a (jpi,jpj), & 354 357 & ehur_a (jpi,jpj) , ehvr_a (jpi,jpj), & -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
r4486 r4488 159 159 CALL dom_vvl_interpol( fse3u_n(:,:,:), fse3uw_n(:,:,:), 'UW' ) 160 160 CALL dom_vvl_interpol( fse3v_n(:,:,:), fse3vw_n(:,:,:), 'VW' ) 161 CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3w_b (:,:,:), 'W' ) 161 162 CALL dom_vvl_interpol( fse3u_b(:,:,:), fse3uw_b(:,:,:), 'UW' ) 162 163 CALL dom_vvl_interpol( fse3v_b(:,:,:), fse3vw_b(:,:,:), 'VW' ) … … 166 167 fsdepw_n(:,:,1) = 0.0_wp 167 168 fsde3w_n(:,:,1) = fsdept_n(:,:,1) - sshn(:,:) 169 fsdept_b(:,:,1) = 0.5_wp * fse3w_b(:,:,1) 170 fsdepw_b(:,:,1) = 0.0_wp 168 171 DO jk = 2, jpk 169 172 fsdept_n(:,:,jk) = fsdept_n(:,:,jk-1) + fse3w_n(:,:,jk) 170 173 fsdepw_n(:,:,jk) = fsdepw_n(:,:,jk-1) + fse3t_n(:,:,jk-1) 171 174 fsde3w_n(:,:,jk) = fsdept_n(:,:,jk ) - sshn (:,:) 175 fsdept_b(:,:,jk) = fsdept_b(:,:,jk-1) + fse3w_b(:,:,jk) 176 fsdepw_b(:,:,jk) = fsdepw_b(:,:,jk-1) + fse3t_b(:,:,jk-1) 172 177 END DO 173 178 ! Reference water column height at t-, u- and v- point … … 600 605 tilde_e3t_n(:,:,:) = tilde_e3t_a(:,:,:) 601 606 ENDIF 607 fsdept_b(:,:,:) = fsdept_n(:,:,:) 608 fsdepw_b(:,:,:) = fsdepw_n(:,:,:) 609 602 610 fse3t_n(:,:,:) = fse3t_a(:,:,:) 603 611 fse3u_n(:,:,:) = fse3u_a(:,:,:) … … 616 624 CALL dom_vvl_interpol( fse3u_n(:,:,:), fse3uw_n(:,:,:), 'UW' ) 617 625 CALL dom_vvl_interpol( fse3v_n(:,:,:), fse3vw_n(:,:,:), 'VW' ) 626 CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3w_b (:,:,:), 'W' ) 618 627 CALL dom_vvl_interpol( fse3u_b(:,:,:), fse3uw_b(:,:,:), 'UW' ) 619 628 CALL dom_vvl_interpol( fse3v_b(:,:,:), fse3vw_b(:,:,:), 'VW' ) -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr_substitute.h90
r4370 r4488 28 28 # define fse3uw_n(i,j,k) e3uw_n(i,j,k) 29 29 # define fse3vw_n(i,j,k) e3vw_n(i,j,k) 30 31 # define fsdept_b(i,j,k) gdept_b(i,j,k) 32 # define fsdepw_b(i,j,k) gdepw_b(i,j,k) 33 # define fse3w_b(i,j,k) e3w_b(i,j,k) 30 34 31 35 # define fse3t_a(i,j,k) e3t_a(i,j,k) … … 78 82 # define fse3vw_n(i,j,k) e3vw_0(i,j,k) 79 83 84 # define fsdept_b(i,j,k) gdept_0(i,j,k) 85 # define fsdepw_b(i,j,k) gdepw_0(i,j,k) 86 # define fse3w_b(i,j,k) e3w_0(i,j,k) 87 80 88 # define fse3t_a(i,j,k) e3t_0(i,j,k) 81 89 # define fse3u_a(i,j,k) e3u_0(i,j,k) -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf.F90
r3294 r4488 20 20 USE dynldf_iso ! lateral mixing (dyn_ldf_iso routine) 21 21 USE dynldf_lap ! lateral mixing (dyn_ldf_lap routine) 22 USE ldftra_oce, ONLY: ln_traldf_hor ! ocean tracers lateral physics 22 23 USE trdmod ! ocean dynamics and tracer trends 23 24 USE trdmod_oce ! ocean variables trends … … 152 153 IF( ioptio > 1 ) CALL ctl_stop( ' use only ONE direction (level/hor/iso)' ) 153 154 155 IF( ln_dynldf_iso .AND. ln_traldf_hor ) CALL ctl_stop( 'Not sensible to use geopotential diffusion for tracers with isoneutral diffusion for dynamics' ) 156 154 157 ! ! Set nldf, the type of lateral diffusion, from ln_dynldf_... logicals 155 158 ierr = 0 -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilapg.F90
r3634 r4488 19 19 USE dom_oce ! ocean space and time domain 20 20 USE ldfdyn_oce ! ocean dynamics lateral physics 21 USE ldftra_oce, ONLY: ln_traldf_iso 21 22 USE zdf_oce ! ocean vertical physics 22 23 USE trdmod ! ocean dynamics trends … … 29 30 USE wrk_nemo ! Memory Allocation 30 31 USE timing ! Timing 31 32 32 33 33 IMPLICIT NONE … … 104 104 IF( dyn_ldf_bilapg_alloc() /= 0 ) CALL ctl_stop('STOP', 'dyn_ldf_bilapg: failed to allocate arrays') 105 105 ENDIF 106 ! 106 107 ! s-coordinate: Iso-level diffusion on tracer, but geopotential level diffusion on momentum 108 IF( ln_dynldf_hor .AND. ln_traldf_iso ) THEN 109 ! 110 DO jk = 1, jpk ! set the slopes of iso-level 111 DO jj = 2, jpjm1 112 DO ji = 2, jpim1 113 uslp (ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept_b(ji+1,jj,jk) - fsdept_b(ji ,jj ,jk) ) * umask(ji,jj,jk) 114 vslp (ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept_b(ji,jj+1,jk) - fsdept_b(ji ,jj ,jk) ) * vmask(ji,jj,jk) 115 wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw_b(ji+1,jj,jk) - fsdepw_b(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5 116 wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw_b(ji,jj+1,jk) - fsdepw_b(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5 117 END DO 118 END DO 119 END DO 120 ! Lateral boundary conditions on the slopes 121 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) 122 CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 123 124 !!bug 125 IF( kt == nit000 ) then 126 IF(lwp) WRITE(numout,*) ' max slop: u', SQRT( MAXVAL(uslp*uslp)), ' v ', SQRT(MAXVAL(vslp)), & 127 & ' wi', sqrt(MAXVAL(wslpi)) , ' wj', sqrt(MAXVAL(wslpj)) 128 endif 129 !!end 130 ENDIF 131 107 132 zwk1(:,:,:) = 0.e0 ; zwk3(:,:,:) = 0.e0 108 133 zwk2(:,:,:) = 0.e0 ; zwk4(:,:,:) = 0.e0 -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_iso.F90
r3294 r4488 131 131 ENDIF 132 132 133 ! s-coordinate: Iso-level diffusion on momentum but not on tracer133 ! s-coordinate: Iso-level diffusion on tracer, but geopotential level diffusion on momentum 134 134 IF( ln_dynldf_hor .AND. ln_traldf_iso ) THEN 135 135 ! 136 136 DO jk = 1, jpk ! set the slopes of iso-level 137 137 DO jj = 2, jpjm1 138 DO ji = fs_2, fs_jpim1 ! vector opt.139 uslp (ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept (ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * umask(ji,jj,jk)140 vslp (ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept (ji,jj+1,jk) - fsdept(ji ,jj ,jk) ) * vmask(ji,jj,jk)141 wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw (ji+1,jj,jk) - fsdepw(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5142 wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw (ji,jj+1,jk) - fsdepw(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5138 DO ji = 2, jpim1 139 uslp (ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept_b(ji+1,jj,jk) - fsdept_b(ji ,jj ,jk) ) * umask(ji,jj,jk) 140 vslp (ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept_b(ji,jj+1,jk) - fsdept_b(ji ,jj ,jk) ) * vmask(ji,jj,jk) 141 wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw_b(ji+1,jj,jk) - fsdepw_b(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5 142 wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw_b(ji,jj+1,jk) - fsdepw_b(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5 143 143 END DO 144 144 END DO -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/LDF/ldfslp.F90
r3848 r4488 117 117 CALL wrk_alloc( jpi,jpj,jpk, zwz, zww, zdzr, zgru, zgrv ) 118 118 119 zeps = 1.e-20_wp !== Local constant initialization ==! 120 z1_16 = 1.0_wp / 16._wp 121 zm1_g = -1.0_wp / grav 122 zm1_2g = -0.5_wp / grav 123 ! 124 zww(:,:,:) = 0._wp 125 zwz(:,:,:) = 0._wp 126 ! 127 DO jk = 1, jpk !== i- & j-gradient of density ==! 128 DO jj = 1, jpjm1 129 DO ji = 1, fs_jpim1 ! vector opt. 130 zgru(ji,jj,jk) = umask(ji,jj,jk) * ( prd(ji+1,jj ,jk) - prd(ji,jj,jk) ) 131 zgrv(ji,jj,jk) = vmask(ji,jj,jk) * ( prd(ji ,jj+1,jk) - prd(ji,jj,jk) ) 132 END DO 133 END DO 134 END DO 135 IF( ln_zps ) THEN ! partial steps correction at the bottom ocean level 119 IF ( ln_traldf_iso .OR. ln_dynldf_iso ) THEN 120 121 zeps = 1.e-20_wp !== Local constant initialization ==! 122 z1_16 = 1.0_wp / 16._wp 123 zm1_g = -1.0_wp / grav 124 zm1_2g = -0.5_wp / grav 125 ! 126 zww(:,:,:) = 0._wp 127 zwz(:,:,:) = 0._wp 128 ! 129 DO jk = 1, jpk !== i- & j-gradient of density ==! 130 DO jj = 1, jpjm1 131 DO ji = 1, fs_jpim1 ! vector opt. 132 zgru(ji,jj,jk) = umask(ji,jj,jk) * ( prd(ji+1,jj ,jk) - prd(ji,jj,jk) ) 133 zgrv(ji,jj,jk) = vmask(ji,jj,jk) * ( prd(ji ,jj+1,jk) - prd(ji,jj,jk) ) 134 END DO 135 END DO 136 END DO 137 IF( ln_zps ) THEN ! partial steps correction at the bottom ocean level 136 138 # if defined key_vectopt_loop 137 DO jj = 1, 1138 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling)139 DO jj = 1, 1 140 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 139 141 # else 140 DO jj = 1, jpjm1141 DO ji = 1, jpim1142 DO jj = 1, jpjm1 143 DO ji = 1, jpim1 142 144 # endif 143 zgru(ji,jj,mbku(ji,jj)) = gru(ji,jj) 144 zgrv(ji,jj,mbkv(ji,jj)) = grv(ji,jj) 145 END DO 146 END DO 145 zgru(ji,jj,mbku(ji,jj)) = gru(ji,jj) 146 zgrv(ji,jj,mbkv(ji,jj)) = grv(ji,jj) 147 END DO 148 END DO 149 ENDIF 150 ! 151 zdzr(:,:,1) = 0._wp !== Local vertical density gradient at T-point == ! (evaluated from N^2) 152 DO jk = 2, jpkm1 153 ! ! zdzr = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point 154 ! ! trick: tmask(ik ) = 0 => all pn2 = 0 => zdzr = 0 155 ! ! else tmask(ik+1) = 0 => pn2(ik+1) = 0 => zdzr divides by 1 156 ! ! umask(ik+1) /= 0 => all pn2 /= 0 => zdzr divides by 2 157 ! ! NB: 1/(tmask+1) = (1-.5*tmask) substitute a / by a * ==> faster 158 zdzr(:,:,jk) = zm1_g * ( prd(:,:,jk) + 1._wp ) & 159 & * ( pn2(:,:,jk) + pn2(:,:,jk+1) ) * ( 1._wp - 0.5_wp * tmask(:,:,jk+1) ) 160 END DO 161 ! 162 ! !== Slopes just below the mixed layer ==! 163 CALL ldf_slp_mxl( prd, pn2, zgru, zgrv, zdzr ) ! output: uslpml, vslpml, wslpiml, wslpjml 164 165 166 ! I. slopes at u and v point | uslp = d/di( prd ) / d/dz( prd ) 167 ! =========================== | vslp = d/dj( prd ) / d/dz( prd ) 168 ! 169 DO jk = 2, jpkm1 !* Slopes at u and v points 170 DO jj = 2, jpjm1 171 DO ji = fs_2, fs_jpim1 ! vector opt. 172 ! ! horizontal and vertical density gradient at u- and v-points 173 zau = zgru(ji,jj,jk) / e1u(ji,jj) 174 zav = zgrv(ji,jj,jk) / e2v(ji,jj) 175 zbu = 0.5_wp * ( zdzr(ji,jj,jk) + zdzr(ji+1,jj ,jk) ) 176 zbv = 0.5_wp * ( zdzr(ji,jj,jk) + zdzr(ji ,jj+1,jk) ) 177 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0 178 ! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 179 zbu = MIN( zbu, -100._wp* ABS( zau ) , -7.e+3_wp/fse3u(ji,jj,jk)* ABS( zau ) ) 180 zbv = MIN( zbv, -100._wp* ABS( zav ) , -7.e+3_wp/fse3v(ji,jj,jk)* ABS( zav ) ) 181 ! ! uslp and vslp output in zwz and zww, resp. 182 zfi = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) ) 183 zfj = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) ) 184 zwz(ji,jj,jk) = ( ( 1. - zfi) * zau / ( zbu - zeps ) & 185 & + zfi * uslpml(ji,jj) & 186 & * 0.5_wp * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk)-fse3u(ji,jj,1) ) & 187 & / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 5._wp ) ) * umask(ji,jj,jk) 188 zww(ji,jj,jk) = ( ( 1. - zfj) * zav / ( zbv - zeps ) & 189 & + zfj * vslpml(ji,jj) & 190 & * 0.5_wp * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk)-fse3v(ji,jj,1) ) & 191 & / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 5. ) ) * vmask(ji,jj,jk) 192 !!gm modif to suppress omlmask.... (as in Griffies case) 193 ! ! ! jk must be >= ML level for zf=1. otherwise zf=0. 194 ! zfi = REAL( 1 - 1/(1 + jk / MAX( nmln(ji+1,jj), nmln(ji,jj) ) ), wp ) 195 ! zfj = REAL( 1 - 1/(1 + jk / MAX( nmln(ji,jj+1), nmln(ji,jj) ) ), wp ) 196 ! zci = 0.5 * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk) ) / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 10. ) ) 197 ! zcj = 0.5 * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk) ) / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 10. ) ) 198 ! zwz(ji,jj,jk) = ( zfi * zai / ( zbi - zeps ) + ( 1._wp - zfi ) * wslpiml(ji,jj) * zci ) * tmask(ji,jj,jk) 199 ! zww(ji,jj,jk) = ( zfj * zaj / ( zbj - zeps ) + ( 1._wp - zfj ) * wslpjml(ji,jj) * zcj ) * tmask(ji,jj,jk) 200 !!gm end modif 201 END DO 202 END DO 203 END DO 204 CALL lbc_lnk( zwz, 'U', -1. ) ; CALL lbc_lnk( zww, 'V', -1. ) ! lateral boundary conditions 205 ! 206 ! !* horizontal Shapiro filter 207 DO jk = 2, jpkm1 208 DO jj = 2, jpjm1, MAX(1, jpj-3) ! rows jj=2 and =jpjm1 only 209 DO ji = 2, jpim1 210 uslp(ji,jj,jk) = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 211 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 212 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 213 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 214 & + 4.* zwz(ji ,jj ,jk) ) 215 vslp(ji,jj,jk) = z1_16 * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 216 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 217 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 218 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 219 & + 4.* zww(ji,jj ,jk) ) 220 END DO 221 END DO 222 DO jj = 3, jpj-2 ! other rows 223 DO ji = fs_2, fs_jpim1 ! vector opt. 224 uslp(ji,jj,jk) = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 225 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 226 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 227 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 228 & + 4.* zwz(ji ,jj ,jk) ) 229 vslp(ji,jj,jk) = z1_16 * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 230 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 231 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 232 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 233 & + 4.* zww(ji,jj ,jk) ) 234 END DO 235 END DO 236 ! !* decrease along coastal boundaries 237 DO jj = 2, jpjm1 238 DO ji = fs_2, fs_jpim1 ! vector opt. 239 uslp(ji,jj,jk) = uslp(ji,jj,jk) * ( umask(ji,jj+1,jk) + umask(ji,jj-1,jk ) ) * 0.5_wp & 240 & * ( umask(ji,jj ,jk) + umask(ji,jj ,jk+1) ) * 0.5_wp 241 vslp(ji,jj,jk) = vslp(ji,jj,jk) * ( vmask(ji+1,jj,jk) + vmask(ji-1,jj,jk ) ) * 0.5_wp & 242 & * ( vmask(ji ,jj,jk) + vmask(ji ,jj,jk+1) ) * 0.5_wp 243 END DO 244 END DO 245 END DO 246 247 248 ! II. slopes at w point | wslpi = mij( d/di( prd ) / d/dz( prd ) 249 ! =========================== | wslpj = mij( d/dj( prd ) / d/dz( prd ) 250 ! 251 DO jk = 2, jpkm1 252 DO jj = 2, jpjm1 253 DO ji = fs_2, fs_jpim1 ! vector opt. 254 ! !* Local vertical density gradient evaluated from N^2 255 zbw = zm1_2g * pn2 (ji,jj,jk) * ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) 256 ! !* Slopes at w point 257 ! ! i- & j-gradient of density at w-points 258 zci = MAX( umask(ji-1,jj,jk ) + umask(ji,jj,jk ) & 259 & + umask(ji-1,jj,jk-1) + umask(ji,jj,jk-1) , zeps ) * e1t(ji,jj) 260 zcj = MAX( vmask(ji,jj-1,jk ) + vmask(ji,jj,jk-1) & 261 & + vmask(ji,jj-1,jk-1) + vmask(ji,jj,jk ) , zeps ) * e2t(ji,jj) 262 zai = ( zgru (ji-1,jj,jk ) + zgru (ji,jj,jk-1) & 263 & + zgru (ji-1,jj,jk-1) + zgru (ji,jj,jk ) ) / zci * tmask (ji,jj,jk) 264 zaj = ( zgrv (ji,jj-1,jk ) + zgrv (ji,jj,jk-1) & 265 & + zgrv (ji,jj-1,jk-1) + zgrv (ji,jj,jk ) ) / zcj * tmask (ji,jj,jk) 266 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0. 267 ! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 268 zbi = MIN( zbw ,- 100._wp* ABS( zai ) , -7.e+3_wp/fse3w(ji,jj,jk)* ABS( zai ) ) 269 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/fse3w(ji,jj,jk)* ABS( zaj ) ) 270 ! ! wslpi and wslpj with ML flattening (output in zwz and zww, resp.) 271 zfk = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) ) ! zfk=1 in the ML otherwise zfk=0 272 zck = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10._wp ) 273 zwz(ji,jj,jk) = ( zai / ( zbi - zeps ) * ( 1._wp - zfk ) + zck * wslpiml(ji,jj) * zfk ) * tmask(ji,jj,jk) 274 zww(ji,jj,jk) = ( zaj / ( zbj - zeps ) * ( 1._wp - zfk ) + zck * wslpjml(ji,jj) * zfk ) * tmask(ji,jj,jk) 275 276 !!gm modif to suppress omlmask.... (as in Griffies operator) 277 ! ! ! jk must be >= ML level for zfk=1. otherwise zfk=0. 278 ! zfk = REAL( 1 - 1/(1 + jk / nmln(ji+1,jj)), wp ) 279 ! zck = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10. ) 280 ! zwz(ji,jj,jk) = ( zfk * zai / ( zbi - zeps ) + ( 1._wp - zfk ) * wslpiml(ji,jj) * zck ) * tmask(ji,jj,jk) 281 ! zww(ji,jj,jk) = ( zfk * zaj / ( zbj - zeps ) + ( 1._wp - zfk ) * wslpjml(ji,jj) * zck ) * tmask(ji,jj,jk) 282 !!gm end modif 283 END DO 284 END DO 285 END DO 286 CALL lbc_lnk( zwz, 'T', -1. ) ; CALL lbc_lnk( zww, 'T', -1. ) ! lateral boundary conditions 287 ! 288 ! !* horizontal Shapiro filter 289 DO jk = 2, jpkm1 290 DO jj = 2, jpjm1, MAX(1, jpj-3) ! rows jj=2 and =jpjm1 only 291 DO ji = 2, jpim1 292 zcofw = tmask(ji,jj,jk) * z1_16 293 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 294 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 295 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 296 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 297 & + 4.* zwz(ji ,jj ,jk) ) * zcofw 298 299 wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 300 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 301 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 302 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 303 & + 4.* zww(ji ,jj ,jk) ) * zcofw 304 END DO 305 END DO 306 DO jj = 3, jpj-2 ! other rows 307 DO ji = fs_2, fs_jpim1 ! vector opt. 308 zcofw = tmask(ji,jj,jk) * z1_16 309 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 310 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 311 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 312 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 313 & + 4.* zwz(ji ,jj ,jk) ) * zcofw 314 315 wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 316 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 317 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 318 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 319 & + 4.* zww(ji ,jj ,jk) ) * zcofw 320 END DO 321 END DO 322 ! !* decrease along coastal boundaries 323 DO jj = 2, jpjm1 324 DO ji = fs_2, fs_jpim1 ! vector opt. 325 zck = ( umask(ji,jj,jk) + umask(ji-1,jj,jk) ) & 326 & * ( vmask(ji,jj,jk) + vmask(ji,jj-1,jk) ) * 0.25 327 wslpi(ji,jj,jk) = wslpi(ji,jj,jk) * zck 328 wslpj(ji,jj,jk) = wslpj(ji,jj,jk) * zck 329 END DO 330 END DO 331 END DO 332 333 ! III. Specific grid points 334 ! =========================== 335 ! 336 IF( cp_cfg == "orca" .AND. jp_cfg == 4 ) THEN ! ORCA_R4 configuration: horizontal diffusion in specific area 337 ! ! Gibraltar Strait 338 ij0 = 50 ; ij1 = 53 339 ii0 = 69 ; ii1 = 71 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 340 ij0 = 51 ; ij1 = 53 341 ii0 = 68 ; ii1 = 71 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 342 ii0 = 69 ; ii1 = 71 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 343 ii0 = 69 ; ii1 = 71 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 344 ! 345 ! ! Mediterrannean Sea 346 ij0 = 49 ; ij1 = 56 347 ii0 = 71 ; ii1 = 90 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 348 ij0 = 50 ; ij1 = 56 349 ii0 = 70 ; ii1 = 90 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 350 ii0 = 71 ; ii1 = 90 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 351 ii0 = 71 ; ii1 = 90 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 352 ENDIF 353 354 355 ! IV. Lateral boundary conditions 356 ! =============================== 357 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) 358 CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 359 360 361 IF(ln_ctl) THEN 362 CALL prt_ctl(tab3d_1=uslp , clinfo1=' slp - u : ', tab3d_2=vslp, clinfo2=' v : ', kdim=jpk) 363 CALL prt_ctl(tab3d_1=wslpi, clinfo1=' slp - wi: ', tab3d_2=wslpj, clinfo2=' wj: ', kdim=jpk) 364 ENDIF 365 ! 366 367 ELSEIF ( lk_vvl ) THEN 368 369 IF(lwp) THEN 370 WRITE(numout,*) ' Horizontal mixing in s-coordinate: slope = slope of s-surfaces' 371 ENDIF 372 373 ! geopotential diffusion in s-coordinates on tracers and/or momentum 374 ! The slopes of s-surfaces are computed at each time step due to vvl 375 ! The slopes for momentum diffusion are i- or j- averaged of those on tracers 376 377 ! set the slope of diffusion to the slope of s-surfaces 378 ! ( c a u t i o n : minus sign as fsdep has positive value ) 379 DO jk = 1, jpk 380 DO jj = 2, jpjm1 381 DO ji = fs_2, fs_jpim1 ! vector opt. 382 uslp(ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept_b(ji+1,jj,jk) - fsdept_b(ji ,jj ,jk) ) * umask(ji,jj,jk) 383 vslp(ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept_b(ji,jj+1,jk) - fsdept_b(ji ,jj ,jk) ) * vmask(ji,jj,jk) 384 wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw_b(ji+1,jj,jk) - fsdepw_b(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5 385 wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw_b(ji,jj+1,jk) - fsdepw_b(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5 386 END DO 387 END DO 388 END DO 389 390 ! Lateral boundary conditions on the slopes 391 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) 392 CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 393 394 if( kt == nit000 ) then 395 IF(lwp) WRITE(numout,*) ' max slop: u',SQRT( MAXVAL(uslp*uslp)), ' v ', SQRT(MAXVAL(vslp)), & 396 & ' wi', sqrt(MAXVAL(wslpi)), ' wj', sqrt(MAXVAL(wslpj)) 397 endif 398 399 IF(ln_ctl) THEN 400 CALL prt_ctl(tab3d_1=uslp , clinfo1=' slp - u : ', tab3d_2=vslp, clinfo2=' v : ', kdim=jpk) 401 CALL prt_ctl(tab3d_1=wslpi, clinfo1=' slp - wi: ', tab3d_2=wslpj, clinfo2=' wj: ', kdim=jpk) 402 ENDIF 403 147 404 ENDIF 148 ! 149 zdzr(:,:,1) = 0._wp !== Local vertical density gradient at T-point == ! (evaluated from N^2) 150 DO jk = 2, jpkm1 151 ! ! zdzr = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point 152 ! ! trick: tmask(ik ) = 0 => all pn2 = 0 => zdzr = 0 153 ! ! else tmask(ik+1) = 0 => pn2(ik+1) = 0 => zdzr divides by 1 154 ! ! umask(ik+1) /= 0 => all pn2 /= 0 => zdzr divides by 2 155 ! ! NB: 1/(tmask+1) = (1-.5*tmask) substitute a / by a * ==> faster 156 zdzr(:,:,jk) = zm1_g * ( prd(:,:,jk) + 1._wp ) & 157 & * ( pn2(:,:,jk) + pn2(:,:,jk+1) ) * ( 1._wp - 0.5_wp * tmask(:,:,jk+1) ) 158 END DO 159 ! 160 ! !== Slopes just below the mixed layer ==! 161 CALL ldf_slp_mxl( prd, pn2, zgru, zgrv, zdzr ) ! output: uslpml, vslpml, wslpiml, wslpjml 162 163 164 ! I. slopes at u and v point | uslp = d/di( prd ) / d/dz( prd ) 165 ! =========================== | vslp = d/dj( prd ) / d/dz( prd ) 166 ! 167 DO jk = 2, jpkm1 !* Slopes at u and v points 168 DO jj = 2, jpjm1 169 DO ji = fs_2, fs_jpim1 ! vector opt. 170 ! ! horizontal and vertical density gradient at u- and v-points 171 zau = zgru(ji,jj,jk) / e1u(ji,jj) 172 zav = zgrv(ji,jj,jk) / e2v(ji,jj) 173 zbu = 0.5_wp * ( zdzr(ji,jj,jk) + zdzr(ji+1,jj ,jk) ) 174 zbv = 0.5_wp * ( zdzr(ji,jj,jk) + zdzr(ji ,jj+1,jk) ) 175 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0 176 ! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 177 zbu = MIN( zbu, -100._wp* ABS( zau ) , -7.e+3_wp/fse3u(ji,jj,jk)* ABS( zau ) ) 178 zbv = MIN( zbv, -100._wp* ABS( zav ) , -7.e+3_wp/fse3v(ji,jj,jk)* ABS( zav ) ) 179 ! ! uslp and vslp output in zwz and zww, resp. 180 zfi = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) ) 181 zfj = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) ) 182 zwz(ji,jj,jk) = ( ( 1. - zfi) * zau / ( zbu - zeps ) & 183 & + zfi * uslpml(ji,jj) & 184 & * 0.5_wp * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk)-fse3u(ji,jj,1) ) & 185 & / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 5._wp ) ) * umask(ji,jj,jk) 186 zww(ji,jj,jk) = ( ( 1. - zfj) * zav / ( zbv - zeps ) & 187 & + zfj * vslpml(ji,jj) & 188 & * 0.5_wp * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk)-fse3v(ji,jj,1) ) & 189 & / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 5. ) ) * vmask(ji,jj,jk) 190 !!gm modif to suppress omlmask.... (as in Griffies case) 191 ! ! ! jk must be >= ML level for zf=1. otherwise zf=0. 192 ! zfi = REAL( 1 - 1/(1 + jk / MAX( nmln(ji+1,jj), nmln(ji,jj) ) ), wp ) 193 ! zfj = REAL( 1 - 1/(1 + jk / MAX( nmln(ji,jj+1), nmln(ji,jj) ) ), wp ) 194 ! zci = 0.5 * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk) ) / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 10. ) ) 195 ! zcj = 0.5 * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk) ) / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 10. ) ) 196 ! zwz(ji,jj,jk) = ( zfi * zai / ( zbi - zeps ) + ( 1._wp - zfi ) * wslpiml(ji,jj) * zci ) * tmask(ji,jj,jk) 197 ! zww(ji,jj,jk) = ( zfj * zaj / ( zbj - zeps ) + ( 1._wp - zfj ) * wslpjml(ji,jj) * zcj ) * tmask(ji,jj,jk) 198 !!gm end modif 199 END DO 200 END DO 201 END DO 202 CALL lbc_lnk( zwz, 'U', -1. ) ; CALL lbc_lnk( zww, 'V', -1. ) ! lateral boundary conditions 203 ! 204 ! !* horizontal Shapiro filter 205 DO jk = 2, jpkm1 206 DO jj = 2, jpjm1, MAX(1, jpj-3) ! rows jj=2 and =jpjm1 only 207 DO ji = 2, jpim1 208 uslp(ji,jj,jk) = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 209 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 210 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 211 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 212 & + 4.* zwz(ji ,jj ,jk) ) 213 vslp(ji,jj,jk) = z1_16 * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 214 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 215 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 216 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 217 & + 4.* zww(ji,jj ,jk) ) 218 END DO 219 END DO 220 DO jj = 3, jpj-2 ! other rows 221 DO ji = fs_2, fs_jpim1 ! vector opt. 222 uslp(ji,jj,jk) = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 223 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 224 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 225 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 226 & + 4.* zwz(ji ,jj ,jk) ) 227 vslp(ji,jj,jk) = z1_16 * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 228 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 229 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 230 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 231 & + 4.* zww(ji,jj ,jk) ) 232 END DO 233 END DO 234 ! !* decrease along coastal boundaries 235 DO jj = 2, jpjm1 236 DO ji = fs_2, fs_jpim1 ! vector opt. 237 uslp(ji,jj,jk) = uslp(ji,jj,jk) * ( umask(ji,jj+1,jk) + umask(ji,jj-1,jk ) ) * 0.5_wp & 238 & * ( umask(ji,jj ,jk) + umask(ji,jj ,jk+1) ) * 0.5_wp 239 vslp(ji,jj,jk) = vslp(ji,jj,jk) * ( vmask(ji+1,jj,jk) + vmask(ji-1,jj,jk ) ) * 0.5_wp & 240 & * ( vmask(ji ,jj,jk) + vmask(ji ,jj,jk+1) ) * 0.5_wp 241 END DO 242 END DO 243 END DO 244 245 246 ! II. slopes at w point | wslpi = mij( d/di( prd ) / d/dz( prd ) 247 ! =========================== | wslpj = mij( d/dj( prd ) / d/dz( prd ) 248 ! 249 DO jk = 2, jpkm1 250 DO jj = 2, jpjm1 251 DO ji = fs_2, fs_jpim1 ! vector opt. 252 ! !* Local vertical density gradient evaluated from N^2 253 zbw = zm1_2g * pn2 (ji,jj,jk) * ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) 254 ! !* Slopes at w point 255 ! ! i- & j-gradient of density at w-points 256 zci = MAX( umask(ji-1,jj,jk ) + umask(ji,jj,jk ) & 257 & + umask(ji-1,jj,jk-1) + umask(ji,jj,jk-1) , zeps ) * e1t(ji,jj) 258 zcj = MAX( vmask(ji,jj-1,jk ) + vmask(ji,jj,jk-1) & 259 & + vmask(ji,jj-1,jk-1) + vmask(ji,jj,jk ) , zeps ) * e2t(ji,jj) 260 zai = ( zgru (ji-1,jj,jk ) + zgru (ji,jj,jk-1) & 261 & + zgru (ji-1,jj,jk-1) + zgru (ji,jj,jk ) ) / zci * tmask (ji,jj,jk) 262 zaj = ( zgrv (ji,jj-1,jk ) + zgrv (ji,jj,jk-1) & 263 & + zgrv (ji,jj-1,jk-1) + zgrv (ji,jj,jk ) ) / zcj * tmask (ji,jj,jk) 264 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0. 265 ! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 266 zbi = MIN( zbw ,- 100._wp* ABS( zai ) , -7.e+3_wp/fse3w(ji,jj,jk)* ABS( zai ) ) 267 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/fse3w(ji,jj,jk)* ABS( zaj ) ) 268 ! ! wslpi and wslpj with ML flattening (output in zwz and zww, resp.) 269 zfk = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) ) ! zfk=1 in the ML otherwise zfk=0 270 zck = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10._wp ) 271 zwz(ji,jj,jk) = ( zai / ( zbi - zeps ) * ( 1._wp - zfk ) + zck * wslpiml(ji,jj) * zfk ) * tmask(ji,jj,jk) 272 zww(ji,jj,jk) = ( zaj / ( zbj - zeps ) * ( 1._wp - zfk ) + zck * wslpjml(ji,jj) * zfk ) * tmask(ji,jj,jk) 273 274 !!gm modif to suppress omlmask.... (as in Griffies operator) 275 ! ! ! jk must be >= ML level for zfk=1. otherwise zfk=0. 276 ! zfk = REAL( 1 - 1/(1 + jk / nmln(ji+1,jj)), wp ) 277 ! zck = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10. ) 278 ! zwz(ji,jj,jk) = ( zfk * zai / ( zbi - zeps ) + ( 1._wp - zfk ) * wslpiml(ji,jj) * zck ) * tmask(ji,jj,jk) 279 ! zww(ji,jj,jk) = ( zfk * zaj / ( zbj - zeps ) + ( 1._wp - zfk ) * wslpjml(ji,jj) * zck ) * tmask(ji,jj,jk) 280 !!gm end modif 281 END DO 282 END DO 283 END DO 284 CALL lbc_lnk( zwz, 'T', -1. ) ; CALL lbc_lnk( zww, 'T', -1. ) ! lateral boundary conditions 285 ! 286 ! !* horizontal Shapiro filter 287 DO jk = 2, jpkm1 288 DO jj = 2, jpjm1, MAX(1, jpj-3) ! rows jj=2 and =jpjm1 only 289 DO ji = 2, jpim1 290 zcofw = tmask(ji,jj,jk) * z1_16 291 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 292 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 293 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 294 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 295 & + 4.* zwz(ji ,jj ,jk) ) * zcofw 296 297 wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 298 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 299 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 300 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 301 & + 4.* zww(ji ,jj ,jk) ) * zcofw 302 END DO 303 END DO 304 DO jj = 3, jpj-2 ! other rows 305 DO ji = fs_2, fs_jpim1 ! vector opt. 306 zcofw = tmask(ji,jj,jk) * z1_16 307 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 308 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 309 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 310 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 311 & + 4.* zwz(ji ,jj ,jk) ) * zcofw 312 313 wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 314 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 315 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 316 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 317 & + 4.* zww(ji ,jj ,jk) ) * zcofw 318 END DO 319 END DO 320 ! !* decrease along coastal boundaries 321 DO jj = 2, jpjm1 322 DO ji = fs_2, fs_jpim1 ! vector opt. 323 zck = ( umask(ji,jj,jk) + umask(ji-1,jj,jk) ) & 324 & * ( vmask(ji,jj,jk) + vmask(ji,jj-1,jk) ) * 0.25 325 wslpi(ji,jj,jk) = wslpi(ji,jj,jk) * zck 326 wslpj(ji,jj,jk) = wslpj(ji,jj,jk) * zck 327 END DO 328 END DO 329 END DO 330 331 ! III. Specific grid points 332 ! =========================== 333 ! 334 IF( cp_cfg == "orca" .AND. jp_cfg == 4 ) THEN ! ORCA_R4 configuration: horizontal diffusion in specific area 335 ! ! Gibraltar Strait 336 ij0 = 50 ; ij1 = 53 337 ii0 = 69 ; ii1 = 71 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 338 ij0 = 51 ; ij1 = 53 339 ii0 = 68 ; ii1 = 71 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 340 ii0 = 69 ; ii1 = 71 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 341 ii0 = 69 ; ii1 = 71 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 342 ! 343 ! ! Mediterrannean Sea 344 ij0 = 49 ; ij1 = 56 345 ii0 = 71 ; ii1 = 90 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 346 ij0 = 50 ; ij1 = 56 347 ii0 = 70 ; ii1 = 90 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 348 ii0 = 71 ; ii1 = 90 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 349 ii0 = 71 ; ii1 = 90 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 350 ENDIF 351 352 353 ! IV. Lateral boundary conditions 354 ! =============================== 355 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) 356 CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 357 358 359 IF(ln_ctl) THEN 360 CALL prt_ctl(tab3d_1=uslp , clinfo1=' slp - u : ', tab3d_2=vslp, clinfo2=' v : ', kdim=jpk) 361 CALL prt_ctl(tab3d_1=wslpi, clinfo1=' slp - wi: ', tab3d_2=wslpj, clinfo2=' wj: ', kdim=jpk) 362 ENDIF 363 ! 405 364 406 CALL wrk_dealloc( jpi,jpj,jpk, zwz, zww, zdzr, zgru, zgrv ) 365 407 ! … … 783 825 wslpj(:,:,:) = 0._wp ; wslpjml(:,:) = 0._wp 784 826 785 !!gm I no longer understand this..... 786 IF( (ln_traldf_hor .OR. ln_dynldf_hor) .AND. .NOT. (lk_vvl .AND. ln_rstart) ) THEN 827 IF( ln_traldf_hor .OR. ln_dynldf_hor ) THEN 787 828 IF(lwp) WRITE(numout,*) ' Horizontal mixing in s-coordinate: slope = slope of s-surfaces' 788 829 … … 796 837 DO jj = 2, jpjm1 797 838 DO ji = fs_2, fs_jpim1 ! vector opt. 798 uslp (ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept (ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * umask(ji,jj,jk)799 vslp (ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept (ji,jj+1,jk) - fsdept(ji ,jj ,jk) ) * vmask(ji,jj,jk)800 wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw (ji+1,jj,jk) - fsdepw(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5801 wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw (ji,jj+1,jk) - fsdepw(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5839 uslp (ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept_b(ji+1,jj,jk) - fsdept_b(ji ,jj ,jk) ) * umask(ji,jj,jk) 840 vslp (ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept_b(ji,jj+1,jk) - fsdept_b(ji ,jj ,jk) ) * vmask(ji,jj,jk) 841 wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw_b(ji+1,jj,jk) - fsdepw_b(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5 842 wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw_b(ji,jj+1,jk) - fsdepw_b(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5 802 843 END DO 803 844 END DO -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/TRA/traldf.F90
r3294 r4488 204 204 ENDIF 205 205 206 IF( nldf == 3 ) CALL ctl_warn( 'geopotential bilaplacian tracer diffusion in s-coords not thoroughly tested' ) 206 207 IF( ierr == 1 ) CALL ctl_stop( ' iso-level in z-coordinate - partial step, not allowed' ) 207 208 IF( ierr == 2 ) CALL ctl_stop( ' isoneutral bilaplacian operator does not exist' )
Note: See TracChangeset
for help on using the changeset viewer.
