Changeset 1515
- Timestamp:
- 2009-07-21T16:35:10+02:00 (15 years ago)
- File:
-
- 1 edited
-
trunk/NEMO/OPA_SRC/LDF/ldfslp.F90 (modified) (30 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/LDF/ldfslp.F90
r1514 r1515 4 4 !! Ocean physics: slopes of neutral surfaces 5 5 !!====================================================================== 6 !! History : OPA ! 1994-12 (G. Madec, M. Imbard) Original code 7 !! 8.0 ! 1997-06 (G. Madec) optimization, lbc 8 !! 8.1 ! 1999-10 (A. Jouzeau) NEW profile in the mixed layer 9 !! NEMO 0.5 ! 2002-10 (G. Madec) Free form, F90 10 !! 1.0 ! 2005-10 (A. Beckmann) correction for s-coordinates 11 !!---------------------------------------------------------------------- 6 12 #if defined key_ldfslp || defined key_esopa 7 13 !!---------------------------------------------------------------------- … … 12 18 !! ldf_slp_init : initialization of the slopes computation 13 19 !!---------------------------------------------------------------------- 14 !! * Modules used15 20 USE oce ! ocean dynamics and tracers 16 21 USE dom_oce ! ocean space and time domain … … 26 31 PRIVATE 27 32 28 !! * Accessibility 29 PUBLIC ldf_slp ! routine called by step.F90 30 31 !! * Share module variables 32 LOGICAL , PUBLIC, PARAMETER :: lk_ldfslp = .TRUE. !: slopes flag 33 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: & !: 34 uslp, wslpi, & !: i_slope at U- and W-points 35 vslp, wslpj !: j-slope at V- and W-points 33 PUBLIC ldf_slp ! routine called by step.F90 34 35 LOGICAL , PUBLIC, PARAMETER :: lk_ldfslp = .TRUE. !: slopes flag 36 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: uslp, wslpi !: i_slope at U- and W-points 37 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: vslp, wslpj !: j-slope at V- and W-points 36 38 37 !! * Module variables 38 REAL(wp), DIMENSION(jpi,jpj,jpk) :: & 39 omlmask ! mask of the surface mixed layer at T-pt 40 REAL(wp), DIMENSION(jpi,jpj) :: & 41 uslpml, wslpiml, & ! i_slope at U- and W-points just below 42 ! ! the surface mixed layer 43 vslpml, wslpjml ! j_slope at V- and W-points just below 44 ! ! the surface mixed layer 39 REAL(wp), DIMENSION(jpi,jpj,jpk) :: omlmask ! mask of the surface mixed layer at T-pt 40 REAL(wp), DIMENSION(jpi,jpj) :: uslpml, wslpiml ! i_slope at U- and W-points just below the mixed layer 41 REAL(wp), DIMENSION(jpi,jpj) :: vslpml, wslpjml ! j_slope at V- and W-points just below the mixed layer 45 42 46 43 !! * Substitutions … … 48 45 # include "vectopt_loop_substitute.h90" 49 46 !!---------------------------------------------------------------------- 50 !! OPA 9.0 , LOCEAN-IPSL (2005)47 !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009) 51 48 !! $Id$ 52 49 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 59 56 !! *** ROUTINE ldf_slp *** 60 57 !! 61 !! ** Purpose : Compute the slopes of neutral surface (slope of iso -62 !! pycnalsurfaces referenced locally) ('key_traldfiso').63 !! 58 !! ** Purpose : Compute the slopes of neutral surface (slope of isopycnal 59 !! surfaces referenced locally) ('key_traldfiso'). 60 !! 64 61 !! ** Method : The slope in the i-direction is computed at U- and 65 62 !! W-points (uslp, wslpi) and the slope in the j-direction is … … 79 76 !! ** Action : - uslp, wslpi, and vslp, wslpj, the i- and j-slopes 80 77 !! of now neutral surfaces at u-, w- and v- w-points, resp. 81 !! 82 !! History : 83 !! 7.0 ! 94-12 (G. Madec, M. Imbard) Original code 84 !! 8.0 ! 97-06 (G. Madec) optimization, lbc 85 !! 8.1 ! 99-10 (A. Jouzeau) NEW profile 86 !! 8.5 ! 99-10 (G. Madec) Free form, F90 87 !! 9.0 ! 05-10 (A. Beckmann) correction for s-coordinates 88 !!---------------------------------------------------------------------- 89 !! * Modules used 90 USE oce , zgru => ua, & ! use ua as workspace 91 zgrv => va, & ! use va as workspace 92 zwy => ta, & ! use ta as workspace 93 zwz => sa ! use sa as workspace 94 !! * Arguments 95 INTEGER, INTENT( in ) :: kt ! ocean time-step index 96 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: & 97 prd, & ! in situ density 98 pn2 ! Brunt-Vaisala frequency (locally ref.) 99 100 !! * Local declarations 101 INTEGER :: ji, jj, jk ! dummy loop indices 102 INTEGER :: ii0, ii1, ij0, ij1, & ! temporary integer 103 & iku, ikv ! " " 104 REAL(wp) :: & 105 zeps, zmg, zm05g, & ! temporary scalars 106 zcoef1, zcoef2, zcoef3, & ! 107 zau, zbu, zav, zbv, & 108 zai, zbi, zaj, zbj, & 109 zcofu, zcofv, zcofw, & 110 z1u, z1v, z1wu, z1wv, & 111 zalpha 112 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zww 113 !!---------------------------------------------------------------------- 114 78 !!---------------------------------------------------------------------- 79 USE oce , zgru => ua ! use ua as workspace 80 USE oce , zgrv => va ! use va as workspace 81 USE oce , zwy => ta ! use ta as workspace 82 USE oce , zwz => sa ! use sa as workspace 83 !! 84 INTEGER , INTENT(in) :: kt ! ocean time-step index 85 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: prd ! in situ density 86 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pn2 ! Brunt-Vaisala frequency (locally ref.) 87 !! 88 INTEGER :: ji , jj , jk ! dummy loop indices 89 INTEGER :: ii0, ii1, iku ! temporary integer 90 INTEGER :: ij0, ij1, ikv ! temporary integer 91 REAL(wp) :: zeps, zmg, zm05g, zalpha ! temporary scalars 92 REAL(wp) :: zcoef1, zcoef2, zcoef3 ! - - 93 REAL(wp) :: zcofu , zcofv , zcofw ! - - 94 REAL(wp) :: zau, zbu, zai, zbi, z1u, z1wu ! - - 95 REAL(wp) :: zav, zbv, zaj, zbj, z1v, z1wv ! 96 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zww ! 3D workspace 97 !!---------------------------------------------------------------------- 115 98 116 ! 0. Initialization (first time-step only) 117 ! -------------- 118 119 IF( kt == nit000 ) CALL ldf_slp_init 120 121 122 ! 0. Local constant initialization 123 ! -------------------------------- 124 125 zeps = 1.e-20 99 IF( kt == nit000 ) CALL ldf_slp_init ! initialization (first time-step only) 100 101 zeps = 1.e-20 ! Local constant initialization 126 102 zmg = -1.0 / grav 127 103 zm05g = -0.5 / grav 128 104 ! 129 105 zww(:,:,:) = 0.e0 130 106 zwz(:,:,:) = 0.e0 131 132 ! horizontal density gradient computation 107 ! ! horizontal density gradient computation 133 108 DO jk = 1, jpk 134 109 DO jj = 1, jpjm1 … … 139 114 END DO 140 115 END DO 141 142 IF( ln_zps ) THEN ! partial steps correction at the bottom ocean level (zps_hde routine) 116 IF( ln_zps ) THEN ! partial steps correction at the bottom ocean level 143 117 # if defined key_vectopt_loop 144 jj =1145 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling)118 DO jj = 1, 1 119 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 146 120 # else 147 121 DO jj = 1, jpjm1 148 122 DO ji = 1, jpim1 149 123 # endif 150 ! last ocean level 151 iku = MIN ( mbathy(ji,jj), mbathy(ji+1,jj) ) - 1 124 iku = MIN ( mbathy(ji,jj), mbathy(ji+1,jj) ) - 1 ! last ocean level 152 125 ikv = MIN ( mbathy(ji,jj), mbathy(ji,jj+1) ) - 1 153 126 zgru(ji,jj,iku) = gru(ji,jj) 154 127 zgrv(ji,jj,ikv) = grv(ji,jj) 155 # if ! defined key_vectopt_loop 156 END DO 157 # endif 128 END DO 158 129 END DO 159 130 ENDIF 160 161 ! Slopes of isopycnal surfaces just below the mixed layer162 ! -------------------------------------------------------163 131 164 CALL ldf_slp_mxl( prd, pn2 ) 132 CALL ldf_slp_mxl( prd, pn2 ) ! Slopes of isopycnal surfaces just below the mixed layer 133 165 134 166 !-------------------synchro--------------------------------------------- 167 168 ! ! =============== 169 DO jk = 2, jpkm1 ! Horizontal slab 170 ! ! =============== 171 172 ! I. slopes at u and v point 173 ! =========================== 174 175 176 ! I.1. Slopes of isopycnal surfaces 177 ! --------------------------------- 178 ! uslp = d/di( prd ) / d/dz( prd ) 179 ! vslp = d/dj( prd ) / d/dz( prd ) 180 181 ! Local vertical density gradient evaluated from N^2 182 ! zwy = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point 183 135 ! I. slopes at u and v point | uslp = d/di( prd ) / d/dz( prd ) 136 ! =========================== | vslp = d/dj( prd ) / d/dz( prd ) 137 ! 138 ! !* Local vertical density gradient evaluated from N^2 139 DO jk = 2, jpkm1 ! zwy = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point 184 140 DO jj = 1, jpj 185 141 DO ji = 1, jpi 186 zwy(ji,jj,jk) = zmg * ( prd(ji,jj,jk) + 1. ) & 187 & * ( pn2(ji,jj,jk) + pn2(ji,jj,jk+1) ) & 188 & / MAX( tmask(ji,jj,jk) + tmask (ji,jj,jk+1), 1. ) 189 END DO 190 END DO 191 192 ! Slope at u and v points 142 zwy(ji,jj,jk) = zmg * ( prd(ji,jj,jk) + 1. ) * ( pn2 (ji,jj,jk) + pn2 (ji,jj,jk+1) ) & 143 & / MAX( tmask(ji,jj,jk) + tmask(ji,jj,jk+1), 1. ) 144 END DO 145 END DO 146 END DO 147 DO jk = 2, jpkm1 !* Slopes at u and v points 193 148 DO jj = 2, jpjm1 194 149 DO ji = fs_2, fs_jpim1 ! vector opt. … … 204 159 ! uslp and vslp output in zwz and zww, resp. 205 160 zalpha = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) ) 206 zwz (ji,jj,jk) = ( zau / ( zbu - zeps ) * ( 1. - zalpha)&207 & + zalpha * uslpml(ji,jj)&208 & * 0.5 * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk)-fse3u(ji,jj,1) ) &209 & / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 5. ) ) * umask(ji,jj,jk)161 zwz (ji,jj,jk) = ( ( 1. - zalpha) * zau / ( zbu - zeps ) & 162 & + zalpha * uslpml(ji,jj) & 163 & * 0.5 * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk)-fse3u(ji,jj,1) ) & 164 & / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 5. ) ) * umask(ji,jj,jk) 210 165 zalpha = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) ) 211 zww (ji,jj,jk) = ( zav / ( zbv - zeps ) * ( 1. - zalpha) & 212 & + zalpha * vslpml(ji,jj) & 213 & * 0.5 * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk)-fse3v(ji,jj,1) ) & 214 & / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 5. ) ) * vmask(ji,jj,jk) 215 END DO 216 END DO 217 ! ! =============== 218 END DO ! end of slab 219 ! ! =============== 220 221 222 ! lateral boundary conditions on zww and zwz 223 CALL lbc_lnk( zwz, 'U', -1. ) 224 CALL lbc_lnk( zww, 'V', -1. ) 225 226 ! ! =============== 227 DO jk = 2, jpkm1 ! Horizontal slab 228 ! ! =============== 229 230 ! Shapiro filter applied in the horizontal direction 231 zcofu = 1. / 16. 232 zcofv = 1. / 16. 233 DO jj = 2, jpjm1, jpj-3 ! row jj=2 and =jpjm1 only 166 zww (ji,jj,jk) = ( ( 1. - zalpha) * zav / ( zbv - zeps ) & 167 & + zalpha * vslpml(ji,jj) & 168 & * 0.5 * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk)-fse3v(ji,jj,1) ) & 169 & / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 5. ) ) * vmask(ji,jj,jk) 170 END DO 171 END DO 172 END DO 173 CALL lbc_lnk( zwz, 'U', -1. ) ; CALL lbc_lnk( zww, 'V', -1. ) ! lateral boundary conditions 174 ! 175 zcofu = 1. / 16. !* horizontal Shapiro filter 176 zcofv = 1. / 16. 177 DO jk = 2, jpkm1 178 DO jj = 2, jpjm1, jpj-3 ! rows jj=2 and =jpjm1 only 234 179 DO ji = 2, jpim1 235 !uslop236 uslp(ji,jj,jk) = zcofu * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) &237 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) &238 & + 2.*(zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) &239 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) &240 & + 4.* zwz(ji ,jj ,jk) )241 ! vslop242 vslp(ji,jj,jk) = zcofv * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) &243 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) &244 & + 2.*(zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) &245 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) &246 & + 4.* zww(ji,jj ,jk) )247 END DO248 END DO249 250 DO jj = 3, jpj-2251 DO ji = fs_2, fs_jpim1 ! vector opt.252 ! uslop253 180 uslp(ji,jj,jk) = zcofu * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 254 181 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & … … 256 183 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 257 184 & + 4.* zwz(ji ,jj ,jk) ) 258 ! vslop259 185 vslp(ji,jj,jk) = zcofv * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 260 186 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & … … 264 190 END DO 265 191 END DO 266 267 ! decrease along coastal boundaries 192 DO jj = 3, jpj-2 ! other rows 193 DO ji = fs_2, fs_jpim1 ! vector opt. 194 uslp(ji,jj,jk) = zcofu * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 195 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 196 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 197 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 198 & + 4.* zwz(ji ,jj ,jk) ) 199 vslp(ji,jj,jk) = zcofv * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 200 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 201 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 202 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 203 & + 4.* zww(ji,jj ,jk) ) 204 END DO 205 END DO 206 ! !* decrease along coastal boundaries 268 207 DO jj = 2, jpjm1 269 208 DO ji = fs_2, fs_jpim1 ! vector opt. … … 276 215 END DO 277 216 END DO 278 279 280 ! II. Computation of slopes at w point 281 ! ==================================== 282 283 284 ! II.1 Slopes of isopycnal surfaces 285 ! --------------------------------- 286 ! wslpi = mij( d/di( prd ) / d/dz( prd ) 287 ! wslpj = mij( d/dj( prd ) / d/dz( prd ) 288 289 290 ! Local vertical density gradient evaluated from N^2 291 ! zwy = d/dz(prd)= - mk ( prd ) / grav * pn2 -- at w point 217 END DO 218 219 220 ! II. slopes at w point | wslpi = mij( d/di( prd ) / d/dz( prd ) 221 ! =========================== | wslpj = mij( d/dj( prd ) / d/dz( prd ) 222 ! 223 ! !* Local vertical density gradient evaluated from N^2 224 DO jk = 2, jpkm1 ! zwy = d/dz(prd)= - mk ( prd ) / grav * pn2 -- at w point 292 225 DO jj = 1, jpj 293 226 DO ji = 1, jpi 294 zwy (ji,jj,jk) = zm05g * pn2 (ji,jj,jk) * & 295 & ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) 296 END DO 297 END DO 298 299 ! Slope at w point 227 zwy(ji,jj,jk) = zm05g * pn2 (ji,jj,jk) * ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) 228 END DO 229 END DO 230 END DO 231 DO jk = 2, jpkm1 !* Slopes at w point 300 232 DO jj = 2, jpjm1 301 233 DO ji = fs_2, fs_jpim1 ! vector opt. 302 ! horizontal density i-gradient at w-points234 ! ! horizontal density i-gradient at w-points 303 235 zcoef1 = MAX( zeps, umask(ji-1,jj,jk )+umask(ji,jj,jk ) & 304 236 & +umask(ji-1,jj,jk-1)+umask(ji,jj,jk-1) ) … … 306 238 zai = zcoef1 * ( zgru(ji ,jj,jk ) + zgru(ji ,jj,jk-1) & 307 239 & + zgru(ji-1,jj,jk-1) + zgru(ji-1,jj,jk ) ) * tmask (ji,jj,jk) 308 ! horizontal density j-gradient at w-points240 ! ! horizontal density j-gradient at w-points 309 241 zcoef2 = MAX( zeps, vmask(ji,jj-1,jk )+vmask(ji,jj,jk-1) & 310 242 & +vmask(ji,jj-1,jk-1)+vmask(ji,jj,jk ) ) … … 312 244 zaj = zcoef2 * ( zgrv(ji,jj ,jk ) + zgrv(ji,jj ,jk-1) & 313 245 & + zgrv(ji,jj-1,jk-1) + zgrv(ji,jj-1,jk ) ) * tmask (ji,jj,jk) 314 ! bound the slopes: abs(zw.)<= 1/100 and zb..<0. 315 ! static instability: 316 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 246 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0. 247 ! ! static instability: kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 317 248 zbi = MIN( zwy (ji,jj,jk),- 100.*ABS(zai), -7.e+3/fse3w(ji,jj,jk)*ABS(zai) ) 318 249 zbj = MIN( zwy (ji,jj,jk), -100.*ABS(zaj), -7.e+3/fse3w(ji,jj,jk)*ABS(zaj) ) 319 ! wslpi and wslpj output in zwz and zww, resp.250 ! ! wslpi and wslpj output in zwz and zww, resp. 320 251 zalpha = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) ) 321 252 zcoef3 = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10. ) … … 326 257 END DO 327 258 END DO 328 ! ! =============== 329 END DO ! end of slab 330 ! ! =============== 331 332 333 ! lateral boundary conditions on zwz and zww 334 CALL lbc_lnk( zwz, 'T', -1. ) 335 CALL lbc_lnk( zww, 'T', -1. ) 336 337 ! ! =============== 338 DO jk = 2, jpkm1 ! Horizontal slab 339 ! ! =============== 340 341 ! Shapiro filter applied in the horizontal direction 342 343 DO jj = 2, jpjm1, jpj-3 ! row jj=2 and =jpjm1 259 END DO 260 CALL lbc_lnk( zwz, 'T', -1. ) ; CALL lbc_lnk( zww, 'T', -1. ) ! lateral boundary conditions on zwz and zww 261 ! 262 ! !* horizontal Shapiro filter 263 DO jk = 2, jpkm1 264 DO jj = 2, jpjm1, jpj-3 ! rows jj=2 and =jpjm1 344 265 DO ji = 2, jpim1 345 zcofw = tmask(ji,jj,jk) /16.266 zcofw = tmask(ji,jj,jk) / 16. 346 267 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 347 268 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & … … 356 277 & + 4.* zww(ji ,jj ,jk) ) * zcofw 357 278 END DO 358 END DO 359 360 DO jj = 3, jpj-2 279 END DO 280 DO jj = 3, jpj-2 ! other rows 361 281 DO ji = fs_2, fs_jpim1 ! vector opt. 362 zcofw = tmask(ji,jj,jk) /16.282 zcofw = tmask(ji,jj,jk) / 16. 363 283 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 364 284 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & … … 374 294 END DO 375 295 END DO 376 377 ! decrease the slope along the boundaries 296 ! !* decrease along coastal boundaries 378 297 DO jj = 2, jpjm1 379 298 DO ji = fs_2, fs_jpim1 ! vector opt. … … 384 303 END DO 385 304 END DO 305 END DO 386 306 387 388 ! III. Specific grid points 389 ! ------------------------- 390 391 IF( cp_cfg == "orca" .AND. jp_cfg == 4 ) THEN 392 ! ! ======================= 393 ! Horizontal diffusion in ! ORCA_R4 configuration 394 ! specific area ! ======================= 395 ! 396 ! ! Gibraltar Strait 397 ij0 = 50 ; ij1 = 53 398 ii0 = 69 ; ii1 = 71 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk ) = 0.e0 399 ij0 = 51 ; ij1 = 53 400 ii0 = 68 ; ii1 = 71 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk ) = 0.e0 401 ii0 = 69 ; ii1 = 71 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk ) = 0.e0 402 ii0 = 69 ; ii1 = 71 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk ) = 0.e0 403 404 ! ! Mediterrannean Sea 405 ij0 = 49 ; ij1 = 56 406 ii0 = 71 ; ii1 = 90 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk ) = 0.e0 407 ij0 = 50 ; ij1 = 56 408 ii0 = 70 ; ii1 = 90 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk ) = 0.e0 409 ii0 = 71 ; ii1 = 90 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk ) = 0.e0 410 ii0 = 71 ; ii1 = 90 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk ) = 0.e0 411 ENDIF 412 ! ! =============== 413 END DO ! end of slab 414 ! ! =============== 307 308 ! III. Specific grid points 309 ! =========================== 310 ! 311 IF( cp_cfg == "orca" .AND. jp_cfg == 4 ) THEN ! ORCA_R4 configuration: horizontal diffusion in specific area 312 ! ! Gibraltar Strait 313 ij0 = 50 ; ij1 = 53 314 ii0 = 69 ; ii1 = 71 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0.e0 315 ij0 = 51 ; ij1 = 53 316 ii0 = 68 ; ii1 = 71 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0.e0 317 ii0 = 69 ; ii1 = 71 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0.e0 318 ii0 = 69 ; ii1 = 71 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0.e0 319 ! 320 ! ! Mediterrannean Sea 321 ij0 = 49 ; ij1 = 56 322 ii0 = 71 ; ii1 = 90 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0.e0 323 ij0 = 50 ; ij1 = 56 324 ii0 = 70 ; ii1 = 90 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0.e0 325 ii0 = 71 ; ii1 = 90 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0.e0 326 ii0 = 71 ; ii1 = 90 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0.e0 327 ENDIF 415 328 416 329 417 ! I II Lateral boundary conditions on all slopes (uslp , vslp,418 ! ------------------------------- wslpi, wslpj )330 ! IV. Lateral boundary conditions 331 ! =============================== 419 332 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) 420 333 CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 334 421 335 422 336 IF(ln_ctl) THEN … … 424 338 CALL prt_ctl(tab3d_1=wslpi, clinfo1=' slp - wi: ', tab3d_2=wslpj, clinfo2=' wj: ', kdim=jpk) 425 339 ENDIF 426 340 ! 427 341 END SUBROUTINE ldf_slp 428 342 … … 431 345 !!---------------------------------------------------------------------- 432 346 !! *** ROUTINE ldf_slp_mxl *** 433 !! ** Purpose :434 !! Compute the slopes of iso-neutral surface (slope of isopycnal435 !! surfaces referenced locally) just abovethe mixed layer.347 !! 348 !! ** Purpose : Compute the slopes of iso-neutral surface just below 349 !! the mixed layer. 436 350 !! 437 351 !! ** Method : … … 445 359 !! of 10cm/s) 446 360 !! 447 !! ** Action : 448 !! Compute uslp, wslpi, and vslp, wslpj, the i- and j-slopes 361 !! ** Action : Compute uslp, wslpi, and vslp, wslpj, the i- and j-slopes 449 362 !! of now neutral surfaces at u-, w- and v- w-points, resp. 450 !! 451 !! History : 452 !! 8.1 ! 99-10 (A. Jouzeau) Original code 453 !! 8.5 ! 99-10 (G. Madec) Free form, F90 454 !!---------------------------------------------------------------------- 455 !! * Modules used 456 USE oce , zgru => ua, & ! ua, va used as workspace and set to hor. 457 zgrv => va ! density gradient in ldf_slp 458 459 !! * Arguments 460 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: & 461 prd, & ! in situ density 462 pn2 ! Brunt-Vaisala frequency (locally ref.) 463 464 !! * Local declarations 465 INTEGER :: ji, jj, jk ! dummy loop indices 466 INTEGER :: ik, ikm1 ! temporary integers 467 REAL(wp), DIMENSION(jpi,jpj) :: & 468 zwy ! temporary workspace 469 REAL(wp) :: & 470 zeps, zmg, zm05g, & ! temporary scalars 471 zcoef1, zcoef2, & ! " " 472 zau, zbu, zav, zbv, & ! " " 473 zai, zbi, zaj, zbj ! " " 474 !!---------------------------------------------------------------------- 475 476 477 ! 0. Local constant initialization 478 ! -------------------------------- 479 480 zeps = 1.e-20 363 !!---------------------------------------------------------------------- 364 USE oce , zgru => ua ! ua, va used as workspace and set to hor. 365 USE oce , zgrv => va ! density gradient in ldf_slp 366 !! 367 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: prd ! in situ density 368 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: pn2 ! Brunt-Vaisala frequency (locally ref.) 369 !! 370 INTEGER :: ji, jj, jk ! dummy loop indices 371 INTEGER :: ik, ikm1 ! temporary integers 372 REAL(wp) :: zeps, zmg, zm05g ! temporary scalars 373 REAL(wp) :: zcoef1, zcoef2 ! - - 374 REAL(wp) :: zau, zbu, zai, zbi ! - - 375 REAL(wp) :: zav, zbv, zaj, zbj ! - - 376 REAL(wp), DIMENSION(jpi,jpj) :: zwy ! 2D workspace 377 !!---------------------------------------------------------------------- 378 379 zeps = 1.e-20 ! Local constant initialization 481 380 zmg = -1.0 / grav 482 381 zm05g = -0.5 / grav 483 484 382 ! 485 383 uslpml (1,:) = 0.e0 ; uslpml (jpi,:) = 0.e0 486 384 vslpml (1,:) = 0.e0 ; vslpml (jpi,:) = 0.e0 … … 488 386 wslpjml(1,:) = 0.e0 ; wslpjml(jpi,:) = 0.e0 489 387 490 ! surface mixed layer mask 491 492 ! mask for mixed layer 493 DO jk = 1, jpk 388 ! ! surface mixed layer mask 389 DO jk = 1, jpk ! =1 inside the mixed layer, =0 otherwise 494 390 # if defined key_vectopt_loop 495 jj =1496 DO ji = 1, jpij ! vector opt. (forced unrolling)391 DO jj = 1, 1 392 DO ji = 1, jpij ! vector opt. (forced unrolling) 497 393 # else 498 394 DO jj = 1, jpj 499 395 DO ji = 1, jpi 500 396 # endif 501 ! mixed layer interior (mask = 1) and exterior (mask = 0)502 397 ik = nmln(ji,jj) - 1 503 IF( jk <= ik ) THEN 504 omlmask(ji,jj,jk) = 1.e0 505 ELSE 506 omlmask(ji,jj,jk) = 0.e0 398 IF( jk <= ik ) THEN ; omlmask(ji,jj,jk) = 1.e0 399 ELSE ; omlmask(ji,jj,jk) = 0.e0 507 400 ENDIF 508 # if ! defined key_vectopt_loop 509 END DO 510 # endif 401 END DO 511 402 END DO 512 403 END DO … … 515 406 ! Slopes of isopycnal surfaces just before bottom of mixed layer 516 407 ! -------------------------------------------------------------- 517 ! uslpml = d/di( prd ) / d/dz( prd )518 ! vslpml = d/dj( prd ) / d/dz( prd )519 520 ! Local vertical density gradient evaluated from N^2521 ! zwy = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point522 408 ! The slope are computed as in the 3D case. 409 ! A key point here is the definition of the mixed layer at u- and v-points. 410 ! It is assumed to be the maximum of the two neighbouring T-point mixed layer depth. 411 ! Otherwise, a n2 value inside the mixed layer can be involved in the computation 412 ! of the slope, resulting in a too steep diagnosed slope and thus a spurious eddy 413 ! induce velocity field near the base of the mixed layer. 523 414 !----------------------------------------------------------------------- 524 zwy(:,jpj) = 0.e0 415 ! 416 zwy(:,jpj) = 0.e0 !* vertical density gradient for u-slope (from N^2) 525 417 zwy(jpi,:) = 0.e0 526 418 # if defined key_vectopt_loop 527 jj =1528 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling)419 DO jj = 1, 1 420 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 529 421 # else 530 422 DO jj = 1, jpjm1 531 423 DO ji = 1, jpim1 532 424 # endif 533 ik = MAX( 1, nmln(ji,jj) , nmln(ji+1,jj) ) 534 ! if ik = jpk take jpkm1 values 535 ik = MIN( ik,jpkm1 ) 536 zwy(ji,jj) = zmg * ( prd(ji,jj,ik) + 1. ) & 537 & * ( pn2(ji,jj,ik) + pn2(ji,jj,ik+1) ) & 538 & / MAX( tmask(ji,jj,ik) + tmask (ji,jj,ik+1), 1. ) 539 # if ! defined key_vectopt_loop 540 END DO 541 # endif 542 END DO 543 CALL lbc_lnk( zwy, 'U', 1. ) ! lateral boundary conditions (NO sign change) 544 545 546 ! Slope at u points 425 ik = MAX( 1, nmln(ji,jj) , nmln(ji+1,jj) ) ! avoid spurious recirculation 426 ik = MIN( ik, jpkm1 ) ! if ik = jpk take jpkm1 values 427 zwy(ji,jj) = zmg * ( prd(ji,jj,ik) + 1. ) * ( pn2 (ji,jj,ik) + pn2 (ji,jj,ik+1) ) & 428 & / MAX( tmask(ji,jj,ik) + tmask(ji,jj,ik+1), 1. ) 429 END DO 430 END DO 431 CALL lbc_lnk( zwy, 'U', 1. ) ! lateral boundary conditions NO sign change 432 433 ! !* Slope at u points 547 434 # if defined key_vectopt_loop 548 jj =1549 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling)435 DO jj = 1, 1 436 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) 550 437 # else 551 438 DO jj = 2, jpjm1 … … 554 441 ! horizontal and vertical density gradient at u-points 555 442 ik = MAX( 1, nmln(ji,jj) , nmln(ji+1,jj) ) 556 ik = MIN( ik, jpkm1 )443 ik = MIN( ik, jpkm1 ) 557 444 zau = 1./ e1u(ji,jj) * zgru(ji,jj,ik) 558 445 zbu = 0.5*( zwy(ji,jj) + zwy(ji+1,jj) ) … … 562 449 ! uslpml 563 450 uslpml (ji,jj) = zau / ( zbu - zeps ) * umask (ji,jj,ik) 564 # if ! defined key_vectopt_loop 565 END DO 566 # endif 567 END DO 568 569 ! lateral boundary conditions on uslpml 570 CALL lbc_lnk( uslpml, 'U', -1. ) 571 572 ! Local vertical density gradient evaluated from N^2 573 ! zwy = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point 574 zwy ( :, jpj) = 0.e0 575 zwy ( jpi, :) = 0.e0 451 END DO 452 END DO 453 CALL lbc_lnk( uslpml, 'U', -1. ) ! lateral boundary conditions (i-gradient => sign change) 454 455 ! !* vertical density gradient for v-slope (from N^2) 576 456 # if defined key_vectopt_loop 577 jj =1578 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling)457 DO jj = 1, 1 458 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 579 459 # else 580 460 DO jj = 1, jpjm1 … … 582 462 # endif 583 463 ik = MAX( 1, nmln(ji,jj) , nmln(ji,jj+1) ) 584 ik = MIN( ik,jpkm1 ) 585 zwy(ji,jj) = zmg * ( prd(ji,jj,ik) + 1. ) & 586 & * ( pn2(ji,jj,ik) + pn2(ji,jj,ik+1) ) & 587 & / MAX( tmask(ji,jj,ik) + tmask (ji,jj,ik+1), 1. ) 588 # if ! defined key_vectopt_loop 589 END DO 590 # endif 591 END DO 592 CALL lbc_lnk( zwy, 'V', 1. ) ! lateral boundary conditions No change in sign 593 594 ! Slope at v points 464 ik = MIN( ik, jpkm1 ) 465 zwy(ji,jj) = zmg * ( prd(ji,jj,ik) + 1. ) * ( pn2 (ji,jj,ik) + pn2 (ji,jj,ik+1) ) & 466 & / MAX( tmask(ji,jj,ik) + tmask(ji,jj,ik+1), 1. ) 467 END DO 468 END DO 469 CALL lbc_lnk( zwy, 'V', 1. ) ! lateral boundary conditions NO sign change 470 471 ! !* Slope at v points 595 472 # if defined key_vectopt_loop 596 jj =1597 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling)473 DO jj = 1, 1 474 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) 598 475 # else 599 476 DO jj = 2, jpjm1 … … 610 487 ! vslpml 611 488 vslpml (ji,jj) = zav / ( zbv - zeps ) * vmask (ji,jj,ik) 612 # if ! defined key_vectopt_loop 613 END DO 614 # endif 615 END DO 616 617 ! lateral boundary conditions on vslpml 618 CALL lbc_lnk( vslpml, 'V', -1. ) 619 620 ! wslpiml = mij( d/di( prd ) / d/dz( prd ) 621 ! wslpjml = mij( d/dj( prd ) / d/dz( prd ) 622 623 624 ! Local vertical density gradient evaluated from N^2 625 ! zwy = d/dz(prd)= - mk ( prd ) / grav * pn2 -- at w point 489 END DO 490 END DO 491 CALL lbc_lnk( vslpml, 'V', -1. ) ! lateral boundary conditions (j-gradient => sign change) 492 493 494 ! !* vertical density gradient for w-slope (from N^2) 626 495 # if defined key_vectopt_loop 627 jj =1628 DO ji = 1, jpij ! vector opt. (forced unrolling)496 DO jj = 1, 1 497 DO ji = 1, jpij ! vector opt. (forced unrolling) 629 498 # else 630 499 DO jj = 1, jpj 631 500 DO ji = 1, jpi 632 501 # endif 633 ik = nmln(ji,jj) +1634 ik = MIN( ik, jpk )502 ik = nmln(ji,jj) + 1 503 ik = MIN( ik, jpk ) 635 504 ikm1 = MAX ( 1, ik-1) 636 505 zwy (ji,jj) = zm05g * pn2 (ji,jj,ik) * & 637 506 & ( prd (ji,jj,ik) + prd (ji,jj,ikm1) + 2. ) 638 # if ! defined key_vectopt_loop 639 END DO 640 # endif 641 END DO 642 643 ! Slope at w point 507 END DO 508 END DO 509 510 ! !* Slopes at w points 644 511 # if defined key_vectopt_loop 645 jj =1646 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling)512 DO jj = 1, 1 513 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) 647 514 # else 648 515 DO jj = 2, jpjm1 649 516 DO ji = 2, jpim1 650 517 # endif 651 ik = nmln(ji,jj) +1652 ik = MIN( ik, jpk )653 ikm1 = MAX ( 1, ik-1 )518 ik = nmln(ji,jj) + 1 519 ik = MIN( ik, jpk ) 520 ikm1 = MAX ( 1, ik-1 ) 654 521 ! horizontal density i-gradient at w-points 655 522 zcoef1 = MAX( zeps, umask(ji-1,jj,ik )+umask(ji,jj,ik ) & … … 672 539 wslpiml (ji,jj) = zai / ( zbi - zeps) * tmask (ji,jj,ik) 673 540 wslpjml (ji,jj) = zaj / ( zbj - zeps) * tmask (ji,jj,ik) 674 # if ! defined key_vectopt_loop 675 END DO 676 # endif 677 END DO 678 679 ! lateral boundary conditions on wslpiml and wslpjml 680 CALL lbc_lnk( wslpiml, 'W', -1. ) 681 CALL lbc_lnk( wslpjml, 'W', -1. ) 682 541 END DO 542 END DO 543 CALL lbc_lnk( wslpiml, 'W', -1. ) ; CALL lbc_lnk( wslpjml, 'W', -1. ) ! lateral boundary conditions 544 ! 683 545 END SUBROUTINE ldf_slp_mxl 684 546 … … 692 554 !! ** Method : read the nammbf namelist and check the parameter 693 555 !! values called by tra_dmp at the first timestep (nit000) 694 !! 695 !! History : 696 !! 8.5 ! 02-06 (G. Madec) original code 697 !!---------------------------------------------------------------------- 698 !! * local declarations 556 !!---------------------------------------------------------------------- 699 557 INTEGER :: ji, jj, jk ! dummy loop indices 700 558 !!---------------------------------------------------------------------- 701 559 702 703 ! Parameter control and print 704 ! --------------------------- 705 IF(lwp) THEN 560 IF(lwp) THEN 706 561 WRITE(numout,*) 707 562 WRITE(numout,*) 'ldf_slp : direction of lateral mixing' … … 744 599 END DO 745 600 END DO 746 747 601 ! Lateral boundary conditions on the slopes 748 602 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) 749 603 CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 750 604 ENDIF 751 605 ! 752 606 END SUBROUTINE ldf_slp_init 753 607 … … 760 614 SUBROUTINE ldf_slp( kt, prd, pn2 ) ! Dummy routine 761 615 INTEGER, INTENT(in) :: kt 762 REAL, DIMENSION(:,:,:), INTENT(in) :: prd, pn2616 REAL, DIMENSION(:,:,:), INTENT(in) :: prd, pn2 763 617 WRITE(*,*) 'ldf_slp: You should not have seen this print! error?', kt, prd(1,1,1), pn2(1,1,1) 764 618 END SUBROUTINE ldf_slp
Note: See TracChangeset
for help on using the changeset viewer.
