- Timestamp:
- 2016-01-08T10:35:19+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OPA_SRC/LDF/ldfslp.F90
r4488 r6225 11 11 !! 3.3 ! 2010-10 (G. Nurser, C. Harris, G. Madec) add Griffies operator 12 12 !! - ! 2010-11 (F. Dupond, G. Madec) bug correction in slopes just below the ML 13 !! 3.7 ! 2013-12 (F. Lemarie, G. Madec) add limiter on triad slopes 13 14 !!---------------------------------------------------------------------- 14 #if defined key_ldfslp || defined key_esopa 15 15 16 !!---------------------------------------------------------------------- 16 !! 'key_ldfslp' Rotation of lateral mixing tensor17 !!----------------------------------------------------------------------18 !! ldf_slp_grif : calculates the triads of isoneutral slopes (Griffies operator)19 17 !! ldf_slp : calculates the slopes of neutral surface (Madec operator) 18 !! ldf_slp_triad : calculates the triads of isoneutral slopes (Griffies operator) 20 19 !! ldf_slp_mxl : calculates the slopes at the base of the mixed layer (Madec operator) 21 20 !! ldf_slp_init : initialization of the slopes computation … … 23 22 USE oce ! ocean dynamics and tracers 24 23 USE dom_oce ! ocean space and time domain 25 USE ldftra_oce ! lateral diffusion: traceur 26 USE ldfdyn_oce ! lateral diffusion: dynamics 24 USE ldfdyn ! lateral diffusion: eddy viscosity coef. 27 25 USE phycst ! physical constants 28 26 USE zdfmxl ! mixed layer depth 29 27 USE eosbn2 ! equation of states 30 USE lbclnk ! ocean lateral boundary conditions (or mpp link)28 ! 31 29 USE in_out_manager ! I/O manager 32 30 USE prtctl ! Print control 31 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 32 USE lib_mpp ! distribued memory computing library 33 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 33 34 USE wrk_nemo ! work arrays 34 35 USE timing ! Timing 35 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)36 36 37 37 IMPLICIT NONE 38 38 PRIVATE 39 39 40 PUBLIC ldf_slp ! routine called by step.F90 41 PUBLIC ldf_slp_grif ! routine called by step.F90 42 PUBLIC ldf_slp_init ! routine called by opa.F90 43 44 LOGICAL , PUBLIC, PARAMETER :: lk_ldfslp = .TRUE. !: slopes flag 45 ! !! Madec operator 46 ! Arrays allocated in ldf_slp_init() routine once we know whether we're using the Griffies or Madec operator 40 PUBLIC ldf_slp ! routine called by step.F90 41 PUBLIC ldf_slp_triad ! routine called by step.F90 42 PUBLIC ldf_slp_init ! routine called by nemogcm.F90 43 44 LOGICAL , PUBLIC :: l_ldfslp = .FALSE. !: slopes flag 45 46 LOGICAL , PUBLIC :: ln_traldf_iso = .TRUE. !: iso-neutral direction 47 LOGICAL , PUBLIC :: ln_traldf_triad = .FALSE. !: griffies triad scheme 48 49 LOGICAL , PUBLIC :: ln_triad_iso = .FALSE. !: pure horizontal mixing in ML 50 LOGICAL , PUBLIC :: ln_botmix_triad = .FALSE. !: mixing on bottom 51 REAL(wp), PUBLIC :: rn_sw_triad = 1._wp !: =1 switching triads ; =0 all four triads used 52 REAL(wp), PUBLIC :: rn_slpmax = 0.01_wp !: slope limit 53 54 LOGICAL , PUBLIC :: l_grad_zps = .FALSE. !: special treatment for Horz Tgradients w partial steps (triad operator) 55 56 ! !! Classic operator (Madec) 47 57 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: uslp, wslpi !: i_slope at U- and W-points 48 58 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: vslp, wslpj !: j-slope at V- and W-points 49 ! !! Griffies operator59 ! !! triad operator (Griffies) 50 60 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: wslp2 !: wslp**2 from Griffies quarter cells 51 61 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:,:) :: triadi_g, triadj_g !: skew flux slopes relative to geopotentials 52 62 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:,:) :: triadi , triadj !: isoneutral slopes relative to model-coordinate 53 54 ! !! Madec operator 55 ! Arrays allocated in ldf_slp_init() routine once we know whether we're using the Griffies or Madec operator 63 ! !! both operators 64 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ah_wslp2 !: ah * slope^2 at w-point 65 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: akz !: stabilizing vertical diffusivity 66 67 ! !! Madec operator 56 68 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: omlmask ! mask of the surface mixed layer at T-pt 57 69 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: uslpml, wslpiml ! i_slope at U- and W-points just below the mixed layer … … 61 73 62 74 !! * Substitutions 63 # include "domzgr_substitute.h90"64 # include "ldftra_substitute.h90"65 # include "ldfeiv_substitute.h90"66 75 # include "vectopt_loop_substitute.h90" 67 76 !!---------------------------------------------------------------------- 68 !! NEMO/OPA 4.0 , NEMO Consortium (201 1)77 !! NEMO/OPA 4.0 , NEMO Consortium (2014) 69 78 !! $Id$ 70 79 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 102 111 !! 103 112 INTEGER :: ji , jj , jk ! dummy loop indices 104 INTEGER :: ii0, ii1 , iku! temporary integer105 INTEGER :: ij0, ij1 , ikv! temporary integer106 REAL(wp) :: zeps, zm1_g, zm1_2g, z1_16, zcofw ! local scalars113 INTEGER :: ii0, ii1 ! temporary integer 114 INTEGER :: ij0, ij1 ! temporary integer 115 REAL(wp) :: zeps, zm1_g, zm1_2g, z1_16, zcofw, z1_slpmax ! local scalars 107 116 REAL(wp) :: zci, zfi, zau, zbu, zai, zbi ! - - 108 117 REAL(wp) :: zcj, zfj, zav, zbv, zaj, zbj ! - - 109 118 REAL(wp) :: zck, zfk, zbw ! - - 119 REAL(wp) :: zdepu, zdepv ! - - 120 REAL(wp), POINTER, DIMENSION(:,: ) :: zslpml_hmlpu, zslpml_hmlpv 110 121 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwz, zww 111 122 REAL(wp), POINTER, DIMENSION(:,:,:) :: zdzr … … 116 127 ! 117 128 CALL wrk_alloc( jpi,jpj,jpk, zwz, zww, zdzr, zgru, zgrv ) 118 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 138 # if defined key_vectopt_loop 139 DO jj = 1, 1 140 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 141 # else 142 DO jj = 1, jpjm1 143 DO ji = 1, jpim1 144 # endif 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) 129 CALL wrk_alloc( jpi,jpj, zslpml_hmlpu, zslpml_hmlpv ) 130 131 zeps = 1.e-20_wp !== Local constant initialization ==! 132 z1_16 = 1.0_wp / 16._wp 133 zm1_g = -1.0_wp / grav 134 zm1_2g = -0.5_wp / grav 135 z1_slpmax = 1._wp / rn_slpmax 136 ! 137 zww(:,:,:) = 0._wp 138 zwz(:,:,:) = 0._wp 139 ! 140 DO jk = 1, jpk !== i- & j-gradient of density ==! 141 DO jj = 1, jpjm1 142 DO ji = 1, fs_jpim1 ! vector opt. 143 zgru(ji,jj,jk) = umask(ji,jj,jk) * ( prd(ji+1,jj ,jk) - prd(ji,jj,jk) ) 144 zgrv(ji,jj,jk) = vmask(ji,jj,jk) * ( prd(ji ,jj+1,jk) - prd(ji,jj,jk) ) 145 END DO 146 END DO 147 END DO 148 IF( ln_zps ) THEN ! partial steps correction at the bottom ocean level 149 DO jj = 1, jpjm1 150 DO ji = 1, jpim1 151 zgru(ji,jj,mbku(ji,jj)) = gru(ji,jj) 152 zgrv(ji,jj,mbkv(ji,jj)) = grv(ji,jj) 153 END DO 154 END DO 155 ENDIF 156 IF( ln_zps .AND. ln_isfcav ) THEN ! partial steps correction at the bottom ocean level 157 DO jj = 1, jpjm1 158 DO ji = 1, jpim1 159 IF ( miku(ji,jj) > 1 ) zgru(ji,jj,miku(ji,jj)) = grui(ji,jj) 160 IF ( mikv(ji,jj) > 1 ) zgrv(ji,jj,mikv(ji,jj)) = grvi(ji,jj) 161 END DO 162 END DO 163 ENDIF 164 ! 165 zdzr(:,:,1) = 0._wp !== Local vertical density gradient at T-point == ! (evaluated from N^2) 166 DO jk = 2, jpkm1 167 ! ! zdzr = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point 168 ! ! trick: tmask(ik ) = 0 => all pn2 = 0 => zdzr = 0 169 ! ! else tmask(ik+1) = 0 => pn2(ik+1) = 0 => zdzr divides by 1 170 ! ! umask(ik+1) /= 0 => all pn2 /= 0 => zdzr divides by 2 171 ! ! NB: 1/(tmask+1) = (1-.5*tmask) substitute a / by a * ==> faster 172 zdzr(:,:,jk) = zm1_g * ( prd(:,:,jk) + 1._wp ) & 173 & * ( pn2(:,:,jk) + pn2(:,:,jk+1) ) * ( 1._wp - 0.5_wp * tmask(:,:,jk+1) ) 174 END DO 175 ! 176 ! !== Slopes just below the mixed layer ==! 177 CALL ldf_slp_mxl( prd, pn2, zgru, zgrv, zdzr ) ! output: uslpml, vslpml, wslpiml, wslpjml 178 179 180 ! I. slopes at u and v point | uslp = d/di( prd ) / d/dz( prd ) 181 ! =========================== | vslp = d/dj( prd ) / d/dz( prd ) 182 ! 183 IF ( ln_isfcav ) THEN 184 DO jj = 2, jpjm1 185 DO ji = fs_2, fs_jpim1 ! vector opt. 186 zslpml_hmlpu(ji,jj) = uslpml(ji,jj) / ( MAX(hmlpt(ji,jj), hmlpt(ji+1,jj ), 5._wp) & 187 & - 0.5_wp * ( risfdep(ji,jj) + risfdep(ji+1,jj ) ) ) 188 zslpml_hmlpv(ji,jj) = vslpml(ji,jj) / ( MAX(hmlpt(ji,jj), hmlpt(ji ,jj+1), 5._wp) & 189 & - 0.5_wp * ( risfdep(ji,jj) + risfdep(ji ,jj+1) ) ) 190 END DO 191 END DO 192 ELSE 193 DO jj = 2, jpjm1 194 DO ji = fs_2, fs_jpim1 ! vector opt. 195 zslpml_hmlpu(ji,jj) = uslpml(ji,jj) / MAX(hmlpt(ji,jj), hmlpt(ji+1,jj ), 5._wp) 196 zslpml_hmlpv(ji,jj) = vslpml(ji,jj) / MAX(hmlpt(ji,jj), hmlpt(ji ,jj+1), 5._wp) 197 END DO 198 END DO 199 END IF 200 201 DO jk = 2, jpkm1 !* Slopes at u and v points 202 DO jj = 2, jpjm1 203 DO ji = fs_2, fs_jpim1 ! vector opt. 204 ! ! horizontal and vertical density gradient at u- and v-points 205 zau = zgru(ji,jj,jk) * r1_e1u(ji,jj) 206 zav = zgrv(ji,jj,jk) * r1_e2v(ji,jj) 207 zbu = 0.5_wp * ( zdzr(ji,jj,jk) + zdzr(ji+1,jj ,jk) ) 208 zbv = 0.5_wp * ( zdzr(ji,jj,jk) + zdzr(ji ,jj+1,jk) ) 209 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0 210 ! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 211 zbu = MIN( zbu, - z1_slpmax * ABS( zau ) , -7.e+3_wp/e3u_n(ji,jj,jk)* ABS( zau ) ) 212 zbv = MIN( zbv, - z1_slpmax * ABS( zav ) , -7.e+3_wp/e3v_n(ji,jj,jk)* ABS( zav ) ) 213 ! ! uslp and vslp output in zwz and zww, resp. 214 zfi = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) ) 215 zfj = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) ) 216 ! thickness of water column between surface and level k at u/v point 217 zdepu = 0.5_wp * ( ( gdept_n (ji,jj,jk) + gdept_n (ji+1,jj ,jk) ) & 218 - ( risfdep(ji,jj) + risfdep(ji+1,jj) ) - e3u_n(ji,jj,miku(ji,jj)) ) 219 zdepv = 0.5_wp * ( ( gdept_n (ji,jj,jk) + gdept_n (ji,jj+1,jk) ) & 220 - ( risfdep(ji,jj) + risfdep(ji,jj+1) ) - e3v_n(ji,jj,mikv(ji,jj)) ) 221 ! 222 zwz(ji,jj,jk) = ( ( 1._wp - zfi) * zau / ( zbu - zeps ) & 223 & + zfi * zdepu * zslpml_hmlpu(ji,jj) ) * umask(ji,jj,jk) 224 zww(ji,jj,jk) = ( ( 1._wp - zfj) * zav / ( zbv - zeps ) & 225 & + zfj * zdepv * zslpml_hmlpv(ji,jj) ) * vmask(ji,jj,jk) 192 226 !!gm modif to suppress omlmask.... (as in Griffies case) 193 ! 194 ! 195 ! 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 ! 199 ! 227 ! ! ! jk must be >= ML level for zf=1. otherwise zf=0. 228 ! zfi = REAL( 1 - 1/(1 + jk / MAX( nmln(ji+1,jj), nmln(ji,jj) ) ), wp ) 229 ! zfj = REAL( 1 - 1/(1 + jk / MAX( nmln(ji,jj+1), nmln(ji,jj) ) ), wp ) 230 ! zci = 0.5 * ( gdept_n(ji+1,jj,jk)+gdept_n(ji,jj,jk) ) / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 10. ) ) 231 ! zcj = 0.5 * ( gdept_n(ji,jj+1,jk)+gdept_n(ji,jj,jk) ) / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 10. ) ) 232 ! zwz(ji,jj,jk) = ( zfi * zai / ( zbi - zeps ) + ( 1._wp - zfi ) * wslpiml(ji,jj) * zci ) * tmask(ji,jj,jk) 233 ! zww(ji,jj,jk) = ( zfj * zaj / ( zbj - zeps ) + ( 1._wp - zfj ) * wslpjml(ji,jj) * zcj ) * tmask(ji,jj,jk) 200 234 !!gm end modif 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk)&232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 & + vmask(ji,jj-1,jk-1) + vmask(ji,jj,jk ) , zeps ) *e2t(ji,jj)262 263 & + zgru (ji-1,jj,jk-1) + zgru (ji,jj,jk ) ) / zci * tmask (ji,jj,jk)264 265 & + zgrv (ji,jj-1,jk-1) + zgrv (ji,jj,jk ) ) / zcj * tmask (ji,jj,jk)266 267 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 271 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)235 END DO 236 END DO 237 END DO 238 CALL lbc_lnk( zwz, 'U', -1. ) ; CALL lbc_lnk( zww, 'V', -1. ) ! lateral boundary conditions 239 ! 240 ! !* horizontal Shapiro filter 241 DO jk = 2, jpkm1 242 DO jj = 2, jpjm1, MAX(1, jpj-3) ! rows jj=2 and =jpjm1 only 243 DO ji = 2, jpim1 244 uslp(ji,jj,jk) = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 245 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 246 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 247 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 248 & + 4.* zwz(ji ,jj ,jk) ) 249 vslp(ji,jj,jk) = z1_16 * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 250 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 251 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 252 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 253 & + 4.* zww(ji,jj ,jk) ) 254 END DO 255 END DO 256 DO jj = 3, jpj-2 ! other rows 257 DO ji = fs_2, fs_jpim1 ! vector opt. 258 uslp(ji,jj,jk) = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 259 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 260 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 261 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 262 & + 4.* zwz(ji ,jj ,jk) ) 263 vslp(ji,jj,jk) = z1_16 * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 264 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 265 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 266 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 267 & + 4.* zww(ji,jj ,jk) ) 268 END DO 269 END DO 270 ! !* decrease along coastal boundaries 271 DO jj = 2, jpjm1 272 DO ji = fs_2, fs_jpim1 ! vector opt. 273 uslp(ji,jj,jk) = uslp(ji,jj,jk) * ( umask(ji,jj+1,jk) + umask(ji,jj-1,jk ) ) * 0.5_wp & 274 & * ( umask(ji,jj ,jk) + umask(ji,jj ,jk+1) ) * 0.5_wp 275 vslp(ji,jj,jk) = vslp(ji,jj,jk) * ( vmask(ji+1,jj,jk) + vmask(ji-1,jj,jk ) ) * 0.5_wp & 276 & * ( vmask(ji ,jj,jk) + vmask(ji ,jj,jk+1) ) * 0.5_wp 277 END DO 278 END DO 279 END DO 280 281 282 ! II. slopes at w point | wslpi = mij( d/di( prd ) / d/dz( prd ) 283 ! =========================== | wslpj = mij( d/dj( prd ) / d/dz( prd ) 284 ! 285 DO jk = 2, jpkm1 286 DO jj = 2, jpjm1 287 DO ji = fs_2, fs_jpim1 ! vector opt. 288 ! !* Local vertical density gradient evaluated from N^2 289 zbw = zm1_2g * pn2 (ji,jj,jk) * ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) 290 ! !* Slopes at w point 291 ! ! i- & j-gradient of density at w-points 292 zci = MAX( umask(ji-1,jj,jk ) + umask(ji,jj,jk ) & 293 & + umask(ji-1,jj,jk-1) + umask(ji,jj,jk-1) , zeps ) * e1t(ji,jj) 294 zcj = MAX( vmask(ji,jj-1,jk ) + vmask(ji,jj,jk-1) & 295 & + vmask(ji,jj-1,jk-1) + vmask(ji,jj,jk ) , zeps ) * e2t(ji,jj) 296 zai = ( zgru (ji-1,jj,jk ) + zgru (ji,jj,jk-1) & 297 & + zgru (ji-1,jj,jk-1) + zgru (ji,jj,jk ) ) / zci * wmask (ji,jj,jk) 298 zaj = ( zgrv (ji,jj-1,jk ) + zgrv (ji,jj,jk-1) & 299 & + zgrv (ji,jj-1,jk-1) + zgrv (ji,jj,jk ) ) / zcj * wmask (ji,jj,jk) 300 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0. 301 ! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 302 zbi = MIN( zbw ,- 100._wp* ABS( zai ) , -7.e+3_wp/e3w_n(ji,jj,jk)* ABS( zai ) ) 303 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/e3w_n(ji,jj,jk)* ABS( zaj ) ) 304 ! ! wslpi and wslpj with ML flattening (output in zwz and zww, resp.) 305 zfk = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) ) ! zfk=1 in the ML otherwise zfk=0 306 zck = ( gdepw_n(ji,jj,jk) - gdepw_n(ji,jj,mikt(ji,jj) ) ) / MAX( hmlp(ji,jj) - gdepw_n(ji,jj,mikt(ji,jj)), 10._wp ) 307 zwz(ji,jj,jk) = ( zai / ( zbi - zeps ) * ( 1._wp - zfk ) + zck * wslpiml(ji,jj) * zfk ) * wmask(ji,jj,jk) 308 zww(ji,jj,jk) = ( zaj / ( zbj - zeps ) * ( 1._wp - zfk ) + zck * wslpjml(ji,jj) * zfk ) * wmask(ji,jj,jk) 275 309 276 310 !!gm modif to suppress omlmask.... (as in Griffies operator) 277 ! 278 ! 279 ! zck = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10. )280 ! 281 ! 311 ! ! ! jk must be >= ML level for zfk=1. otherwise zfk=0. 312 ! zfk = REAL( 1 - 1/(1 + jk / nmln(ji+1,jj)), wp ) 313 ! zck = gdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10. ) 314 ! zwz(ji,jj,jk) = ( zfk * zai / ( zbi - zeps ) + ( 1._wp - zfk ) * wslpiml(ji,jj) * zck ) * tmask(ji,jj,jk) 315 ! zww(ji,jj,jk) = ( zfk * zaj / ( zbj - zeps ) + ( 1._wp - zfk ) * wslpjml(ji,jj) * zck ) * tmask(ji,jj,jk) 282 316 !!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 317 END DO 318 END DO 319 END DO 320 CALL lbc_lnk( zwz, 'T', -1. ) ; CALL lbc_lnk( zww, 'T', -1. ) ! lateral boundary conditions 321 ! 322 ! !* horizontal Shapiro filter 323 DO jk = 2, jpkm1 324 DO jj = 2, jpjm1, MAX(1, jpj-3) ! rows jj=2 and =jpjm1 only 325 DO ji = 2, jpim1 326 zcofw = wmask(ji,jj,jk) * z1_16 327 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 328 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 329 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 330 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 331 & + 4.* zwz(ji ,jj ,jk) ) * zcofw 332 333 wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 334 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 335 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 336 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 337 & + 4.* zww(ji ,jj ,jk) ) * zcofw 338 END DO 339 END DO 340 DO jj = 3, jpj-2 ! other rows 341 DO ji = fs_2, fs_jpim1 ! vector opt. 342 zcofw = wmask(ji,jj,jk) * z1_16 343 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 344 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 345 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 346 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 347 & + 4.* zwz(ji ,jj ,jk) ) * zcofw 348 349 wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 350 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 351 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 352 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 353 & + 4.* zww(ji ,jj ,jk) ) * zcofw 354 END DO 355 END DO 356 ! !* decrease in vicinity of topography 357 DO jj = 2, jpjm1 358 DO ji = fs_2, fs_jpim1 ! vector opt. 359 zck = ( umask(ji,jj,jk) + umask(ji-1,jj,jk) ) & 360 & * ( vmask(ji,jj,jk) + vmask(ji,jj-1,jk) ) * 0.25 361 wslpi(ji,jj,jk) = wslpi(ji,jj,jk) * zck 362 wslpj(ji,jj,jk) = wslpj(ji,jj,jk) * zck 363 END DO 364 END DO 365 END DO 366 367 ! IV. Lateral boundary conditions 368 ! =============================== 369 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) 370 CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 371 372 IF(ln_ctl) THEN 373 CALL prt_ctl(tab3d_1=uslp , clinfo1=' slp - u : ', tab3d_2=vslp, clinfo2=' v : ', kdim=jpk) 374 CALL prt_ctl(tab3d_1=wslpi, clinfo1=' slp - wi: ', tab3d_2=wslpj, clinfo2=' wj: ', kdim=jpk) 404 375 ENDIF 405 376 ! 406 377 CALL wrk_dealloc( jpi,jpj,jpk, zwz, zww, zdzr, zgru, zgrv ) 378 CALL wrk_dealloc( jpi,jpj, zslpml_hmlpu, zslpml_hmlpv ) 407 379 ! 408 380 IF( nn_timing == 1 ) CALL timing_stop('ldf_slp') … … 411 383 412 384 413 SUBROUTINE ldf_slp_ grif( kt )414 !!---------------------------------------------------------------------- 415 !! *** ROUTINE ldf_slp_ grif***385 SUBROUTINE ldf_slp_triad ( kt ) 386 !!---------------------------------------------------------------------- 387 !! *** ROUTINE ldf_slp_triad *** 416 388 !! 417 389 !! ** Purpose : Compute the squared slopes of neutral surfaces (slope 418 !! of iso-pycnal surfaces referenced locally) (ln_traldf_ grif=T)390 !! of iso-pycnal surfaces referenced locally) (ln_traldf_triad=T) 419 391 !! at W-points using the Griffies quarter-cells. 420 392 !! … … 431 403 REAL(wp) :: zfacti, zfactj ! local scalars 432 404 REAL(wp) :: znot_thru_surface ! local scalars 433 REAL(wp) :: zdit, zdis, zdjt, zdjs, zdkt, zdks, zbu, zbv, zbti, zbtj 405 REAL(wp) :: zdit, zdis, zdkt, zbu, zbti, zisw 406 REAL(wp) :: zdjt, zdjs, zdks, zbv, zbtj, zjsw 434 407 REAL(wp) :: zdxrho_raw, zti_coord, zti_raw, zti_lim, zti_g_raw, zti_g_lim 435 408 REAL(wp) :: zdyrho_raw, ztj_coord, ztj_raw, ztj_lim, ztj_g_raw, ztj_g_lim 436 409 REAL(wp) :: zdzrho_raw 437 REAL(wp) :: zbeta0 410 REAL(wp) :: zbeta0, ze3_e1, ze3_e2 438 411 REAL(wp), POINTER, DIMENSION(:,:) :: z1_mlbw 439 412 REAL(wp), POINTER, DIMENSION(:,:,:) :: zalbet … … 442 415 !!---------------------------------------------------------------------- 443 416 ! 444 IF( nn_timing == 1 ) CALL timing_start('ldf_slp_ grif')417 IF( nn_timing == 1 ) CALL timing_start('ldf_slp_triad') 445 418 ! 446 419 CALL wrk_alloc( jpi,jpj, z1_mlbw ) … … 452 425 ! Some preliminary calculation ! 453 426 !--------------------------------! 454 !455 CALL eos_alpbet( tsb, zalbet, zbeta0 ) !== before local thermal/haline expension ratio at T-points ==!456 427 ! 457 428 DO jl = 0, 1 !== unmasked before density i- j-, k-gradients ==! … … 465 436 zdjt = ( tsb(ji,jj+1,jk,jp_tem) - tsb(ji,jj,jk,jp_tem) ) ! j-gradient of T & S at v-point 466 437 zdjs = ( tsb(ji,jj+1,jk,jp_sal) - tsb(ji,jj,jk,jp_sal) ) 467 zdxrho_raw = ( - zalbet(ji+ip,jj ,jk) * zdit + zbeta0*zdis ) /e1u(ji,jj)468 zdyrho_raw = ( - zalbet(ji ,jj+jp,jk) * zdjt + zbeta0*zdjs ) /e2v(ji,jj)469 zdxrho(ji+ip,jj ,jk,1-ip) = SIGN( MAX( repsln, ABS( zdxrho_raw ) ), zdxrho_raw) ! keep the sign470 zdyrho(ji ,jj+jp,jk,1-jp) = SIGN( MAX( repsln, ABS( zdyrho_raw ) ), zdyrho_raw)438 zdxrho_raw = ( - rab_b(ji+ip,jj ,jk,jp_tem) * zdit + rab_b(ji+ip,jj ,jk,jp_sal) * zdis ) * r1_e1u(ji,jj) 439 zdyrho_raw = ( - rab_b(ji ,jj+jp,jk,jp_tem) * zdjt + rab_b(ji ,jj+jp,jk,jp_sal) * zdjs ) * r1_e2v(ji,jj) 440 zdxrho(ji+ip,jj ,jk,1-ip) = SIGN( MAX( repsln, ABS( zdxrho_raw ) ), zdxrho_raw ) ! keep the sign 441 zdyrho(ji ,jj+jp,jk,1-jp) = SIGN( MAX( repsln, ABS( zdyrho_raw ) ), zdyrho_raw ) 471 442 END DO 472 443 END DO 473 444 END DO 474 445 ! 475 IF( ln_zps.and.l_grad_zps ) THEN ! partial steps: correction of i- & j-grad on bottom 476 # if defined key_vectopt_loop 477 DO jj = 1, 1 478 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 479 # else 446 IF( ln_zps .AND. l_grad_zps ) THEN ! partial steps: correction of i- & j-grad on bottom 480 447 DO jj = 1, jpjm1 481 448 DO ji = 1, jpim1 482 # endif483 449 iku = mbku(ji,jj) ; ikv = mbkv(ji,jj) ! last ocean level (u- & v-points) 484 450 zdit = gtsu(ji,jj,jp_tem) ; zdjt = gtsv(ji,jj,jp_tem) ! i- & j-gradient of Temperature 485 451 zdis = gtsu(ji,jj,jp_sal) ; zdjs = gtsv(ji,jj,jp_sal) ! i- & j-gradient of Salinity 486 zdxrho_raw = ( - zalbet(ji+ip,jj ,iku) * zdit + zbeta0*zdis ) /e1u(ji,jj)487 zdyrho_raw = ( - zalbet(ji ,jj+jp,ikv) * zdjt + zbeta0*zdjs ) /e2v(ji,jj)452 zdxrho_raw = ( - rab_b(ji+ip,jj ,iku,jp_tem) * zdit + rab_b(ji+ip,jj ,iku,jp_sal) * zdis ) * r1_e1u(ji,jj) 453 zdyrho_raw = ( - rab_b(ji ,jj+jp,ikv,jp_tem) * zdjt + rab_b(ji ,jj+jp,ikv,jp_sal) * zdjs ) * r1_e2v(ji,jj) 488 454 zdxrho(ji+ip,jj ,iku,1-ip) = SIGN( MAX( repsln, ABS( zdxrho_raw ) ), zdxrho_raw ) ! keep the sign 489 455 zdyrho(ji ,jj+jp,ikv,1-jp) = SIGN( MAX( repsln, ABS( zdyrho_raw ) ), zdyrho_raw ) … … 505 471 zdks = 0._wp 506 472 ENDIF 507 zdzrho_raw = ( - zalbet(ji ,jj ,jk) * zdkt + zbeta0*zdks ) / fse3w(ji,jj,jk+kp)508 zdzrho(ji ,jj ,jk, kp) = - MIN( - repsln,zdzrho_raw ) ! force zdzrho >= repsln473 zdzrho_raw = ( - zalbet(ji,jj,jk) * zdkt + zbeta0*zdks ) / e3w_n(ji,jj,jk+kp) 474 zdzrho(ji,jj,jk,kp) = - MIN( - repsln , zdzrho_raw ) ! force zdzrho >= repsln 509 475 END DO 510 476 END DO … … 515 481 DO ji = 1, jpi 516 482 jk = MIN( nmln(ji,jj), mbkt(ji,jj) ) + 1 ! MIN in case ML depth is the ocean depth 517 z1_mlbw(ji,jj) = 1._wp / fsdepw(ji,jj,jk)483 z1_mlbw(ji,jj) = 1._wp / gdepw_n(ji,jj,jk) 518 484 END DO 519 485 END DO … … 539 505 ! 540 506 jk = nmln(ji+ip,jj) + 1 541 IF( jk .GT. mbkt(ji+ip,jj) ) THEN !ML reaches bottom 542 zti_mlb(ji+ip,jj ,1-ip,kp) = 0.0_wp 543 ELSE 544 ! Add s-coordinate slope at t-points (do this by *subtracting* gradient of depth) 545 zti_g_raw = ( zdxrho(ji+ip,jj,jk-kp,1-ip) / zdzrho(ji+ip,jj,jk-kp,kp) & 546 & - ( fsdept(ji+1,jj,jk-kp) - fsdept(ji,jj,jk-kp) ) / e1u(ji,jj) ) * umask(ji,jj,jk) 547 zti_mlb(ji+ip,jj ,1-ip,kp) = SIGN( MIN( rn_slpmax, ABS( zti_g_raw ) ), zti_g_raw ) 507 IF( jk > mbkt(ji+ip,jj) ) THEN ! ML reaches bottom 508 zti_mlb(ji+ip,jj ,1-ip,kp) = 0.0_wp 509 ELSE 510 ! Add s-coordinate slope at t-points (do this by *subtracting* gradient of depth) 511 zti_g_raw = ( zdxrho(ji+ip,jj,jk-kp,1-ip) / zdzrho(ji+ip,jj,jk-kp,kp) & 512 & - ( gdept_n(ji+1,jj,jk-kp) - gdept_n(ji,jj,jk-kp) ) * r1_e1u(ji,jj) ) * umask(ji,jj,jk) 513 ze3_e1 = e3w_n(ji+ip,jj,jk-kp) * r1_e1u(ji,jj) 514 zti_mlb(ji+ip,jj ,1-ip,kp) = SIGN( MIN( rn_slpmax, 5.0_wp * ze3_e1 , ABS( zti_g_raw ) ), zti_g_raw ) 548 515 ENDIF 549 516 ! 550 517 jk = nmln(ji,jj+jp) + 1 551 IF( jk .GT.mbkt(ji,jj+jp) ) THEN !ML reaches bottom552 ztj_mlb(ji ,jj+jp,1-jp,kp) = 0.0_wp518 IF( jk > mbkt(ji,jj+jp) ) THEN !ML reaches bottom 519 ztj_mlb(ji ,jj+jp,1-jp,kp) = 0.0_wp 553 520 ELSE 554 ztj_g_raw = ( zdyrho(ji,jj+jp,jk-kp,1-jp) / zdzrho(ji,jj+jp,jk-kp,kp) & 555 & - ( fsdept(ji,jj+1,jk-kp) - fsdept(ji,jj,jk-kp) ) / e2v(ji,jj) ) * vmask(ji,jj,jk) 556 ztj_mlb(ji ,jj+jp,1-jp,kp) = SIGN( MIN( rn_slpmax, ABS( ztj_g_raw ) ), ztj_g_raw ) 521 ztj_g_raw = ( zdyrho(ji,jj+jp,jk-kp,1-jp) / zdzrho(ji,jj+jp,jk-kp,kp) & 522 & - ( gdept_n(ji,jj+1,jk-kp) - gdept_n(ji,jj,jk-kp) ) / e2v(ji,jj) ) * vmask(ji,jj,jk) 523 ze3_e2 = e3w_n(ji,jj+jp,jk-kp) / e2v(ji,jj) 524 ztj_mlb(ji ,jj+jp,1-jp,kp) = SIGN( MIN( rn_slpmax, 5.0_wp * ze3_e2 , ABS( ztj_g_raw ) ), ztj_g_raw ) 557 525 ENDIF 558 526 END DO … … 583 551 zti_raw = zdxrho(ji+ip,jj ,jk,1-ip) / zdzrho(ji+ip,jj ,jk,kp) ! unmasked 584 552 ztj_raw = zdyrho(ji ,jj+jp,jk,1-jp) / zdzrho(ji ,jj+jp,jk,kp) 585 553 ! 586 554 ! Must mask contribution to slope for triad jk=1,kp=0 that poke up though ocean surface 587 zti_coord = znot_thru_surface * ( fsdept(ji+1,jj ,jk) - fsdept(ji,jj,jk) ) /e1u(ji,jj)588 ztj_coord = znot_thru_surface * ( fsdept(ji ,jj+1,jk) - fsdept(ji,jj,jk) ) / e2v(ji,jj)! unmasked555 zti_coord = znot_thru_surface * ( gdept_n(ji+1,jj ,jk) - gdept_n(ji,jj,jk) ) * r1_e1u(ji,jj) 556 ztj_coord = znot_thru_surface * ( gdept_n(ji ,jj+1,jk) - gdept_n(ji,jj,jk) ) * r1_e2v(ji,jj) ! unmasked 589 557 zti_g_raw = zti_raw - zti_coord ! ref to geopot surfaces 590 558 ztj_g_raw = ztj_raw - ztj_coord 591 zti_g_lim = SIGN( MIN( rn_slpmax, ABS( zti_g_raw ) ), zti_g_raw ) 592 ztj_g_lim = SIGN( MIN( rn_slpmax, ABS( ztj_g_raw ) ), ztj_g_raw ) 559 ! additional limit required in bilaplacian case 560 ze3_e1 = e3w_n(ji+ip,jj ,jk+kp) * r1_e1u(ji,jj) 561 ze3_e2 = e3w_n(ji ,jj+jp,jk+kp) * r1_e2v(ji,jj) 562 ! NB: hard coded factor 5 (can be a namelist parameter...) 563 zti_g_lim = SIGN( MIN( rn_slpmax, 5.0_wp * ze3_e1, ABS( zti_g_raw ) ), zti_g_raw ) 564 ztj_g_lim = SIGN( MIN( rn_slpmax, 5.0_wp * ze3_e2, ABS( ztj_g_raw ) ), ztj_g_raw ) 593 565 ! 594 566 ! Below ML use limited zti_g as is & mask … … 600 572 zti_g_lim = ( zfacti * zti_g_lim & 601 573 & + ( 1._wp - zfacti ) * zti_mlb(ji+ip,jj,1-ip,kp) & 602 & * fsdepw(ji+ip,jj,jk+kp) * z1_mlbw(ji+ip,jj) ) * umask(ji,jj,jk+kp)574 & * gdepw_n(ji+ip,jj,jk+kp) * z1_mlbw(ji+ip,jj) ) * umask(ji,jj,jk+kp) 603 575 ztj_g_lim = ( zfactj * ztj_g_lim & 604 576 & + ( 1._wp - zfactj ) * ztj_mlb(ji,jj+jp,1-jp,kp) & 605 & * fsdepw(ji,jj+jp,jk+kp) * z1_mlbw(ji,jj+jp) ) * vmask(ji,jj,jk+kp)577 & * gdepw_n(ji,jj+jp,jk+kp) * z1_mlbw(ji,jj+jp) ) * vmask(ji,jj,jk+kp) 606 578 ! 607 579 triadi_g(ji+ip,jj ,jk,1-ip,kp) = zti_g_lim … … 619 591 ! 620 592 IF( ln_triad_iso ) THEN 621 zti_raw = zti_lim* *2/ zti_raw622 ztj_raw = ztj_lim* *2/ ztj_raw593 zti_raw = zti_lim*zti_lim / zti_raw 594 ztj_raw = ztj_lim*ztj_lim / ztj_raw 623 595 zti_raw = SIGN( MIN( ABS(zti_lim), ABS( zti_raw ) ), zti_raw ) 624 596 ztj_raw = SIGN( MIN( ABS(ztj_lim), ABS( ztj_raw ) ), ztj_raw ) 625 zti_lim = zfacti * zti_lim & 626 & + ( 1._wp - zfacti ) * zti_raw 627 ztj_lim = zfactj * ztj_lim & 628 & + ( 1._wp - zfactj ) * ztj_raw 597 zti_lim = zfacti * zti_lim + ( 1._wp - zfacti ) * zti_raw 598 ztj_lim = zfactj * ztj_lim + ( 1._wp - zfactj ) * ztj_raw 629 599 ENDIF 630 triadi(ji+ip,jj ,jk,1-ip,kp) = zti_lim 631 triadj(ji ,jj+jp,jk,1-jp,kp) = ztj_lim 632 ! 633 zbu = e1u(ji ,jj) * e2u(ji ,jj) * fse3u(ji ,jj,jk ) 634 zbv = e1v(ji ,jj) * e2v(ji ,jj) * fse3v(ji ,jj,jk ) 635 zbti = e1t(ji+ip,jj) * e2t(ji+ip,jj) * fse3w(ji+ip,jj,jk+kp) 636 zbtj = e1t(ji,jj+jp) * e2t(ji,jj+jp) * fse3w(ji,jj+jp,jk+kp) 637 ! 638 !!gm this may inhibit vectorization on Vect Computers, and even on scalar computers.... ==> to be checked 639 wslp2 (ji+ip,jj,jk+kp) = wslp2(ji+ip,jj,jk+kp) + 0.25_wp * zbu / zbti * zti_g_lim**2 ! masked 640 wslp2 (ji,jj+jp,jk+kp) = wslp2(ji,jj+jp,jk+kp) + 0.25_wp * zbv / zbtj * ztj_g_lim**2 600 ! ! switching triad scheme 601 zisw = (rn_sw_triad - 1._wp ) + rn_sw_triad & 602 & * 2._wp * ABS( 0.5_wp - kp - ( 0.5_wp - ip ) * SIGN( 1._wp , zdxrho(ji+ip,jj,jk,1-ip) ) ) 603 zjsw = (rn_sw_triad - 1._wp ) + rn_sw_triad & 604 & * 2._wp * ABS( 0.5_wp - kp - ( 0.5_wp - jp ) * SIGN( 1._wp , zdyrho(ji,jj+jp,jk,1-jp) ) ) 605 ! 606 triadi(ji+ip,jj ,jk,1-ip,kp) = zti_lim * zisw 607 triadj(ji ,jj+jp,jk,1-jp,kp) = ztj_lim * zjsw 608 ! 609 zbu = e1e2u(ji ,jj ) * e3u_n(ji ,jj ,jk ) 610 zbv = e1e2v(ji ,jj ) * e3v_n(ji ,jj ,jk ) 611 zbti = e1e2t(ji+ip,jj ) * e3w_n(ji+ip,jj ,jk+kp) 612 zbtj = e1e2t(ji ,jj+jp) * e3w_n(ji ,jj+jp,jk+kp) 613 ! 614 wslp2(ji+ip,jj,jk+kp) = wslp2(ji+ip,jj,jk+kp) + 0.25_wp * zbu / zbti * zti_g_lim*zti_g_lim ! masked 615 wslp2(ji,jj+jp,jk+kp) = wslp2(ji,jj+jp,jk+kp) + 0.25_wp * zbv / zbtj * ztj_g_lim*ztj_g_lim 641 616 END DO 642 617 END DO … … 654 629 CALL wrk_dealloc( jpi,jpj, 2,2, zti_mlb, ztj_mlb, kkstart = 0, klstart = 0 ) 655 630 ! 656 IF( nn_timing == 1 ) CALL timing_stop('ldf_slp_ grif')657 ! 658 END SUBROUTINE ldf_slp_ grif631 IF( nn_timing == 1 ) CALL timing_stop('ldf_slp_triad') 632 ! 633 END SUBROUTINE ldf_slp_triad 659 634 660 635 … … 682 657 INTEGER :: ji , jj , jk ! dummy loop indices 683 658 INTEGER :: iku, ikv, ik, ikm1 ! local integers 684 REAL(wp) :: zeps, zm1_g, zm1_2g 659 REAL(wp) :: zeps, zm1_g, zm1_2g, z1_slpmax ! local scalars 685 660 REAL(wp) :: zci, zfi, zau, zbu, zai, zbi ! - - 686 661 REAL(wp) :: zcj, zfj, zav, zbv, zaj, zbj ! - - … … 693 668 zm1_g = -1.0_wp / grav 694 669 zm1_2g = -0.5_wp / grav 670 z1_slpmax = 1._wp / rn_slpmax 695 671 ! 696 672 uslpml (1,:) = 0._wp ; uslpml (jpi,:) = 0._wp … … 701 677 ! !== surface mixed layer mask ! 702 678 DO jk = 1, jpk ! =1 inside the mixed layer, =0 otherwise 703 # if defined key_vectopt_loop704 DO jj = 1, 1705 DO ji = 1, jpij ! vector opt. (forced unrolling)706 # else707 679 DO jj = 1, jpj 708 680 DO ji = 1, jpi 709 # endif710 681 ik = nmln(ji,jj) - 1 711 682 IF( jk <= ik ) THEN ; omlmask(ji,jj,jk) = 1._wp … … 727 698 !----------------------------------------------------------------------- 728 699 ! 729 # if defined key_vectopt_loop730 DO jj = 1, 1731 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling)732 # else733 700 DO jj = 2, jpjm1 734 701 DO ji = 2, jpim1 735 # endif736 702 ! !== Slope at u- & v-points just below the Mixed Layer ==! 737 703 ! … … 742 708 zbv = 0.5_wp * ( p_dzr(ji,jj,ikv) + p_dzr(ji ,jj+1,ikv) ) 743 709 ! !- horizontal density gradient at u- & v-points 744 zau = p_gru(ji,jj,iku) /e1u(ji,jj)745 zav = p_grv(ji,jj,ikv) /e2v(ji,jj)710 zau = p_gru(ji,jj,iku) * r1_e1u(ji,jj) 711 zav = p_grv(ji,jj,ikv) * r1_e2v(ji,jj) 746 712 ! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0 747 713 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 748 zbu = MIN( zbu , - 100._wp* ABS( zau ) , -7.e+3_wp/fse3u(ji,jj,iku)* ABS( zau ) )749 zbv = MIN( zbv , - 100._wp* ABS( zav ) , -7.e+3_wp/fse3v(ji,jj,ikv)* ABS( zav ) )714 zbu = MIN( zbu , - z1_slpmax * ABS( zau ) , -7.e+3_wp/e3u_n(ji,jj,iku)* ABS( zau ) ) 715 zbv = MIN( zbv , - z1_slpmax * ABS( zav ) , -7.e+3_wp/e3v_n(ji,jj,ikv)* ABS( zav ) ) 750 716 ! !- Slope at u- & v-points (uslpml, vslpml) 751 717 uslpml(ji,jj) = zau / ( zbu - zeps ) * umask(ji,jj,iku) … … 763 729 zcj = MAX( vmask(ji,jj-1,ik ) + vmask(ji,jj,ik ) & 764 730 & + vmask(ji,jj-1,ikm1) + vmask(ji,jj,ikm1) , zeps ) * e2t(ji,jj) 765 zai = ( p_gru(ji-1,jj,ik ) + p_gru(ji,jj,ik) &731 zai = ( p_gru(ji-1,jj,ik ) + p_gru(ji,jj,ik) & 766 732 & + p_gru(ji-1,jj,ikm1) + p_gru(ji,jj,ikm1 ) ) / zci * tmask(ji,jj,ik) 767 733 zaj = ( p_grv(ji,jj-1,ik ) + p_grv(ji,jj,ik ) & … … 769 735 ! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0. 770 736 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 771 zbi = MIN( zbw , -100._wp* ABS( zai ) , -7.e+3_wp/ fse3w(ji,jj,ik)* ABS( zai ) )772 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/ fse3w(ji,jj,ik)* ABS( zaj ) )737 zbi = MIN( zbw , -100._wp* ABS( zai ) , -7.e+3_wp/e3w_n(ji,jj,ik)* ABS( zai ) ) 738 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/e3w_n(ji,jj,ik)* ABS( zaj ) ) 773 739 ! !- i- & j-slope at w-points (wslpiml, wslpjml) 774 740 wslpiml(ji,jj) = zai / ( zbi - zeps ) * tmask (ji,jj,ik) … … 791 757 !! ** Purpose : Initialization for the isopycnal slopes computation 792 758 !! 793 !! ** Method : read the nammbf namelist and check the parameter 794 !! values called by tra_dmp at the first timestep (nit000) 759 !! ** Method : 795 760 !!---------------------------------------------------------------------- 796 761 INTEGER :: ji, jj, jk ! dummy loop indices … … 805 770 WRITE(numout,*) '~~~~~~~~~~~~' 806 771 ENDIF 807 808 IF( ln_traldf_grif ) THEN ! Griffies operator : triad of slopes 809 ALLOCATE( triadi_g(jpi,jpj,jpk,0:1,0:1) , triadj_g(jpi,jpj,jpk,0:1,0:1) , wslp2(jpi,jpj,jpk) , STAT=ierr ) 810 ALLOCATE( triadi (jpi,jpj,jpk,0:1,0:1) , triadj (jpi,jpj,jpk,0:1,0:1) , STAT=ierr ) 811 IF( ierr > 0 ) CALL ctl_stop( 'STOP', 'ldf_slp_init : unable to allocate Griffies operator slope' ) 812 ! 772 ! 773 ALLOCATE( ah_wslp2(jpi,jpj,jpk) , akz(jpi,jpj,jpk) , STAT=ierr ) 774 IF( ierr > 0 ) CALL ctl_stop( 'STOP', 'ldf_slp_init : unable to allocate ah_slp2 or akz' ) 775 ! 776 IF( ln_traldf_triad ) THEN ! Griffies operator : triad of slopes 777 IF(lwp) WRITE(numout,*) ' Griffies (triad) operator initialisation' 778 ALLOCATE( triadi_g(jpi,jpj,jpk,0:1,0:1) , triadj_g(jpi,jpj,jpk,0:1,0:1) , & 779 & triadi (jpi,jpj,jpk,0:1,0:1) , triadj (jpi,jpj,jpk,0:1,0:1) , & 780 & wslp2 (jpi,jpj,jpk) , STAT=ierr ) 781 IF( ierr > 0 ) CALL ctl_stop( 'STOP', 'ldf_slp_init : unable to allocate Griffies operator slope' ) 813 782 IF( ln_dynldf_iso ) CALL ctl_stop( 'ldf_slp_init: Griffies operator on momentum not supported' ) 814 783 ! 815 784 ELSE ! Madec operator : slopes at u-, v-, and w-points 816 ALLOCATE( uslp(jpi,jpj,jpk) , vslp(jpi,jpj,jpk) , wslpi(jpi,jpj,jpk) , wslpj(jpi,jpj,jpk) , & 817 & omlmask(jpi,jpj,jpk) , uslpml(jpi,jpj) , vslpml(jpi,jpj) , wslpiml(jpi,jpj) , wslpjml(jpi,jpj) , STAT=ierr ) 785 IF(lwp) WRITE(numout,*) ' Madec operator initialisation' 786 ALLOCATE( omlmask(jpi,jpj,jpk) , & 787 & uslp(jpi,jpj,jpk) , uslpml(jpi,jpj) , wslpi(jpi,jpj,jpk) , wslpiml(jpi,jpj) , & 788 & vslp(jpi,jpj,jpk) , vslpml(jpi,jpj) , wslpj(jpi,jpj,jpk) , wslpjml(jpi,jpj) , STAT=ierr ) 818 789 IF( ierr > 0 ) CALL ctl_stop( 'STOP', 'ldf_slp_init : unable to allocate Madec operator slope ' ) 819 790 … … 825 796 wslpj(:,:,:) = 0._wp ; wslpjml(:,:) = 0._wp 826 797 827 IF( ln_traldf_hor .OR. ln_dynldf_hor ) THEN 828 IF(lwp) WRITE(numout,*) ' Horizontal mixing in s-coordinate: slope = slope of s-surfaces' 829 830 ! geopotential diffusion in s-coordinates on tracers and/or momentum 831 ! The slopes of s-surfaces are computed once (no call to ldfslp in step) 832 ! The slopes for momentum diffusion are i- or j- averaged of those on tracers 833 834 ! set the slope of diffusion to the slope of s-surfaces 835 ! ( c a u t i o n : minus sign as fsdep has positive value ) 836 DO jk = 1, jpk 837 DO jj = 2, jpjm1 838 DO ji = fs_2, fs_jpim1 ! vector opt. 839 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 843 END DO 844 END DO 845 END DO 846 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) ! Lateral boundary conditions 847 CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 848 ENDIF 798 !!gm I no longer understand this..... 799 !!gm IF( (ln_traldf_hor .OR. ln_dynldf_hor) .AND. .NOT. (.NOT.ln_linssh .AND. ln_rstart) ) THEN 800 ! IF(lwp) WRITE(numout,*) ' Horizontal mixing in s-coordinate: slope = slope of s-surfaces' 801 ! 802 ! ! geopotential diffusion in s-coordinates on tracers and/or momentum 803 ! ! The slopes of s-surfaces are computed once (no call to ldfslp in step) 804 ! ! The slopes for momentum diffusion are i- or j- averaged of those on tracers 805 ! 806 ! ! set the slope of diffusion to the slope of s-surfaces 807 ! ! ( c a u t i o n : minus sign as dep has positive value ) 808 ! DO jk = 1, jpk 809 ! DO jj = 2, jpjm1 810 ! DO ji = fs_2, fs_jpim1 ! vector opt. 811 ! uslp (ji,jj,jk) = - ( gdept_n(ji+1,jj,jk) - gdept_n(ji ,jj ,jk) ) * r1_e1u(ji,jj) * umask(ji,jj,jk) 812 ! vslp (ji,jj,jk) = - ( gdept_n(ji,jj+1,jk) - gdept_n(ji ,jj ,jk) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk) 813 ! wslpi(ji,jj,jk) = - ( gdepw_n(ji+1,jj,jk) - gdepw_n(ji-1,jj,jk) ) * r1_e1t(ji,jj) * wmask(ji,jj,jk) * 0.5 814 ! wslpj(ji,jj,jk) = - ( gdepw_n(ji,jj+1,jk) - gdepw_n(ji,jj-1,jk) ) * r1_e2t(ji,jj) * wmask(ji,jj,jk) * 0.5 815 ! END DO 816 ! END DO 817 ! END DO 818 ! CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) ! Lateral boundary conditions 819 ! CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 820 !!gm ENDIF 849 821 ENDIF 850 822 ! … … 852 824 ! 853 825 END SUBROUTINE ldf_slp_init 854 855 #else856 !!------------------------------------------------------------------------857 !! Dummy module : NO Rotation of lateral mixing tensor858 !!------------------------------------------------------------------------859 LOGICAL, PUBLIC, PARAMETER :: lk_ldfslp = .FALSE. !: slopes flag860 CONTAINS861 SUBROUTINE ldf_slp( kt, prd, pn2 ) ! Dummy routine862 INTEGER, INTENT(in) :: kt863 REAL, DIMENSION(:,:,:), INTENT(in) :: prd, pn2864 WRITE(*,*) 'ldf_slp: You should not have seen this print! error?', kt, prd(1,1,1), pn2(1,1,1)865 END SUBROUTINE ldf_slp866 SUBROUTINE ldf_slp_grif( kt ) ! Dummy routine867 INTEGER, INTENT(in) :: kt868 WRITE(*,*) 'ldf_slp_grif: You should not have seen this print! error?', kt869 END SUBROUTINE ldf_slp_grif870 SUBROUTINE ldf_slp_init ! Dummy routine871 END SUBROUTINE ldf_slp_init872 #endif873 826 874 827 !!======================================================================
Note: See TracChangeset
for help on using the changeset viewer.
