Changeset 2528 for trunk/NEMOGCM/NEMO/OPA_SRC/LDF
- Timestamp:
- 2010-12-27T18:33:53+01:00 (14 years ago)
- Location:
- trunk/NEMOGCM/NEMO/OPA_SRC/LDF
- Files:
-
- 17 edited
-
ldfdyn.F90 (modified) (4 diffs, 2 props)
-
ldfdyn_antarctic.h90 (modified) (1 diff, 1 prop)
-
ldfdyn_arctic.h90 (modified) (1 diff, 1 prop)
-
ldfdyn_c1d.h90 (modified) (6 diffs, 1 prop)
-
ldfdyn_c2d.h90 (modified) (5 diffs, 2 props)
-
ldfdyn_c3d.h90 (modified) (11 diffs, 2 props)
-
ldfdyn_oce.F90 (modified) (2 diffs, 1 prop)
-
ldfdyn_substitute.h90 (modified) (2 diffs, 1 prop)
-
ldfeiv.F90 (modified) (10 diffs, 1 prop)
-
ldfeiv_substitute.h90 (modified) (1 diff, 1 prop)
-
ldfslp.F90 (modified) (24 diffs, 1 prop)
-
ldftra.F90 (modified) (3 diffs, 2 props)
-
ldftra_c1d.h90 (modified) (1 diff, 1 prop)
-
ldftra_c2d.h90 (modified) (1 diff, 1 prop)
-
ldftra_c3d.h90 (modified) (1 diff, 1 prop)
-
ldftra_oce.F90 (modified) (3 diffs, 1 prop)
-
ldftra_substitute.h90 (modified) (1 diff, 1 prop)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn.F90
- Property svn:eol-style deleted
- Property svn:executable deleted
r1601 r2528 36 36 # include "domzgr_substitute.h90" 37 37 !!---------------------------------------------------------------------- 38 !! NEMO/OPA 3. 2 , LOCEAN-IPSL (2009)38 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 39 39 !! $Id$ 40 !! Software governed by the CeCILL licence ( modipsl/doc/NEMO_CeCILL.txt)40 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 41 41 !!---------------------------------------------------------------------- 42 42 … … 67 67 NAMELIST/namdyn_ldf/ ln_dynldf_lap , ln_dynldf_bilap, & 68 68 & ln_dynldf_level, ln_dynldf_hor , ln_dynldf_iso, & 69 & rn_ahm_0 , rn_ahmb_069 & rn_ahm_0_lap , rn_ahmb_0 , rn_ahm_0_blp 70 70 !!---------------------------------------------------------------------- 71 71 … … 78 78 WRITE(numout,*) '~~~~~~~' 79 79 WRITE(numout,*) ' Namelist nam_dynldf : set lateral mixing parameters' 80 WRITE(numout,*) ' laplacian operator ln_dynldf_lap = ', ln_dynldf_lap 81 WRITE(numout,*) ' bilaplacian operator ln_dynldf_bilap = ', ln_dynldf_bilap 82 WRITE(numout,*) ' iso-level ln_dynldf_level = ', ln_dynldf_level 83 WRITE(numout,*) ' horizontal (geopotential) ln_dynldf_hor = ', ln_dynldf_hor 84 WRITE(numout,*) ' iso-neutral ln_dynldf_iso = ', ln_dynldf_iso 85 WRITE(numout,*) ' horizontal eddy viscosity rn_ahm_0 = ', rn_ahm_0 86 WRITE(numout,*) ' background viscosity rn_ahmb_0 = ', rn_ahmb_0 87 ENDIF 88 89 ahm0 = rn_ahm_0 ! OLD namelist variables defined from DOCTOR namelist variables 90 ahmb0 = rn_ahmb_0 80 WRITE(numout,*) ' laplacian operator ln_dynldf_lap = ', ln_dynldf_lap 81 WRITE(numout,*) ' bilaplacian operator ln_dynldf_bilap = ', ln_dynldf_bilap 82 WRITE(numout,*) ' iso-level ln_dynldf_level = ', ln_dynldf_level 83 WRITE(numout,*) ' horizontal (geopotential) ln_dynldf_hor = ', ln_dynldf_hor 84 WRITE(numout,*) ' iso-neutral ln_dynldf_iso = ', ln_dynldf_iso 85 WRITE(numout,*) ' horizontal laplacian eddy viscosity rn_ahm_0_lap = ', rn_ahm_0_lap 86 WRITE(numout,*) ' background viscosity rn_ahmb_0 = ', rn_ahmb_0 87 WRITE(numout,*) ' horizontal bilaplacian eddy viscosity rn_ahm_0_blp = ', rn_ahm_0_blp 88 ENDIF 89 90 ahm0 = rn_ahm_0_lap ! OLD namelist variables defined from DOCTOR namelist variables 91 ahmb0 = rn_ahmb_0 92 ahm0_blp = rn_ahm_0_blp 91 93 92 94 ! ... check of lateral diffusive operator on tracers … … 117 119 IF( ln_dynldf_bilap ) THEN 118 120 IF(lwp) WRITE(numout,*) ' biharmonic momentum diffusion' 119 IF( ahm0 > 0 .AND. .NOT. lk_esopa ) CALL ctl_stop( 'The horizontal viscosity coef. ahm0 must be negative' ) 121 IF( .NOT. ln_dynldf_lap ) ahm0 = ahm0_blp ! Allow spatially varying coefs, which use ahm0 as input 122 IF( ahm0_blp > 0 .AND. .NOT. lk_esopa ) CALL ctl_stop( 'The horizontal viscosity coef. ahm0 must be negative' ) 120 123 ELSE 121 124 IF(lwp) WRITE(numout,*) ' harmonic momentum diff. (default)' -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_antarctic.h90
- Property svn:eol-style deleted
r1152 r2528 3 3 !!---------------------------------------------------------------------- 4 4 !!---------------------------------------------------------------------- 5 !! OPA 9.0 , LOCEAN-IPSL (2005)5 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 6 6 !! $Id$ 7 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt7 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 8 8 !!---------------------------------------------------------------------- 9 9 ! -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_arctic.h90
- Property svn:eol-style deleted
r1152 r2528 3 3 !!---------------------------------------------------------------------- 4 4 !!---------------------------------------------------------------------- 5 !! OPA 9.0 , LOCEAN-IPSL (2005)5 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 6 6 !! $Id$ 7 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt7 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 8 8 !!---------------------------------------------------------------------- 9 9 -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_c1d.h90
- Property svn:eol-style deleted
r2146 r2528 4 4 5 5 !!---------------------------------------------------------------------- 6 !! OPA 9.0 , LOCEAN-IPSL (2005)6 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 7 7 !! $Id$ 8 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt8 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 9 9 !!---------------------------------------------------------------------- 10 10 … … 17 17 !! ** Method : 1D eddy viscosity coefficients ( depth ) 18 18 !! ahm3, ahm4 never used 19 !! biharmonic or harmonic operator : ahm1=ahm2 defined at T-level 19 !! harmonic operator : ahm1 defined at T-level 20 !! biharmonic operator : ahm2 defined at T-level 20 21 !! isopycnal or geopotential harmonic operator 21 22 !! : ahm1 defined at T-level … … 27 28 28 29 !! * Local variables 29 INTEGER :: jk ! dummy loop indice 30 REAL(wp) :: zdam, zwam, zm00, zm01, zmhf, zmhs 30 INTEGER :: jk ! dummy loop indice 31 REAL(wp) :: zdam, zwam, zm00, zm01, zmhf, zmhs 32 REAL(wp) :: zdam2, zwam2, zm200, zm201, zmh2f, zmh2s 31 33 REAL(wp) :: zahmf, zahms 32 34 !!---------------------------------------------------------------------- … … 37 39 IF(lwp) WRITE(numout,*) 38 40 39 ! Set ahm1 =ahm2(always at t-level)41 ! Set ahm1 for laplacian (always at t-level) 40 42 ! ============= 41 43 ! (USER: modify ahm1 following your desiderata) … … 54 56 zmhs = zahms-zmhf * zm00 55 57 56 ! set ahm1=ahm2 at T-level 58 ! Set ahm2 for bilaplacian (always at t-level) 59 ! ============= 60 ! (USER: modify ahm2 following your desiderata) 61 62 ! initialization of the profile 63 ! ahms, ahmf: surface and bottom values 64 zahm2s = ahm0_blp 65 zahm2f = ahm0_blp/4. 66 ! zdam, zwam: depth of the inflection pt and width of inflection 67 zdam2 = -300. 68 zwam2 = 300. 69 ! computation coefficients 70 zm200 = TANH( (0-zdam2)/zwam2 ) 71 zm201 = TANH( (-fsdept(1,1,jpk)-zdam2)/zwam2 ) 72 zmh2f = (zahm2s-zahm2f)/(zm200-zm201) 73 zmh2s = zahm2s-zmh2f * zm00 74 75 76 ! set ahm1 and ahm2 at T-level 57 77 DO jk = 1, jpk 58 ahm1(jk) = zmhs + zmhf * TANH( (-fsdept(1,1,jk)-zdam) / zwam)59 ahm2(jk) = ahm1(jk)78 ahm1(jk) = zmhs + zmhf * TANH( (-fsdept(1,1,jk)-zdam ) / zwam ) 79 ahm2(jk) = zmh2s + zmh2f * TANH( (-fsdept(1,1,jk)-zdam2) / zwam2 ) 60 80 END DO 61 81 … … 63 83 IF(lwp .AND. ld_print ) THEN 64 84 WRITE(numout,*) 65 WRITE(numout,*) ' ahm profile : '85 WRITE(numout,*) ' ahm profile (laplacian): ' 66 86 WRITE(numout,*) 67 87 WRITE(numout,9100) 68 88 DO jk = 1, jpk 69 89 WRITE(numout,9110) jk, ahm1(jk), fsdept(1,1,jk) 90 END DO 91 WRITE(numout,*) 92 WRITE(numout,*) ' ahm profile (bilaplacian): ' 93 WRITE(numout,*) 94 WRITE(numout,9100) 95 DO jk = 1, jpk 96 WRITE(numout,9110) jk, ahm2(jk), fsdept(1,1,jk) 70 97 END DO 71 98 ENDIF -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_c2d.h90
- Property svn:eol-style deleted
- Property svn:executable deleted
r1694 r2528 7 7 8 8 !!---------------------------------------------------------------------- 9 !! OPA 9.0 , LOCEAN-IPSL (2005)9 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 10 10 !! $Id$ 11 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt11 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 12 12 !!---------------------------------------------------------------------- 13 13 … … 27 27 !! iso-model level : ahm3, ahm4 not used 28 28 !! 29 !! biharmonic operator : ahm 1is defined at u-point30 !! ahm 2is defined at v-point31 !! : ahm 3, ahm4not used29 !! biharmonic operator : ahm3 is defined at u-point 30 !! ahm4 is defined at v-point 31 !! : ahm1, ahm2 not used 32 32 !! 33 33 !!---------------------------------------------------------------------- … … 74 74 ENDIF 75 75 76 ! Special case for ORCA R 2 and R4 configurations (overwrite the value of ahm1 ahm2)76 ! Special case for ORCA R1, R2 and R4 configurations (overwrite the value of ahm1 ahm2) 77 77 ! ============================================== 78 78 IF( cp_cfg == "orca" .AND. ( jp_cfg == 2 .OR. jp_cfg == 4 ) ) CALL ldf_dyn_c2d_orca( ld_print ) 79 IF( cp_cfg == "orca" .AND. jp_cfg == 1) CALL ldf_dyn_c2d_orca_R1( ld_print ) 79 80 80 81 ! Control print … … 102 103 IF(lwp) WRITE(numout,*) ' maximum grid-spacing = ', zd_max, ' maximum value for ahm = ', ahm0 103 104 104 za00 = ahm0 / ( zd_max * zd_max * zd_max )105 za00 = ahm0_blp / ( zd_max * zd_max * zd_max ) 105 106 DO jj = 1, jpj 106 107 DO ji = 1, jpi … … 288 289 289 290 END SUBROUTINE ldf_dyn_c2d_orca 291 292 SUBROUTINE ldf_dyn_c2d_orca_R1( ld_print ) 293 !!---------------------------------------------------------------------- 294 !! *** ROUTINE ldf_dyn_c2d *** 295 !! 296 !! **** W A R N I N G **** 297 !! 298 !! ORCA R1 configuration 299 !! 300 !! **** W A R N I N G **** 301 !! 302 !! ** Purpose : initializations of the lateral viscosity for orca R1 303 !! 304 !! ** Method : blah blah blah... 305 !! 306 !!---------------------------------------------------------------------- 307 !! * Modules used 308 USE ldftra_oce, ONLY : aht0 309 310 !! * Arguments 311 LOGICAL, INTENT (in) :: ld_print ! If true, output arrays on numout 312 313 !! * Local variables 314 INTEGER :: ji, jj, jn ! dummy loop indices 315 INTEGER :: inum ! temporary logical unit 316 INTEGER :: iim, ijm 317 INTEGER :: ifreq, il1, il2, ij, ii 318 INTEGER, DIMENSION(jpidta,jpidta) :: idata 319 INTEGER, DIMENSION(jpi ,jpj ) :: icof 320 321 REAL(wp) :: zahmeq, zcoft, zcoff, zmsk, zam20s 322 323 CHARACTER (len=15) :: clexp 324 !!---------------------------------------------------------------------- 325 326 IF(lwp) WRITE(numout,*) 327 IF(lwp) WRITE(numout,*) 'inildf: 2d eddy viscosity coefficient' 328 IF(lwp) WRITE(numout,*) '~~~~~~ --' 329 IF(lwp) WRITE(numout,*) 330 IF(lwp) WRITE(numout,*) ' orca_r1 ocean model' 331 IF(lwp) WRITE(numout,*) 332 333 #if defined key_antarctic 334 # include "ldfdyn_antarctic.h90" 335 #elif defined key_arctic 336 # include "ldfdyn_arctic.h90" 337 #else 338 ! Read 2d integer array to specify western boundary increase in the 339 ! ===================== equatorial strip (20N-20S) defined at t-points 340 341 CALL ctl_opn( inum, 'ahmcoef', 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL', & 342 & 1, numout, lwp ) 343 REWIND inum 344 READ(inum,9101) clexp, iim, ijm 345 READ(inum,'(/)') 346 ifreq = 40 347 il1 = 1 348 DO jn = 1, jpidta/ifreq+1 349 READ(inum,'(/)') 350 il2 = MIN( jpidta, il1+ifreq-1 ) 351 READ(inum,9201) ( ii, ji = il1, il2, 5 ) 352 READ(inum,'(/)') 353 DO jj = jpjdta, 1, -1 354 READ(inum,9202) ij, ( idata(ji,jj), ji = il1, il2 ) 355 END DO 356 il1 = il1 + ifreq 357 END DO 358 359 DO jj = 1, nlcj 360 DO ji = 1, nlci 361 icof(ji,jj) = idata( mig(ji), mjg(jj) ) 362 END DO 363 END DO 364 DO jj = nlcj+1, jpj 365 DO ji = 1, nlci 366 icof(ji,jj) = icof(ji,nlcj) 367 END DO 368 END DO 369 DO jj = 1, jpj 370 DO ji = nlci+1, jpi 371 icof(ji,jj) = icof(nlci,jj) 372 END DO 373 END DO 374 375 9101 FORMAT(1x,a15,2i8) 376 9201 FORMAT(3x,13(i3,12x)) 377 9202 FORMAT(i3,41i3) 378 379 380 ! Set ahm1 and ahm2 ( T- and F- points) (used for laplacian operator) 381 ! ================= 382 ! define ahm1 and ahm2 at the right grid point position 383 ! (USER: modify ahm1 and ahm2 following your desiderata) 384 385 386 ! Decrease ahm to zahmeq m2/s in the tropics 387 ! (from 90 to 20 degrees: ahm = scaled by local metrics 388 ! from 20 to 2.5 degrees: ahm = decrease in (1-cos)/2 389 ! from 2.5 to 0 degrees: ahm = constant 390 ! symmetric in the south hemisphere) 391 392 zahmeq = aht0 393 zam20s = ahm0*COS( rad * 20. ) 394 395 DO jj = 1, jpj 396 DO ji = 1, jpi 397 IF( ABS( gphif(ji,jj) ) >= 20. ) THEN 398 ! leave as set in ldf_dyn_c2d 399 ELSEIF( ABS( gphif(ji,jj) ) <= 2.5 ) THEN 400 ahm2(ji,jj) = zahmeq 401 ELSE 402 ahm2(ji,jj) = zahmeq + (zam20s-zahmeq)/2. & 403 * ( 1. - COS( rad * ( ABS(gphif(ji,jj))-2.5 ) * 180. / 17.5 ) ) 404 ENDIF 405 IF( ABS( gphit(ji,jj) ) >= 20. ) THEN 406 ! leave as set in ldf_dyn_c2d 407 ELSEIF( ABS( gphit(ji,jj) ) <= 2.5 ) THEN 408 ahm1(ji,jj) = zahmeq 409 ELSE 410 ahm1(ji,jj) = zahmeq + (zam20s-zahmeq)/2. & 411 * ( 1. - COS( rad * ( ABS(gphit(ji,jj))-2.5 ) * 180. / 17.5 ) ) 412 ENDIF 413 END DO 414 END DO 415 416 ! increase along western boundaries of equatorial strip 417 ! t-point 418 DO jj = 1, jpjm1 419 DO ji = 1, jpim1 420 IF( ABS( gphit(ji,jj) ) < 20. ) THEN 421 zcoft = FLOAT( icof(ji,jj) ) / 100. 422 ahm1(ji,jj) = zcoft * ahm0 + (1.-zcoft) * ahm1(ji,jj) 423 ENDIF 424 END DO 425 END DO 426 ! f-point 427 icof(:,:) = icof(:,:) * tmask(:,:,1) 428 DO jj = 1, jpjm1 429 DO ji = 1, jpim1 430 IF( ABS( gphif(ji,jj) ) < 20. ) THEN 431 zmsk = tmask(ji,jj+1,1) + tmask(ji+1,jj+1,1) + tmask(ji,jj,1) + tmask(ji,jj+1,1) 432 IF( zmsk == 0. ) THEN 433 zcoff = 1. 434 ELSE 435 zcoff = FLOAT( icof(ji,jj+1) + icof(ji+1,jj+1) + icof(ji,jj) + icof(ji,jj+1) ) & 436 / (zmsk * 100.) 437 ENDIF 438 ahm2(ji,jj) = zcoff * ahm0 + (1.-zcoff) * ahm2(ji,jj) 439 ENDIF 440 END DO 441 END DO 442 #endif 443 444 ! Lateral boundary conditions on ( ahm1, ahm2 ) 445 ! ============== 446 CALL lbc_lnk( ahm1, 'T', 1. ) ! T-point, unchanged sign 447 CALL lbc_lnk( ahm2, 'F', 1. ) ! F-point, unchanged sign 448 449 ! Control print 450 IF( lwp .AND. ld_print ) THEN 451 WRITE(numout,*) 452 WRITE(numout,*) 'inildf: 2D ahm1 array' 453 CALL prihre(ahm1,jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 454 WRITE(numout,*) 455 WRITE(numout,*) 'inildf: 2D ahm2 array' 456 CALL prihre(ahm2,jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 457 ENDIF 458 459 END SUBROUTINE ldf_dyn_c2d_orca_R1 -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_c3d.h90
- Property svn:eol-style deleted
- Property svn:executable deleted
r1694 r2528 4 4 5 5 !!---------------------------------------------------------------------- 6 !! OPA 9.0 , LOCEAN-IPSL (2005)6 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 7 7 !! $Id$ 8 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt8 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 9 9 !!---------------------------------------------------------------------- 10 10 … … 26 26 !! ??? explanation of the default is missing 27 27 !!---------------------------------------------------------------------- 28 !! * Modules used29 28 USE ldftra_oce, ONLY : aht0 30 31 !! * Arguments 29 !! 32 30 LOGICAL, INTENT (in) :: ld_print ! If true, output arrays on numout 33 34 !! * local variables 31 !! 35 32 INTEGER :: ji, jj, jk ! dummy loop indices 36 33 REAL(wp) :: & … … 86 83 87 84 88 ! Special case for ORCA R 2 and R4 configurations (overwrite the value of ahm1 ahm2)85 ! Special case for ORCA R1, R2 and R4 configurations (overwrite the value of ahm1 ahm2) 89 86 ! ============================================== 90 IF( cp_cfg == "orca" .AND. ( jp_cfg == 2 .OR. jp_cfg == 4 ) ) THEN87 IF( cp_cfg == "orca" .AND. ( jp_cfg == 1 .OR. jp_cfg == 2 .OR. jp_cfg == 4 ) ) THEN 91 88 IF(lwp) WRITE(numout,*) 92 IF(lwp) WRITE(numout,*) ' ORCA R 2 or R4: overwrite the previous definition of ahm'93 IF(lwp) WRITE(numout,*) ' ============= '89 IF(lwp) WRITE(numout,*) ' ORCA R1, R2 or R4: overwrite the previous definition of ahm' 90 IF(lwp) WRITE(numout,*) ' =================' 94 91 CALL ldf_dyn_c3d_orca( ld_print ) 95 92 ENDIF … … 122 119 IF(lwp) WRITE(numout,*) ' maximum grid-spacing = ', zd_max, ' maximum value for ahm = ', ahm0 123 120 124 za00 = ahm0 / ( zd_max * zd_max * zd_max )121 za00 = ahm0_blp / ( zd_max * zd_max * zd_max ) 125 122 DO jj = 1, jpj 126 123 DO ji = 1, jpi … … 148 145 END DO 149 146 ELSE ! partial steps or s-ccordinate 150 # if defined key_zco151 zc = gdept_0(jpkm1)152 # else153 147 zc = MAXVAL( fsdept(:,:,jpkm1) ) 154 # endif155 148 IF( lk_mpp ) CALL mpp_max( zc ) ! max over the global domain 156 149 … … 196 189 !! *** ROUTINE ldf_dyn_c3d *** 197 190 !! 198 !! ** Purpose : ORCA R 2 anR4 only191 !! ** Purpose : ORCA R1, R2 and R4 only 199 192 !! 200 193 !! ** Method : blah blah blah .... 201 194 !!---------------------------------------------------------------------- 202 !! * Modules used203 195 USE ldftra_oce, ONLY : aht0 204 205 !! * Arguments 196 !! 206 197 LOGICAL, INTENT (in) :: ld_print ! If true, output arrays on numout 207 208 !! * local variables 198 !! 209 199 INTEGER :: ji, jj, jk, jn ! dummy loop indices 210 200 INTEGER :: ii0, ii1, ij0, ij1 ! temporary integers … … 228 218 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~' 229 219 IF(lwp) WRITE(numout,*) 230 IF(lwp) WRITE(numout,*) ' orca R 2 or R4 ocean model'220 IF(lwp) WRITE(numout,*) ' orca R1, R2 or R4 ocean model' 231 221 IF(lwp) WRITE(numout,*) ' reduced in the surface Eq. strip ' 232 222 IF(lwp) WRITE(numout,*) … … 318 308 ENDIF 319 309 310 IF( jp_cfg == 1 ) THEN 311 zahmeq = aht0 ! reduced to aht0 on equator; set to ahm0 if no tropical reduction is required 312 zahmm = ahm0 313 zahm0(:,:) = ahm0 314 ENDIF 315 320 316 DO jj = 1, jpj 321 317 DO ji = 1, jpi … … 363 359 364 360 ! other level: re-increase the coef in the deep ocean 365 366 #if defined key_orca_lev10 367 DO jk = 1, 210 368 zcoef(jk) = 1. 369 END DO 370 DO jk= 211, 230 371 zcoef(jk) = 1. + 0.1 * FLOAT(jk-210) 372 END DO 373 DO jk= 231, 260 374 zcoef(jk) = 3. + 0.2 * FLOAT(jk-230) 375 END DO 376 DO jk= 261, 270 377 zcoef(jk) = 9. + 0.1 * FLOAT(jk-260) 378 END DO 379 DO jk= 271, jpk 380 zcoef(jk) = 10. 381 END DO 382 DO jk= 1, jpk 383 IF(lwp) WRITE(numout,*) 'k= ',jk, 'cof ', zcoef(jk) 384 END DO 385 #else 386 DO jk = 1, 21 387 zcoef(jk) = 1. 388 END DO 389 zcoef(22) = 2. 390 zcoef(23) = 3. 391 zcoef(24) = 5. 392 zcoef(25) = 7. 393 zcoef(26) = 9. 394 DO jk = 27, jpk 395 zcoef(jk) = 10. 396 END DO 397 #endif 361 !================================================================== 362 ! Prior to v3.3, zcoeff was hardwired according to k-index jk. 363 ! 364 ! From v3.3 onwards this has been generalised to a function of 365 ! depth so that it can be used with any number of levels. 366 ! 367 ! The function has been chosen to match the original values (shown 368 ! in the following comments) when using the standard 31 ORCA levels. 369 ! DO jk = 1, 21 370 ! zcoef(jk) = 1._wp 371 ! END DO 372 ! zcoef(22) = 2._wp 373 ! zcoef(23) = 3._wp 374 ! zcoef(24) = 5._wp 375 ! zcoef(25) = 7._wp 376 ! zcoef(26) = 9._wp 377 ! DO jk = 27, jpk 378 ! zcoef(jk) = 10._wp 379 ! END DO 380 !================================================================== 381 382 IF(lwp) THEN 383 WRITE(numout,*) 384 WRITE(numout,*) ' 1D zcoef array ' 385 WRITE(numout,*) ' ~~~~~~~~~~~~~~ ' 386 WRITE(numout,*) 387 WRITE(numout,*) ' jk zcoef ' 388 ENDIF 389 390 DO jk=1, jpk 391 zcoef(jk) = 1.0_wp + NINT(9.0_wp*(gdept_0(jk)-800.0_wp)/(3000.0_wp-800.0_wp)) 392 zcoef(jk) = MIN(10.0_wp, MAX(1.0_wp, zcoef(jk))) 393 IF(lwp) WRITE(numout,'(4x,i3,6x,f7.3)') jk,zcoef(jk) 394 END DO 398 395 399 396 DO jk = 2, jpk … … 402 399 END DO 403 400 404 IF( jp_cfg == 4 ) THEN 405 ! Limit AHM in Gibraltar strait 401 IF( jp_cfg == 4 ) THEN ! Limit AHM in Gibraltar strait 406 402 ij0 = 50 ; ij1 = 53 407 403 ii0 = 69 ; ii1 = 71 408 404 DO jk = 1, jpk 409 ahm1(mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk) = min( zahmm, ahm1 (mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk) ) 410 ahm2(mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk) = min( zahmm, ahm2 (mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , jk) ) 411 END DO 412 ENDIF 413 414 ! Lateral boundary conditions on ( ahm1, ahm2 ) 415 ! ============== 416 CALL lbc_lnk( ahm1, 'T', 1. ) ! T-point, unchanged sign 417 CALL lbc_lnk( ahm2, 'F', 1. ) ! F-point, unchanged sign 418 419 ! Control print 420 421 IF(lwp) THEN 405 ahm1(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),jk) = MIN( zahmm, ahm1(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),jk) ) 406 ahm2(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),jk) = MIN( zahmm, ahm2(mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),jk) ) 407 END DO 408 ENDIF 409 CALL lbc_lnk( ahm1, 'T', 1. ) ! Lateral boundary conditions (unchanged sign) 410 CALL lbc_lnk( ahm2, 'F', 1. ) 411 412 413 IF(lwp) THEN ! Control print 422 414 WRITE(numout,*) 423 415 WRITE(numout,*) ' 3D ahm1 array (k=1)' … … 451 443 ahm4 ( :, 1, :) = ahm4 ( :, 2, :) 452 444 453 ! Lateral boundary conditions on ( ahm3, ahm4 ) 454 ! ============== 455 CALL lbc_lnk( ahm3, 'U', 1. ) ! U-point, unchanged sign 456 CALL lbc_lnk( ahm4, 'V', 1. ) ! V-point, unchanged sign 445 CALL lbc_lnk( ahm3, 'U', 1. ) ! Lateral boundary conditions (unchanged sign) 446 CALL lbc_lnk( ahm4, 'V', 1. ) 457 447 458 448 ! Control print -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_oce.F90
- Property svn:eol-style deleted
r1601 r2528 17 17 LOGICAL , PUBLIC :: ln_dynldf_hor = .TRUE. !: horizontal (geopotential) direction 18 18 LOGICAL , PUBLIC :: ln_dynldf_iso = .FALSE. !: iso-neutral direction 19 REAL(wp), PUBLIC :: rn_ahm_0 = 40000._wp !: lateraleddy viscosity (m2/s)20 REAL(wp), PUBLIC :: rn_ahmb_0 = 0._wp !: lateral background eddy viscosity (m2/s)21 22 REAL(wp), PUBLIC :: ahm0, ahmb0 19 REAL(wp), PUBLIC :: rn_ahm_0_lap = 40000._wp !: lateral laplacian eddy viscosity (m2/s) 20 REAL(wp), PUBLIC :: rn_ahmb_0 = 0._wp !: lateral laplacian background eddy viscosity (m2/s) 21 REAL(wp), PUBLIC :: rn_ahm_0_blp = 0._wp !: lateral bilaplacian eddy viscosity (m4/s) 22 REAL(wp), PUBLIC :: ahm0, ahmb0, ahm0_blp ! OLD namelist names 23 23 24 24 #if defined key_dynldf_c3d … … 33 33 34 34 !!---------------------------------------------------------------------- 35 !! NEMO/OPA 3. 2 , LOCEAN-IPSL (2009)35 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 36 36 !! $Id$ 37 !! Software governed by the CeCILL licence ( modipsl/doc/NEMO_CeCILL.txt)37 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 38 38 !!====================================================================== 39 39 END MODULE ldfdyn_oce -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn_substitute.h90
- Property svn:eol-style deleted
r1152 r2528 6 6 !!---------------------------------------------------------------------- 7 7 !!---------------------------------------------------------------------- 8 !! OPA 9.0 , LOCEAN-IPSL (2005)8 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 9 9 !! $Id$ 10 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt10 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 11 11 !!---------------------------------------------------------------------- 12 !! 13 !! fsahmt, fsahmf - used for laplaian operator only 14 !! fsahmu, fsahmv - used for bilaplacian operator only 15 !! 12 16 #if defined key_dynldf_c3d 13 17 ! ' key_dynldf_c3d' : 3D coefficient … … 26 30 # define fsahmt(i,j,k) ahm1(k) 27 31 # define fsahmf(i,j,k) ahm1(k) 28 # define fsahmu(i,j,k) ahm 1(k)29 # define fsahmv(i,j,k) ahm 1(k)32 # define fsahmu(i,j,k) ahm2(k) 33 # define fsahmv(i,j,k) ahm2(k) 30 34 #else 31 35 ! default option : Constant coefficient 32 36 # define fsahmt(i,j,k) ahm0 33 37 # define fsahmf(i,j,k) ahm0 34 # define fsahmu(i,j,k) ahm0 35 # define fsahmv(i,j,k) ahm0 38 # define fsahmu(i,j,k) ahm0_blp 39 # define fsahmv(i,j,k) ahm0_blp 36 40 #endif -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfeiv.F90
- Property svn:eol-style deleted
r1482 r2528 4 4 !! Ocean physics: variable eddy induced velocity coefficients 5 5 !!====================================================================== 6 !! History : OPA ! 1999-03 (G. Madec, A. Jouzeau) Original code 7 !! NEMO 1.0 ! 2002-06 (G. Madec) Free form, F90 8 !!---------------------------------------------------------------------- 6 9 #if defined key_traldf_eiv && defined key_traldf_c2d 7 10 !!---------------------------------------------------------------------- … … 11 14 !! ldf_eiv : compute the eddy induced velocity coefficients 12 15 !!---------------------------------------------------------------------- 13 !! * Modules used14 16 USE oce ! ocean dynamics and tracers 15 17 USE dom_oce ! ocean space and time domain … … 22 24 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 23 25 USE prtctl ! Print control 24 USE iom 26 USE iom ! I/O library 25 27 26 28 IMPLICIT NONE 27 29 PRIVATE 28 30 29 !! * Routine accessibility 30 PUBLIC ldf_eiv ! routine called by step.F90 31 !!---------------------------------------------------------------------- 32 !! OPA 9.0 , LOCEAN-IPSL (2005) 33 !! $Id$ 34 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 35 !!---------------------------------------------------------------------- 31 PUBLIC ldf_eiv ! routine called by step.F90 32 36 33 !! * Substitutions 37 34 # include "domzgr_substitute.h90" 38 35 # include "vectopt_loop_substitute.h90" 39 36 !!---------------------------------------------------------------------- 40 37 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 38 !! $Id$ 39 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 40 !!---------------------------------------------------------------------- 41 41 CONTAINS 42 42 … … 46 46 !! 47 47 !! ** Purpose : Compute the eddy induced velocity coefficient from the 48 !! growth rate of baroclinic instability.48 !! growth rate of baroclinic instability. 49 49 !! 50 50 !! ** Method : 51 51 !! 52 !! ** Action : - uslp(), : i- and j-slopes of neutral surfaces 53 !! - vslp() at u- and v-points, resp. 54 !! - wslpi(), : i- and j-slopes of neutral surfaces 55 !! - wslpj() at w-points. 56 !! 57 !! History : 58 !! 8.1 ! 99-03 (G. Madec, A. Jouzeau) Original code 59 !! 8.5 ! 02-06 (G. Madec) Free form, F90 52 !! ** Action : - uslp , vslp : i- and j-slopes of neutral surfaces at u- & v-points 53 !! - wslpi, wslpj : i- and j-slopes of neutral surfaces at w-points. 60 54 !!---------------------------------------------------------------------- 61 !! * Arguments 62 INTEGER, INTENT( in ) :: kt ! ocean time-step inedx 63 64 !! * Local declarations 65 INTEGER :: ji, jj, jk ! dummy loop indices 66 REAL(wp) :: & 67 zfw, ze3w, zn2, zf20, & ! temporary scalars 68 zaht, zaht_min 69 REAL(wp), DIMENSION(jpi,jpj) :: & 70 zn, zah, zhw, zross ! workspace 55 INTEGER, INTENT(in) :: kt ! ocean time-step inedx 56 !! 57 INTEGER :: ji, jj, jk ! dummy loop indices 58 REAL(wp) :: zfw, ze3w, zn2, zf20, zaht, zaht_min ! temporary scalars 59 REAL(wp), DIMENSION(jpi,jpj) :: zn, zah, zhw, zross ! 2D workspace 71 60 !!---------------------------------------------------------------------- 72 61 … … 79 68 ! 0. Local initialization 80 69 ! ----------------------- 81 zn (:,:) = 0. e082 zhw (:,:) = 5. e083 zah (:,:) = 0. e084 zross(:,:) = 0. e070 zn (:,:) = 0._wp 71 zhw (:,:) = 5._wp 72 zah (:,:) = 0._wp 73 zross(:,:) = 0._wp 85 74 86 75 87 76 ! 1. Compute lateral diffusive coefficient 88 77 ! ---------------------------------------- 89 90 DO jk = 1, jpk78 IF( ln_traldf_grif ) THEN 79 DO jk = 1, jpk 91 80 # if defined key_vectopt_loop 92 81 !CDIR NOVERRCHK 93 DO ji = 1, jpij ! vector opt. 94 ! Take the max of N^2 and zero then take the vertical sum 95 ! of the square root of the resulting N^2 ( required to compute 96 ! internal Rossby radius Ro = .5 * sum_jpk(N) / f 97 zn2 = MAX( rn2b(ji,1,jk), 0.e0 ) 98 zn(ji,1) = zn(ji,1) + SQRT( zn2 ) * fse3w(ji,1,jk) 99 ! Compute elements required for the inverse time scale of baroclinic 100 ! eddies using the isopycnal slopes calculated in ldfslp.F : 101 ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w)) 102 ze3w = fse3w(ji,1,jk) * tmask(ji,1,jk) 103 zah(ji,1) = zah(ji,1) + zn2 & 104 * ( wslpi(ji,1,jk) * wslpi(ji,1,jk) & 105 + wslpj(ji,1,jk) * wslpj(ji,1,jk) ) & 106 * ze3w 107 zhw(ji,1) = zhw(ji,1) + ze3w 108 END DO 82 DO ji = 1, jpij ! vector opt. 83 ! Take the max of N^2 and zero then take the vertical sum 84 ! of the square root of the resulting N^2 ( required to compute 85 ! internal Rossby radius Ro = .5 * sum_jpk(N) / f 86 zn2 = MAX( rn2b(ji,1,jk), 0._wp ) 87 zn(ji,1) = zn(ji,1) + SQRT( zn2 ) * fse3w(ji,1,jk) 88 ! Compute elements required for the inverse time scale of baroclinic 89 ! eddies using the isopycnal slopes calculated in ldfslp.F : 90 ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w)) 91 ze3w = fse3w(ji,1,jk) * tmask(ji,1,jk) 92 zah(ji,1) = zah(ji,1) + zn2 * wslp2(ji,1,jk) * ze3w 93 zhw(ji,1) = zhw(ji,1) + ze3w 94 END DO 109 95 # else 110 DO jj = 2, jpjm1 111 !CDIR NOVERRCHK 112 DO ji = 2, jpim1 113 ! Take the max of N^2 and zero then take the vertical sum 114 ! of the square root of the resulting N^2 ( required to compute 115 ! internal Rossby radius Ro = .5 * sum_jpk(N) / f 116 zn2 = MAX( rn2b(ji,jj,jk), 0.e0 ) 117 zn(ji,jj) = zn(ji,jj) + SQRT( zn2 ) * fse3w(ji,jj,jk) 96 DO jj = 2, jpjm1 97 !CDIR NOVERRCHK 98 DO ji = 2, jpim1 99 ! Take the max of N^2 and zero then take the vertical sum 100 ! of the square root of the resulting N^2 ( required to compute 101 ! internal Rossby radius Ro = .5 * sum_jpk(N) / f 102 zn2 = MAX( rn2b(ji,jj,jk), 0._wp ) 103 zn(ji,jj) = zn(ji,jj) + SQRT( zn2 ) * fse3w(ji,jj,jk) 104 ! Compute elements required for the inverse time scale of baroclinic 105 ! eddies using the isopycnal slopes calculated in ldfslp.F : 106 ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w)) 107 ze3w = fse3w(ji,jj,jk) * tmask(ji,jj,jk) 108 zah(ji,jj) = zah(ji,jj) + zn2 * wslp2(ji,jj,jk) * ze3w 109 zhw(ji,jj) = zhw(ji,jj) + ze3w 110 END DO 111 END DO 112 # endif 113 END DO 114 ELSE 115 DO jk = 1, jpk 116 # if defined key_vectopt_loop 117 !CDIR NOVERRCHK 118 DO ji = 1, jpij ! vector opt. 119 ! Take the max of N^2 and zero then take the vertical sum 120 ! of the square root of the resulting N^2 ( required to compute 121 ! internal Rossby radius Ro = .5 * sum_jpk(N) / f 122 zn2 = MAX( rn2b(ji,1,jk), 0._wp ) 123 zn(ji,1) = zn(ji,1) + SQRT( zn2 ) * fse3w(ji,1,jk) 118 124 ! Compute elements required for the inverse time scale of baroclinic 119 ! eddies using the isopycnal slopes calculated in ldfslp.F : 125 ! eddies using the isopycnal slopes calculated in ldfslp.F : 120 126 ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w)) 121 ze3w = fse3w(ji,jj,jk) * tmask(ji,jj,jk) 122 zah(ji,jj) = zah(ji,jj) + zn2 & 123 * ( wslpi(ji,jj,jk) * wslpi(ji,jj,jk) & 124 + wslpj(ji,jj,jk) * wslpj(ji,jj,jk) ) & 125 * ze3w 126 zhw(ji,jj) = zhw(ji,jj) + ze3w 127 END DO 128 END DO 127 ze3w = fse3w(ji,1,jk) * tmask(ji,1,jk) 128 zah(ji,1) = zah(ji,1) + zn2 * ( wslpi(ji,1,jk) * wslpi(ji,1,jk) & 129 & + wslpj(ji,1,jk) * wslpj(ji,1,jk) ) * ze3w 130 zhw(ji,1) = zhw(ji,1) + ze3w 131 END DO 132 # else 133 DO jj = 2, jpjm1 134 !CDIR NOVERRCHK 135 DO ji = 2, jpim1 136 ! Take the max of N^2 and zero then take the vertical sum 137 ! of the square root of the resulting N^2 ( required to compute 138 ! internal Rossby radius Ro = .5 * sum_jpk(N) / f 139 zn2 = MAX( rn2b(ji,jj,jk), 0._wp ) 140 zn(ji,jj) = zn(ji,jj) + SQRT( zn2 ) * fse3w(ji,jj,jk) 141 ! Compute elements required for the inverse time scale of baroclinic 142 ! eddies using the isopycnal slopes calculated in ldfslp.F : 143 ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w)) 144 ze3w = fse3w(ji,jj,jk) * tmask(ji,jj,jk) 145 zah(ji,jj) = zah(ji,jj) + zn2 * ( wslpi(ji,jj,jk) * wslpi(ji,jj,jk) & 146 & + wslpj(ji,jj,jk) * wslpj(ji,jj,jk) ) * ze3w 147 zhw(ji,jj) = zhw(ji,jj) + ze3w 148 END DO 149 END DO 129 150 # endif 130 END DO 151 END DO 152 END IF 131 153 132 154 DO jj = 2, jpjm1 … … 141 163 END DO 142 164 143 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN ! ORCA R 02165 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN ! ORCA R2 144 166 DO jj = 2, jpjm1 145 167 !CDIR NOVERRCHK 146 168 DO ji = fs_2, fs_jpim1 ! vector opt. 147 ! Take the minimum between aeiw and aeiv0 for depth levels 148 ! lower than 20 (21 in w- point) 149 IF( mbathy(ji,jj) <= 21. ) aeiw(ji,jj) = MIN( aeiw(ji,jj), 1000. ) 169 ! Take the minimum between aeiw and 1000 m2/s over shelves (depth shallower than 650 m) 170 IF( mbkt(ji,jj) <= 20 ) aeiw(ji,jj) = MIN( aeiw(ji,jj), 1000. ) 150 171 END DO 151 172 END DO … … 153 174 154 175 ! Decrease the coefficient in the tropics (20N-20S) 155 zf20 = 2. * omega * sin( rad * 20.)176 zf20 = 2._wp * omega * sin( rad * 20._wp ) 156 177 DO jj = 2, jpjm1 157 178 DO ji = fs_2, fs_jpim1 ! vector opt. … … 168 189 END DO 169 190 ENDIF 170 171 ! lateral boundary condition on aeiw 172 CALL lbc_lnk( aeiw, 'W', 1. ) 191 CALL lbc_lnk( aeiw, 'W', 1. ) ! lateral boundary condition on aeiw 192 173 193 174 194 ! Average the diffusive coefficient at u- v- points 175 195 DO jj = 2, jpjm1 176 196 DO ji = fs_2, fs_jpim1 ! vector opt. 177 aeiu(ji,jj) = .5* ( aeiw(ji,jj) + aeiw(ji+1,jj ) )178 aeiv(ji,jj) = .5* ( aeiw(ji,jj) + aeiw(ji ,jj+1) )197 aeiu(ji,jj) = 0.5_wp * ( aeiw(ji,jj) + aeiw(ji+1,jj ) ) 198 aeiv(ji,jj) = 0.5_wp * ( aeiw(ji,jj) + aeiw(ji ,jj+1) ) 179 199 END DO 180 200 END DO 181 182 ! lateral boundary condition on aeiu, aeiv 183 CALL lbc_lnk( aeiu, 'U', 1. ) 184 CALL lbc_lnk( aeiv, 'V', 1. ) 185 186 IF(ln_ctl) THEN 201 CALL lbc_lnk( aeiu, 'U', 1. ) ; CALL lbc_lnk( aeiv, 'V', 1. ) ! lateral boundary condition 202 203 204 IF(ln_ctl) THEN 187 205 CALL prt_ctl(tab2d_1=aeiu, clinfo1=' eiv - u: ', ovlap=1) 188 206 CALL prt_ctl(tab2d_1=aeiv, clinfo1=' eiv - v: ', ovlap=1) … … 191 209 ! ORCA R05: add a space variation on aht (=aeiv except at the equator and river mouth) 192 210 IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN 193 zf20 = 2. * omega * SIN( rad * 20.)194 zaht_min = 100. 211 zf20 = 2._wp * omega * SIN( rad * 20._wp ) 212 zaht_min = 100._wp ! minimum value for aht 195 213 DO jj = 1, jpj 196 214 DO ji = 1, jpi 197 zaht = ( 1. - MIN( 1., ABS( ff(ji,jj) / zf20 ) ) ) * ( aht0 - zaht_min ) &215 zaht = ( 1._wp - MIN( 1._wp , ABS( ff(ji,jj) / zf20 ) ) ) * ( aht0 - zaht_min ) & 198 216 & + aht0 * rnfmsk(ji,jj) ! enhanced near river mouths 199 217 ahtu(ji,jj) = MAX( MAX( zaht_min, aeiu(ji,jj) ) + zaht, aht0 ) … … 209 227 ENDIF 210 228 211 IF( aeiv0 == 0. e0) THEN212 aeiu(:,:) = 0. e0213 aeiv(:,:) = 0. e0214 aeiw(:,:) = 0. e0229 IF( aeiv0 == 0._wp ) THEN 230 aeiu(:,:) = 0._wp 231 aeiv(:,:) = 0._wp 232 aeiw(:,:) = 0._wp 215 233 ENDIF 216 234 217 235 CALL iom_put( "aht2d" , ahtw ) ! lateral eddy diffusivity 218 236 CALL iom_put( "aht2d_eiv", aeiw ) ! EIV lateral eddy diffusivity 219 237 ! 220 238 END SUBROUTINE ldf_eiv 221 239 -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfeiv_substitute.h90
- Property svn:eol-style deleted
r1152 r2528 7 7 !!---------------------------------------------------------------------- 8 8 !!---------------------------------------------------------------------- 9 !! OPA 9.0 , LOCEAN-IPSL (2005)9 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 10 10 !! $Id$ 11 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt11 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 12 12 !!---------------------------------------------------------------------- 13 13 # if defined key_traldf_c3d -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldfslp.F90
- Property svn:eol-style deleted
r2389 r2528 7 7 !! 8.0 ! 1997-06 (G. Madec) optimization, lbc 8 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 !! 3.2 ! 2010-11 (F. Dupond, G. Madec) bug correction in slopes just below the ML 9 !! NEMO 1.0 ! 2002-10 (G. Madec) Free form, F90 10 !! - ! 2005-10 (A. Beckmann) correction for s-coordinates 11 !! 3.3 ! 2010-10 (G. Nurser, C. Harris, G. Madec) add Griffies operator 12 !! - ! 2010-11 (F. Dupond, G. Madec) bug correction in slopes just below the ML 12 13 !!---------------------------------------------------------------------- 13 14 #if defined key_ldfslp || defined key_esopa … … 15 16 !! 'key_ldfslp' Rotation of lateral mixing tensor 16 17 !!---------------------------------------------------------------------- 17 !! ldf_slp : compute the slopes of iso-neutral surface 18 !! ldf_slp_mxl : compute the slopes of iso-neutral surface just below the Mixed Layer 18 !! ldf_slp_grif : calculates the triads of isoneutral slopes (Griffies operator) 19 !! ldf_slp : calculates the slopes of neutral surface (Madec operator) 20 !! ldf_slp_mxl : calculates the slopes at the base of the mixed layer (Madec operator) 19 21 !! ldf_slp_init : initialization of the slopes computation 20 22 !!---------------------------------------------------------------------- 21 23 USE oce ! ocean dynamics and tracers 22 24 USE dom_oce ! ocean space and time domain 23 USE ldftra_oce 24 USE ldfdyn_oce 25 USE ldftra_oce ! lateral diffusion: traceur 26 USE ldfdyn_oce ! lateral diffusion: dynamics 25 27 USE phycst ! physical constants 26 28 USE zdfmxl ! mixed layer depth 29 USE eosbn2 ! equation of states 27 30 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 28 31 USE in_out_manager ! I/O manager … … 32 35 PRIVATE 33 36 34 PUBLIC ldf_slp ! routine called by step.F90 35 36 LOGICAL , PUBLIC, PARAMETER :: lk_ldfslp = .TRUE. !: slopes flag 37 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: uslp, wslpi !: i_slope at U- and W-points 38 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: vslp, wslpj !: j-slope at V- and W-points 39 40 REAL(wp), DIMENSION(jpi,jpj,jpk) :: omlmask ! mask of the surface mixed layer at T-pt 41 REAL(wp), DIMENSION(jpi,jpj) :: uslpml, wslpiml ! i_slope at U- and W-points just below the mixed layer 42 REAL(wp), DIMENSION(jpi,jpj) :: vslpml, wslpjml ! j_slope at V- and W-points just below the mixed layer 37 PUBLIC ldf_slp ! routine called by step.F90 38 PUBLIC ldf_slp_grif ! routine called by step.F90 39 PUBLIC ldf_slp_init ! routine called by opa.F90 40 41 LOGICAL , PUBLIC, PARAMETER :: lk_ldfslp = .TRUE. !: slopes flag 42 ! !! Madec operator 43 REAL(wp), PUBLIC, DIMENSION(:,:,:) , ALLOCATABLE :: uslp, wslpi !: i_slope at U- and W-points 44 REAL(wp), PUBLIC, DIMENSION(:,:,:) , ALLOCATABLE :: vslp, wslpj !: j-slope at V- and W-points 45 ! !! Griffies operator 46 REAL(wp), PUBLIC, DIMENSION(:,:,:) , ALLOCATABLE :: wslp2 !: wslp**2 from Griffies quarter cells 47 REAL(wp), PUBLIC, DIMENSION(:,:,:,:,:), ALLOCATABLE :: triadi_g, triadj_g !: skew flux slopes relative to geopotentials 48 REAL(wp), PUBLIC, DIMENSION(:,:,:,:,:), ALLOCATABLE :: triadi , triadj !: isoneutral slopes relative to model-coordinate 49 50 ! !! Madec operator 51 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: omlmask ! mask of the surface mixed layer at T-pt 52 REAL(wp), DIMENSION(:,:) , ALLOCATABLE :: uslpml, wslpiml ! i_slope at U- and W-points just below the mixed layer 53 REAL(wp), DIMENSION(:,:) , ALLOCATABLE :: vslpml, wslpjml ! j_slope at V- and W-points just below the mixed layer 54 55 REAL(wp) :: repsln = 1.e-25_wp ! tiny value used as minium of di(rho), dj(rho) and dk(rho) 43 56 44 57 !! * Substitutions 45 58 # include "domzgr_substitute.h90" 59 # include "ldftra_substitute.h90" 60 # include "ldfeiv_substitute.h90" 46 61 # include "vectopt_loop_substitute.h90" 47 62 !!---------------------------------------------------------------------- 48 !! NEMO/OPA 3. 2 , LOCEAN-IPSL (2009)63 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 49 64 !! $Id$ 50 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)65 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 51 66 !!---------------------------------------------------------------------- 52 53 67 CONTAINS 54 68 … … 58 72 !! 59 73 !! ** Purpose : Compute the slopes of neutral surface (slope of isopycnal 60 !! surfaces referenced locally) ( 'key_traldfiso').74 !! surfaces referenced locally) (ln_traldf_iso=T). 61 75 !! 62 76 !! ** Method : The slope in the i-direction is computed at U- and … … 80 94 USE oce , zgru => ua ! use ua as workspace 81 95 USE oce , zgrv => va ! use va as workspace 82 USE oce , zw y=> ta ! use ta as workspace96 USE oce , zww => ta ! use ta as workspace 83 97 USE oce , zwz => sa ! use sa as workspace 84 98 !! … … 90 104 INTEGER :: ii0, ii1, iku ! temporary integer 91 105 INTEGER :: ij0, ij1, ikv ! temporary integer 92 REAL(wp) :: zeps, zmg, zm05g, zalpha ! temporary scalars 93 REAL(wp) :: zcoef1, zcoef2, zcoef3 ! - - 94 REAL(wp) :: zcofu , zcofv , zcofw ! - - 95 REAL(wp) :: zau, zbu, zai, zbi, z1u, z1wu ! - - 96 REAL(wp) :: zav, zbv, zaj, zbj, z1v, z1wv ! 97 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zww ! 3D workspace 106 REAL(wp) :: zeps, zm1_g, zm1_2g, z1_16 ! local scalars 107 REAL(wp) :: zci, zfi, zau, zbu, zai, zbi ! - - 108 REAL(wp) :: zcj, zfj, zav, zbv, zaj, zbj ! - - 109 REAL(wp) :: zck, zfk, zbw ! - - 110 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdzr ! 3D workspace 98 111 !!---------------------------------------------------------------------- 99 112 100 IF( kt == nit000 ) CALL ldf_slp_init ! initialization (first time-step only) 101 102 zeps = 1.e-20 ! Local constant initialization 103 zmg = -1.0 / grav 104 zm05g = -0.5 / grav 105 ! 106 zww(:,:,:) = 0.e0 107 zwz(:,:,:) = 0.e0 108 ! ! horizontal density gradient computation 109 DO jk = 1, jpk 113 zeps = 1.e-20_wp !== Local constant initialization ==! 114 z1_16 = 1.0_wp / 16._wp 115 zm1_g = -1.0_wp / grav 116 zm1_2g = -0.5_wp / grav 117 ! 118 zww(:,:,:) = 0._wp 119 zwz(:,:,:) = 0._wp 120 ! 121 DO jk = 1, jpk !== i- & j-gradient of density ==! 110 122 DO jj = 1, jpjm1 111 123 DO ji = 1, fs_jpim1 ! vector opt. … … 123 135 DO ji = 1, jpim1 124 136 # endif 125 iku = MIN ( mbathy(ji,jj), mbathy(ji+1,jj) ) - 1 ! last ocean level 126 ikv = MIN ( mbathy(ji,jj), mbathy(ji,jj+1) ) - 1 127 zgru(ji,jj,iku) = gru(ji,jj) 128 zgrv(ji,jj,ikv) = grv(ji,jj) 137 zgru(ji,jj,mbku(ji,jj)) = gru(ji,jj) 138 zgrv(ji,jj,mbkv(ji,jj)) = grv(ji,jj) 129 139 END DO 130 140 END DO 131 141 ENDIF 132 133 134 ! !* Local vertical density gradient evaluated from N^2 135 DO jk = 2, jpkm1 ! zwy = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point 136 zwy(:,:,jk) = zmg * ( prd(:,:,jk) + 1. ) * ( pn2 (:,:,jk) + pn2 (:,:,jk+1) ) & 137 & / MAX( tmask(:,:,jk) + tmask(:,:,jk+1), 1. ) 138 END DO 139 zwy(:,:,1) = 0.e0 ! surface value set to zero 140 141 142 CALL ldf_slp_mxl( prd, pn2 ) !* Slopes of isopycnal surfaces just below the mixed layer 143 142 ! 143 zdzr(:,:,1) = 0._wp !== Local vertical density gradient at T-point == ! (evaluated from N^2) 144 DO jk = 2, jpkm1 145 ! ! zdzr = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point 146 ! ! trick: tmask(ik ) = 0 => all pn2 = 0 => zdzr = 0 147 ! ! else tmask(ik+1) = 0 => pn2(ik+1) = 0 => zdzr divides by 1 148 ! ! umask(ik+1) /= 0 => all pn2 /= 0 => zdzr divides by 2 149 ! ! NB: 1/(tmask+1) = (1-.5*tmask) substitute a / by a * ==> faster 150 zdzr(:,:,jk) = zm1_g * ( prd(:,:,jk) + 1._wp ) & 151 & * ( pn2(:,:,jk) + pn2(:,:,jk+1) ) * ( 1._wp - 0.5_wp * tmask(:,:,jk+1) ) 152 END DO 153 ! 154 ! !== Slopes just below the mixed layer ==! 155 CALL ldf_slp_mxl( prd, pn2, zgru, zgrv, zdzr ) ! output: uslpml, vslpml, wslpiml, wslpjml 156 144 157 145 158 ! I. slopes at u and v point | uslp = d/di( prd ) / d/dz( prd ) … … 149 162 DO jj = 2, jpjm1 150 163 DO ji = fs_2, fs_jpim1 ! vector opt. 151 ! horizontal and vertical density gradient at u- and v-points 152 zau = 1. / e1u(ji,jj) * zgru(ji,jj,jk) 153 zav = 1. / e2v(ji,jj) * zgrv(ji,jj,jk) 154 zbu = 0.5 * ( zwy(ji,jj,jk) + zwy(ji+1,jj ,jk) ) 155 zbv = 0.5 * ( zwy(ji,jj,jk) + zwy(ji ,jj+1,jk) ) 156 ! bound the slopes: abs(zw.)<= 1/100 and zb..<0 157 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 158 zbu = MIN( zbu, -100.*ABS( zau ), -7.e+3/fse3u(ji,jj,jk)*ABS( zau ) ) 159 zbv = MIN( zbv, -100.*ABS( zav ), -7.e+3/fse3v(ji,jj,jk)*ABS( zav ) ) 160 ! uslp and vslp output in zwz and zww, resp. 161 zalpha = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) ) 162 zwz (ji,jj,jk) = ( ( 1. - zalpha) * zau / ( zbu - zeps ) & 163 & + zalpha * uslpml(ji,jj) & 164 & * 0.5 * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk)-fse3u(ji,jj,1) ) & 165 & / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 5. ) ) * umask(ji,jj,jk) 166 zalpha = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) ) 167 zww (ji,jj,jk) = ( ( 1. - zalpha) * zav / ( zbv - zeps ) & 168 & + zalpha * vslpml(ji,jj) & 169 & * 0.5 * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk)-fse3v(ji,jj,1) ) & 170 & / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 5. ) ) * vmask(ji,jj,jk) 164 ! ! horizontal and vertical density gradient at u- and v-points 165 zau = zgru(ji,jj,jk) / e1u(ji,jj) 166 zav = zgrv(ji,jj,jk) / e2v(ji,jj) 167 zbu = 0.5_wp * ( zdzr(ji,jj,jk) + zdzr(ji+1,jj ,jk) ) 168 zbv = 0.5_wp * ( zdzr(ji,jj,jk) + zdzr(ji ,jj+1,jk) ) 169 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0 170 ! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 171 zbu = MIN( zbu, -100._wp* ABS( zau ) , -7.e+3_wp/fse3u(ji,jj,jk)* ABS( zau ) ) 172 zbv = MIN( zbv, -100._wp* ABS( zav ) , -7.e+3_wp/fse3v(ji,jj,jk)* ABS( zav ) ) 173 ! ! uslp and vslp output in zwz and zww, resp. 174 zfi = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) ) 175 zfj = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) ) 176 zwz(ji,jj,jk) = ( ( 1. - zfi) * zau / ( zbu - zeps ) & 177 & + zfi * uslpml(ji,jj) & 178 & * 0.5_wp * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk)-fse3u(ji,jj,1) ) & 179 & / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 5._wp ) ) * umask(ji,jj,jk) 180 zww(ji,jj,jk) = ( ( 1. - zfj) * zav / ( zbv - zeps ) & 181 & + zfj * vslpml(ji,jj) & 182 & * 0.5_wp * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk)-fse3v(ji,jj,1) ) & 183 & / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 5. ) ) * vmask(ji,jj,jk) 184 !!gm modif to suppress omlmask.... (as in Griffies case) 185 ! ! ! jk must be >= ML level for zf=1. otherwise zf=0. 186 ! zfi = REAL( 1 - 1/(1 + jk / MAX( nmln(ji+1,jj), nmln(ji,jj) ) ), wp ) 187 ! zfj = REAL( 1 - 1/(1 + jk / MAX( nmln(ji,jj+1), nmln(ji,jj) ) ), wp ) 188 ! zci = 0.5 * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk) ) / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 10. ) ) 189 ! zcj = 0.5 * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk) ) / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 10. ) ) 190 ! zwz(ji,jj,jk) = ( zfi * zai / ( zbi - zeps ) + ( 1._wp - zfi ) * wslpiml(ji,jj) * zci ) * tmask(ji,jj,jk) 191 ! zww(ji,jj,jk) = ( zfj * zaj / ( zbj - zeps ) + ( 1._wp - zfj ) * wslpjml(ji,jj) * zcj ) * tmask(ji,jj,jk) 192 !!gm end modif 171 193 END DO 172 194 END DO … … 174 196 CALL lbc_lnk( zwz, 'U', -1. ) ; CALL lbc_lnk( zww, 'V', -1. ) ! lateral boundary conditions 175 197 ! 176 zcofu = 1. / 16. !* horizontal Shapiro filter 177 zcofv = 1. / 16. 198 ! !* horizontal Shapiro filter 178 199 DO jk = 2, jpkm1 179 DO jj = 2, jpjm1, jpj-3! rows jj=2 and =jpjm1 only200 DO jj = 2, jpjm1, MAX(1, jpj-3) ! rows jj=2 and =jpjm1 only 180 201 DO ji = 2, jpim1 181 uslp(ji,jj,jk) = z cofu* ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) &202 uslp(ji,jj,jk) = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 182 203 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 183 204 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 184 205 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 185 206 & + 4.* zwz(ji ,jj ,jk) ) 186 vslp(ji,jj,jk) = z cofv* ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) &207 vslp(ji,jj,jk) = z1_16 * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 187 208 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 188 209 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & … … 193 214 DO jj = 3, jpj-2 ! other rows 194 215 DO ji = fs_2, fs_jpim1 ! vector opt. 195 uslp(ji,jj,jk) = z cofu* ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) &216 uslp(ji,jj,jk) = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 196 217 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 197 218 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 198 219 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 199 220 & + 4.* zwz(ji ,jj ,jk) ) 200 vslp(ji,jj,jk) = z cofv* ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) &221 vslp(ji,jj,jk) = z1_16 * ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 201 222 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 202 223 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & … … 208 229 DO jj = 2, jpjm1 209 230 DO ji = fs_2, fs_jpim1 ! vector opt. 210 z1u = ( umask(ji,jj+1,jk) + umask(ji,jj-1,jk) )*.5 211 z1v = ( vmask(ji+1,jj,jk) + vmask(ji-1,jj,jk) )*.5 212 z1wu = ( umask(ji,jj,jk) + umask(ji,jj,jk+1) )*.5 213 z1wv = ( vmask(ji,jj,jk) + vmask(ji,jj,jk+1) )*.5 214 uslp(ji,jj,jk) = uslp(ji,jj,jk) * z1u * z1wu 215 vslp(ji,jj,jk) = vslp(ji,jj,jk) * z1v * z1wv 231 uslp(ji,jj,jk) = uslp(ji,jj,jk) * ( umask(ji,jj+1,jk) + umask(ji,jj-1,jk ) ) * 0.5_wp & 232 & * ( umask(ji,jj ,jk) + umask(ji,jj ,jk+1) ) * 0.5_wp 233 vslp(ji,jj,jk) = vslp(ji,jj,jk) * ( vmask(ji+1,jj,jk) + vmask(ji-1,jj,jk ) ) * 0.5_wp & 234 & * ( vmask(ji ,jj,jk) + vmask(ji ,jj,jk+1) ) * 0.5_wp 216 235 END DO 217 236 END DO … … 222 241 ! =========================== | wslpj = mij( d/dj( prd ) / d/dz( prd ) 223 242 ! 224 ! !* Local vertical density gradient evaluated from N^2 225 DO jk = 2, jpkm1 ! zwy = d/dz(prd)= - mk ( prd ) / grav * pn2 -- at w point 226 DO jj = 1, jpj 227 DO ji = 1, jpi 228 zwy(ji,jj,jk) = zm05g * pn2 (ji,jj,jk) * ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) 229 END DO 230 END DO 231 END DO 232 DO jk = 2, jpkm1 !* Slopes at w point 243 DO jk = 2, jpkm1 233 244 DO jj = 2, jpjm1 234 245 DO ji = fs_2, fs_jpim1 ! vector opt. 235 ! ! horizontal density i-gradient at w-points236 z coef1 = MAX( zeps, umask(ji-1,jj,jk )+umask(ji,jj,jk ) &237 & +umask(ji-1,jj,jk-1)+umask(ji,jj,jk-1) )238 zcoef1 = 1. / ( zcoef1 * e1t (ji,jj) )239 z ai = zcoef1 * ( zgru(ji ,jj,jk ) + zgru(ji ,jj,jk-1)&240 & + zgru(ji-1,jj,jk-1) + zgru(ji-1,jj,jk ) ) * tmask (ji,jj,jk)241 ! ! horizontal density j-gradient at w-points242 zcoef2 = MAX( zeps, vmask(ji,jj-1,jk )+vmask(ji,jj,jk-1) &243 & +vmask(ji,jj-1,jk-1)+vmask(ji,jj,jk ) )244 zcoef2 = 1.0 / ( zcoef2 * e2t (ji,jj))245 zaj = zcoef2 * ( zgrv(ji,jj ,jk ) + zgrv(ji,jj ,jk-1)&246 & + zgrv(ji,jj-1,jk-1) + zgrv(ji,jj-1,jk ) )* tmask (ji,jj,jk)246 ! !* Local vertical density gradient evaluated from N^2 247 zbw = zm1_2g * pn2 (ji,jj,jk) * ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) 248 ! !* Slopes at w point 249 ! ! i- & j-gradient of density at w-points 250 zci = MAX( umask(ji-1,jj,jk ) + umask(ji,jj,jk ) & 251 & + umask(ji-1,jj,jk-1) + umask(ji,jj,jk-1) , zeps ) * e1t(ji,jj) 252 zcj = MAX( vmask(ji,jj-1,jk ) + vmask(ji,jj,jk-1) & 253 & + vmask(ji,jj-1,jk-1) + vmask(ji,jj,jk ) , zeps ) * e2t(ji,jj) 254 zai = ( zgru (ji-1,jj,jk ) + zgru (ji,jj,jk-1) & 255 & + zgru (ji-1,jj,jk-1) + zgru (ji,jj,jk ) ) / zci * tmask (ji,jj,jk) 256 zaj = ( zgrv (ji,jj-1,jk ) + zgrv (ji,jj,jk-1) & 257 & + zgrv (ji,jj-1,jk-1) + zgrv (ji,jj,jk ) ) / zcj * tmask (ji,jj,jk) 247 258 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0. 248 ! ! static instability: kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 249 zbi = MIN( zwy (ji,jj,jk),- 100.*ABS(zai), -7.e+3/fse3w(ji,jj,jk)*ABS(zai) ) 250 zbj = MIN( zwy (ji,jj,jk), -100.*ABS(zaj), -7.e+3/fse3w(ji,jj,jk)*ABS(zaj) ) 251 ! ! wslpi and wslpj output in zwz and zww, resp. 252 zalpha = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) ) 253 zcoef3 = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10. ) 254 zwz(ji,jj,jk) = ( zai / ( zbi - zeps) * ( 1. - zalpha ) & 255 & + zcoef3 * wslpiml(ji,jj) * zalpha ) * tmask (ji,jj,jk) 256 zww(ji,jj,jk) = ( zaj / ( zbj - zeps) * ( 1. - zalpha ) & 257 & + zcoef3 * wslpjml(ji,jj) * zalpha ) * tmask (ji,jj,jk) 258 END DO 259 END DO 260 END DO 261 CALL lbc_lnk( zwz, 'T', -1. ) ; CALL lbc_lnk( zww, 'T', -1. ) ! lateral boundary conditions on zwz and zww 259 ! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 260 zbi = MIN( zbw ,- 100._wp* ABS( zai ) , -7.e+3_wp/fse3w(ji,jj,jk)* ABS( zai ) ) 261 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/fse3w(ji,jj,jk)* ABS( zaj ) ) 262 ! ! wslpi and wslpj with ML flattening (output in zwz and zww, resp.) 263 zfk = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) ) ! zfk=1 in the ML otherwise zfk=0 264 zck = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10._wp ) 265 zwz(ji,jj,jk) = ( zai / ( zbi - zeps ) * ( 1._wp - zfk ) + zck * wslpiml(ji,jj) * zfk ) * tmask(ji,jj,jk) 266 zww(ji,jj,jk) = ( zaj / ( zbj - zeps ) * ( 1._wp - zfk ) + zck * wslpjml(ji,jj) * zfk ) * tmask(ji,jj,jk) 267 268 !!gm modif to suppress omlmask.... (as in Griffies operator) 269 ! ! ! jk must be >= ML level for zfk=1. otherwise zfk=0. 270 ! zfk = REAL( 1 - 1/(1 + jk / nmln(ji+1,jj)), wp ) 271 ! zck = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10. ) 272 ! zwz(ji,jj,jk) = ( zfk * zai / ( zbi - zeps ) + ( 1._wp - zfk ) * wslpiml(ji,jj) * zck ) * tmask(ji,jj,jk) 273 ! zww(ji,jj,jk) = ( zfk * zaj / ( zbj - zeps ) + ( 1._wp - zfk ) * wslpjml(ji,jj) * zck ) * tmask(ji,jj,jk) 274 !!gm end modif 275 END DO 276 END DO 277 END DO 278 CALL lbc_lnk( zwz, 'T', -1. ) ; CALL lbc_lnk( zww, 'T', -1. ) ! lateral boundary conditions 262 279 ! 263 280 ! !* horizontal Shapiro filter 264 281 DO jk = 2, jpkm1 265 DO jj = 2, jpjm1, jpj-3 ! rows jj=2 and =jpjm1282 DO jj = 2, jpjm1, MAX(1, jpj-3) ! rows jj=2 and =jpjm1 only 266 283 DO ji = 2, jpim1 267 zcofw = tmask(ji,jj,jk) / 16.268 284 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 269 285 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 270 286 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 271 287 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 272 & + 4.* zwz(ji ,jj ,jk) ) * zcofw288 & + 4.* zwz(ji ,jj ,jk) ) * z1_16 * tmask(ji,jj,jk) 273 289 274 290 wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & … … 276 292 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 277 293 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 278 & + 4.* zww(ji ,jj ,jk) ) * zcofw294 & + 4.* zww(ji ,jj ,jk) ) * z1_16 * tmask(ji,jj,jk) 279 295 END DO 280 296 END DO 281 297 DO jj = 3, jpj-2 ! other rows 282 298 DO ji = fs_2, fs_jpim1 ! vector opt. 283 zcofw = tmask(ji,jj,jk) / 16.284 299 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 285 300 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 286 301 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 287 302 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 288 & + 4.* zwz(ji ,jj ,jk) ) * zcofw303 & + 4.* zwz(ji ,jj ,jk) ) * z1_16 * tmask(ji,jj,jk) 289 304 290 305 wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & … … 292 307 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 293 308 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 294 & + 4.* zww(ji ,jj ,jk) ) * zcofw309 & + 4.* zww(ji ,jj ,jk) ) * z1_16 * tmask(ji,jj,jk) 295 310 END DO 296 311 END DO … … 298 313 DO jj = 2, jpjm1 299 314 DO ji = fs_2, fs_jpim1 ! vector opt. 300 z 1u = ( umask(ji,jj,jk) + umask(ji-1,jj,jk) ) *.5301 z1v = ( vmask(ji,jj,jk) + vmask(ji,jj-1,jk) ) *.5302 wslpi(ji,jj,jk) = wslpi(ji,jj,jk) * z 1u * z1v303 wslpj(ji,jj,jk) = wslpj(ji,jj,jk) * z 1u * z1v315 zck = ( umask(ji,jj,jk) + umask(ji-1,jj,jk) ) & 316 & * ( vmask(ji,jj,jk) + vmask(ji,jj-1,jk) ) * 0.25 317 wslpi(ji,jj,jk) = wslpi(ji,jj,jk) * zck 318 wslpj(ji,jj,jk) = wslpj(ji,jj,jk) * zck 304 319 END DO 305 320 END DO 306 321 END DO 307 322 308 309 323 ! III. Specific grid points 310 324 ! =========================== … … 313 327 ! ! Gibraltar Strait 314 328 ij0 = 50 ; ij1 = 53 315 ii0 = 69 ; ii1 = 71 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0. e0329 ii0 = 69 ; ii1 = 71 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 316 330 ij0 = 51 ; ij1 = 53 317 ii0 = 68 ; ii1 = 71 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0. e0318 ii0 = 69 ; ii1 = 71 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0. e0319 ii0 = 69 ; ii1 = 71 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0. e0331 ii0 = 68 ; ii1 = 71 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 332 ii0 = 69 ; ii1 = 71 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 333 ii0 = 69 ; ii1 = 71 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 320 334 ! 321 335 ! ! Mediterrannean Sea 322 336 ij0 = 49 ; ij1 = 56 323 ii0 = 71 ; ii1 = 90 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0. e0337 ii0 = 71 ; ii1 = 90 ; uslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 324 338 ij0 = 50 ; ij1 = 56 325 ii0 = 70 ; ii1 = 90 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0. e0326 ii0 = 71 ; ii1 = 90 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0. e0327 ii0 = 71 ; ii1 = 90 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0. e0339 ii0 = 70 ; ii1 = 90 ; vslp ( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 340 ii0 = 71 ; ii1 = 90 ; wslpi( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 341 ii0 = 71 ; ii1 = 90 ; wslpj( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , : ) = 0._wp 328 342 ENDIF 329 343 … … 341 355 ! 342 356 END SUBROUTINE ldf_slp 343 344 345 SUBROUTINE ldf_slp_mxl( prd, pn2 ) 357 358 359 SUBROUTINE ldf_slp_grif ( kt ) 360 !!---------------------------------------------------------------------- 361 !! *** ROUTINE ldf_slp_grif *** 362 !! 363 !! ** Purpose : Compute the squared slopes of neutral surfaces (slope 364 !! of iso-pycnal surfaces referenced locally) (ln_traldf_grif=T) 365 !! at W-points using the Griffies quarter-cells. 366 !! 367 !! ** Method : calculates alpha and beta at T-points 368 !! 369 !! ** Action : - triadi_g, triadj_g T-pts i- and j-slope triads relative to geopot. (used for eiv) 370 !! - triadi , triadj T-pts i- and j-slope triads relative to model-coordinate 371 !! - wslp2 squared slope of neutral surfaces at w-points. 372 !!---------------------------------------------------------------------- 373 USE oce, zdit => ua ! use ua as workspace 374 USE oce, zdis => va ! use va as workspace 375 USE oce, zdjt => ta ! use ta as workspace 376 USE oce, zdjs => sa ! use sa as workspace 377 !! 378 INTEGER, INTENT( in ) :: kt ! ocean time-step index 379 !! 380 INTEGER :: ji, jj, jk, jl, ip, jp, kp ! dummy loop indices 381 INTEGER :: iku, ikv ! temporary integer 382 REAL(wp) :: zfacti, zfactj, zatempw,zatempu,zatempv ! local scalars 383 REAL(wp) :: zbu, zbv, zbti, zbtj 384 REAL(wp) :: zdxrho_raw, zti_coord, zti_raw, zti_lim, zti_lim2, zti_g_raw, zti_g_lim 385 REAL(wp) :: zdyrho_raw, ztj_coord, ztj_raw, ztj_lim, ztj_lim2, ztj_g_raw, ztj_g_lim 386 REAL(wp) :: zdzrho_raw 387 REAL(wp), DIMENSION(jpi,jpj,jpk,0:1) :: zdzrho, zdyrho, zdxrho ! Horizontal and vertical density gradients 388 REAL(wp), DIMENSION(jpi,jpj,0:1,0:1) :: zti_mlb, ztj_mlb 389 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdkt, zdks 390 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zalpha, zbeta ! alpha, beta at T points, at depth fsgdept 391 REAL(wp), DIMENSION(jpi,jpj) :: z1_mlbw 392 !!---------------------------------------------------------------------- 393 394 !--------------------------------! 395 ! Some preliminary calculation ! 396 !--------------------------------! 397 ! 398 CALL eos_alpbet( tsb, zalpha, zbeta ) !== before thermal and haline expension coeff. at T-points ==! 399 ! 400 DO jk = 1, jpkm1 !== before lateral T & S gradients at T-level jk ==! 401 DO jj = 1, jpjm1 402 DO ji = 1, fs_jpim1 ! vector opt. 403 zdit(ji,jj,jk) = ( tb(ji+1,jj,jk) - tb(ji,jj,jk) ) * umask(ji,jj,jk) ! i-gradient of T and S at jj 404 zdis(ji,jj,jk) = ( sb(ji+1,jj,jk) - sb(ji,jj,jk) ) * umask(ji,jj,jk) 405 zdjt(ji,jj,jk) = ( tb(ji,jj+1,jk) - tb(ji,jj,jk) ) * vmask(ji,jj,jk) ! j-gradient of T and S at jj 406 zdjs(ji,jj,jk) = ( sb(ji,jj+1,jk) - sb(ji,jj,jk) ) * vmask(ji,jj,jk) 407 END DO 408 END DO 409 END DO 410 IF( ln_zps ) THEN ! partial steps: correction at the last level 411 # if defined key_vectopt_loop 412 DO jj = 1, 1 413 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 414 # else 415 DO jj = 1, jpjm1 416 DO ji = 1, jpim1 417 # endif 418 zdit(ji,jj,mbku(ji,jj)) = gtsu(ji,jj,jp_tem) ! i-gradient of T and S 419 zdis(ji,jj,mbku(ji,jj)) = gtsu(ji,jj,jp_sal) 420 zdjt(ji,jj,mbkv(ji,jj)) = gtsv(ji,jj,jp_tem) ! j-gradient of T and S 421 zdjs(ji,jj,mbkv(ji,jj)) = gtsv(ji,jj,jp_sal) 422 END DO 423 END DO 424 ENDIF 425 ! 426 zdkt(:,:,1) = 0._wp !== before vertical T & S gradient at w-level ==! 427 zdks(:,:,1) = 0._wp 428 DO jk = 2, jpk 429 zdkt(:,:,jk) = ( tb(:,:,jk-1) - tb(:,:,jk) ) * tmask(:,:,jk) 430 zdks(:,:,jk) = ( sb(:,:,jk-1) - sb(:,:,jk) ) * tmask(:,:,jk) 431 END DO 432 ! 433 ! 434 DO jl = 0, 1 !== density i-, j-, and k-gradients ==! 435 ip = jl ; jp = jl ! guaranteed nonzero gradients ( absolute value larger than repsln) 436 DO jk = 1, jpkm1 ! done each pair of triad 437 DO jj = 1, jpjm1 ! NB: not masked due to the minimum value set 438 DO ji = 1, fs_jpim1 ! vector opt. 439 zdxrho_raw = ( zalpha(ji+ip,jj ,jk) * zdit(ji,jj,jk) + zbeta(ji+ip,jj ,jk) * zdis(ji,jj,jk) ) / e1u(ji,jj) 440 zdyrho_raw = ( zalpha(ji ,jj+jp,jk) * zdjt(ji,jj,jk) + zbeta(ji ,jj+jp,jk) * zdjs(ji,jj,jk) ) / e2v(ji,jj) 441 zdxrho(ji+ip,jj ,jk,1-ip) = SIGN( MAX( repsln, ABS( zdxrho_raw ) ), zdxrho_raw ) ! keep the sign 442 zdyrho(ji ,jj+jp,jk,1-jp) = SIGN( MAX( repsln, ABS( zdyrho_raw ) ), zdyrho_raw ) 443 END DO 444 END DO 445 END DO 446 END DO 447 DO kp = 0, 1 !== density i-, j-, and k-gradients ==! 448 DO jk = 1, jpkm1 ! done each pair of triad 449 DO jj = 1, jpj ! NB: not masked due to the minimum value set 450 DO ji = 1, jpi ! vector opt. 451 zdzrho_raw = ( zalpha(ji,jj,jk) * zdkt(ji,jj,jk+kp) + zbeta(ji,jj,jk) * zdks(ji,jj,jk+kp) ) & 452 & / fse3w(ji,jj,jk+kp) 453 zdzrho(ji ,jj ,jk, kp) = - MIN( - repsln, zdzrho_raw ) ! force zdzrho >= repsln 454 END DO 455 END DO 456 END DO 457 END DO 458 ! 459 DO jj = 1, jpj !== Reciprocal depth of the w-point below ML base ==! 460 DO ji = 1, jpi 461 jk = MIN( nmln(ji,jj), mbkt(ji,jj) ) + 1 ! MIN in case ML depth is the ocean depth 462 z1_mlbw(ji,jj) = 1._wp / fsdepw(ji,jj,jk) 463 END DO 464 END DO 465 ! 466 ! !== intialisations to zero ==! 467 ! 468 wslp2 (:,:,:) = 0._wp ! wslp2 will be cumulated 3D field set to zero 469 triadi_g(:,:,1,:,:) = 0._wp ; triadi_g(:,:,jpk,:,:) = 0._wp ! set surface and bottom slope to zero 470 triadj_g(:,:,1,:,:) = 0._wp ; triadj_g(:,:,jpk,:,:) = 0._wp 471 !!gm _iso set to zero missing 472 triadi (:,:,1,:,:) = 0._wp ; triadj (:,:,jpk,:,:) = 0._wp ! set surface and bottom slope to zero 473 triadj (:,:,1,:,:) = 0._wp ; triadj (:,:,jpk,:,:) = 0._wp 474 475 !-------------------------------------! 476 ! Triads just below the Mixed Layer ! 477 !-------------------------------------! 478 ! 479 DO jl = 0, 1 ! calculate slope of the 4 triads immediately ONE level below mixed-layer base 480 DO kp = 0, 1 ! with only the slope-max limit and MASKED 481 DO jj = 1, jpjm1 482 DO ji = 1, fs_jpim1 483 ip = jl ; jp = jl 484 jk = MIN( nmln(ji+ip,jj) , mbkt(ji+ip,jj) ) + 1 ! ML level+1 (MIN in case ML depth is the ocean depth) 485 zti_g_raw = ( zdxrho(ji+ip,jj,jk-kp,1-ip) / zdzrho(ji+ip,jj,jk-kp,kp) & 486 & + ( fsdept(ji+1,jj,jk-kp) - fsdept(ji,jj,jk-kp) ) / e1u(ji,jj) ) * umask(ji,jj,jk) 487 jk = MIN( nmln(ji,jj+jp) , mbkt(ji,jj+jp) ) + 1 488 ztj_g_raw = ( zdyrho(ji,jj+jp,jk-kp,1-jp) / zdzrho(ji,jj+jp,jk-kp,kp) & 489 & + ( fsdept(ji,jj+1,jk-kp) - fsdept(ji,jj,jk-kp) ) / e2v(ji,jj) ) * vmask(ji,jj,jk) 490 zti_mlb(ji+ip,jj ,1-ip,kp) = SIGN( MIN( rn_slpmax, ABS( zti_g_raw ) ), zti_g_raw ) 491 ztj_mlb(ji ,jj+jp,1-jp,kp) = SIGN( MIN( rn_slpmax, ABS( ztj_g_raw ) ), ztj_g_raw ) 492 END DO 493 END DO 494 END DO 495 END DO 496 497 !-------------------------------------! 498 ! Triads with surface limits ! 499 !-------------------------------------! 500 ! 501 DO kp = 0, 1 ! k-index of triads 502 DO jl = 0, 1 503 ip = jl ; jp = jl ! i- and j-indices of triads (i-k and j-k planes) 504 DO jk = 1, jpkm1 505 DO jj = 1, jpjm1 506 DO ji = 1, fs_jpim1 ! vector opt. 507 ! 508 ! Calculate slope relative to geopotentials used for GM skew fluxes 509 ! For s-coordinate, subtract slope at t-points (equivalent to *adding* gradient of depth) 510 ! Limit by slope *relative to geopotentials* by rn_slpmax, and mask by psi-point 511 ! masked by umask taken at the level of dz(rho) 512 ! 513 ! raw slopes: unmasked unbounded slopes (relative to geopotential (zti_g) and model surface (zti) 514 ! 515 zti_raw = zdxrho(ji+ip,jj ,jk,1-ip) / zdzrho(ji+ip,jj ,jk,kp) ! unmasked 516 ztj_raw = zdyrho(ji ,jj+jp,jk,1-jp) / zdzrho(ji ,jj+jp,jk,kp) 517 zti_coord = ( fsdept(ji+1,jj ,jk) - fsdept(ji,jj,jk) ) / e1u(ji,jj) 518 ztj_coord = ( fsdept(ji ,jj+1,jk) - fsdept(ji,jj,jk) ) / e2v(ji,jj) 519 ! unmasked 520 zti_g_raw = zti_raw + zti_coord ! ref to geopot surfaces 521 ztj_g_raw = ztj_raw + ztj_coord 522 zti_g_lim = SIGN( MIN( rn_slpmax, ABS( zti_g_raw ) ), zti_g_raw ) 523 ztj_g_lim = SIGN( MIN( rn_slpmax, ABS( ztj_g_raw ) ), ztj_g_raw ) 524 ! 525 ! Below ML use limited zti_g as is 526 ! Inside ML replace by linearly reducing sx_mlb towards surface 527 ! 528 zfacti = REAL( 1 - 1/(1 + (jk+kp-1)/nmln(ji+ip,jj)), wp ) ! k index of uppermost point(s) of triad is jk+kp-1 529 zfactj = REAL( 1 - 1/(1 + (jk+kp-1)/nmln(ji,jj+jp)), wp ) ! must be .ge. nmln(ji,jj) for zfact=1 530 ! ! otherwise zfact=0 531 zti_g_lim = zfacti * zti_g_lim & 532 & + ( 1._wp - zfacti ) * zti_mlb(ji+ip,jj,1-ip,kp) & 533 & * fsdepw(ji+ip,jj,jk+kp) * z1_mlbw(ji+ip,jj) 534 ztj_g_lim = zfactj * ztj_g_lim & 535 & + ( 1._wp - zfactj ) * ztj_mlb(ji,jj+jp,1-jp,kp) & 536 & * fsdepw(ji,jj+jp,jk+kp) * z1_mlbw(ji,jj+jp) ! masked 537 ! 538 triadi_g(ji+ip,jj ,jk,1-ip,kp) = zti_g_lim * umask(ji,jj,jk+kp) 539 triadj_g(ji ,jj+jp,jk,1-jp,kp) = ztj_g_lim * vmask(ji,jj,jk+kp) 540 ! 541 ! Get coefficients of isoneutral diffusion tensor 542 ! 1. Utilise gradients *relative* to s-coordinate, so add t-point slopes (*subtract* depth gradients) 543 ! 2. We require that isoneutral diffusion gives no vertical buoyancy flux 544 ! i.e. 33 term = (real slope* 31, 13 terms) 545 ! To do this, retain limited sx**2 in vertical flux, but divide by real slope for 13/31 terms 546 ! Equivalent to tapering A_iso = sx_limited**2/(real slope)**2 547 ! 548 zti_lim = zti_g_lim - zti_coord ! remove the coordinate slope ==> relative to coordinate surfaces 549 ztj_lim = ztj_g_lim - ztj_coord 550 zti_lim2 = zti_lim * zti_lim * umask(ji,jj,jk+kp) ! square of limited slopes ! masked <<== 551 ztj_lim2 = ztj_lim * ztj_lim * vmask(ji,jj,jk+kp) 552 ! 553 zbu = e1u(ji ,jj) * e2u(ji ,jj) * fse3u(ji ,jj,jk ) 554 zbv = e1v(ji ,jj) * e2v(ji ,jj) * fse3v(ji ,jj,jk ) 555 zbti = e1t(ji+ip,jj) * e2t(ji+ip,jj) * fse3w(ji+ip,jj,jk+kp) 556 zbtj = e1t(ji,jj+jp) * e2t(ji,jj+jp) * fse3w(ji,jj+jp,jk+kp) 557 ! 558 triadi(ji+ip,jj ,jk,1-ip,kp) = zti_lim2 / zti_raw ! masked 559 triadj(ji ,jj+jp,jk,1-jp,kp) = ztj_lim2 / ztj_raw 560 ! 561 !!gm this may inhibit vectorization on Vect Computers, and even on scalar computers.... ==> to be checked 562 wslp2 (ji+ip,jj,jk+kp) = wslp2(ji+ip,jj,jk+kp) + 0.25_wp * zbu / zbti * zti_lim2 ! masked 563 wslp2 (ji,jj+jp,jk+kp) = wslp2(ji,jj+jp,jk+kp) + 0.25_wp * zbv / zbtj * ztj_lim2 564 END DO 565 END DO 566 END DO 567 END DO 568 END DO 569 ! 570 wslp2(:,:,1) = 0._wp ! force the surface wslp to zero 571 572 CALL lbc_lnk( wslp2, 'W', 1. ) ! lateral boundary confition on wslp2 only ==>>> gm : necessary ? to be checked 573 ! 574 END SUBROUTINE ldf_slp_grif 575 576 577 SUBROUTINE ldf_slp_mxl( prd, pn2, p_gru, p_grv, p_dzr ) 346 578 !!---------------------------------------------------------------------- 347 579 !! *** ROUTINE ldf_slp_mxl *** … … 350 582 !! the mixed layer. 351 583 !! 352 !! ** Method : caution: use zgru, zgrv and zwy computed in ldf_slp 584 !! ** Method : The slope in the i-direction is computed at u- & w-points 585 !! (uslpml, wslpiml) and the slope in the j-direction is computed 586 !! at v- and w-points (vslpml, wslpjml) with the same bounds as 587 !! in ldf_slp. 353 588 !! 354 589 !! ** Action : uslpml, wslpiml : i- & j-slopes of neutral surfaces … … 356 591 !! omlmask : mixed layer mask 357 592 !!---------------------------------------------------------------------- 358 USE oce , zgru => ua ! zgru, zgrv (ua, va used as workspace) 359 USE oce , zgrv => va ! set to i- & j-density gradient in ldf_slp 360 USE oce , zwy => ta ! set to vertical density gradient at T-point in ldfslp 361 !! 362 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: prd ! in situ density 363 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: pn2 ! Brunt-Vaisala frequency (locally ref.) 364 !! 365 INTEGER :: ji, jj, jk ! dummy loop indices 366 INTEGER :: ik, ikm1 ! temporary integers 367 REAL(wp) :: zeps, zmg, zm05g ! temporary scalars 368 REAL(wp) :: zcoef1, zcoef2 ! - - 369 REAL(wp) :: zau, zbu, zai, zbi ! - - 370 REAL(wp) :: zav, zbv, zaj, zbj ! - - 371 REAL(wp) :: zbw ! - - 372 !!---------------------------------------------------------------------- 373 374 zeps = 1.e-20 ! Local constant initialization 375 zmg = -1.0 / grav 376 zm05g = -0.5 / grav 377 ! 378 uslpml (1,:) = 0.e0 ; uslpml (jpi,:) = 0.e0 379 vslpml (1,:) = 0.e0 ; vslpml (jpi,:) = 0.e0 380 wslpiml(1,:) = 0.e0 ; wslpiml(jpi,:) = 0.e0 381 wslpjml(1,:) = 0.e0 ; wslpjml(jpi,:) = 0.e0 382 383 ! ! surface mixed layer mask 384 DO jk = 1, jpk ! =1 inside the mixed layer, =0 otherwise 593 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: prd ! in situ density 594 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: pn2 ! Brunt-Vaisala frequency (locally ref.) 595 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: p_gru, p_grv ! i- & j-gradient of density (u- & v-pts) 596 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: p_dzr ! z-gradient of density (T-point) 597 !! 598 INTEGER :: ji , jj , jk ! dummy loop indices 599 INTEGER :: iku, ikv, ik, ikm1 ! local integers 600 REAL(wp) :: zeps, zm1_g, zm1_2g ! local scalars 601 REAL(wp) :: zci, zfi, zau, zbu, zai, zbi ! - - 602 REAL(wp) :: zcj, zfj, zav, zbv, zaj, zbj ! - - 603 REAL(wp) :: zck, zfk, zbw ! - - 604 !!---------------------------------------------------------------------- 605 606 zeps = 1.e-20_wp !== Local constant initialization ==! 607 zm1_g = -1.0_wp / grav 608 zm1_2g = -0.5_wp / grav 609 ! 610 uslpml (1,:) = 0._wp ; uslpml (jpi,:) = 0._wp 611 vslpml (1,:) = 0._wp ; vslpml (jpi,:) = 0._wp 612 wslpiml(1,:) = 0._wp ; wslpiml(jpi,:) = 0._wp 613 wslpjml(1,:) = 0._wp ; wslpjml(jpi,:) = 0._wp 614 ! 615 ! !== surface mixed layer mask ! 616 DO jk = 1, jpk ! =1 inside the mixed layer, =0 otherwise 385 617 # if defined key_vectopt_loop 386 618 DO jj = 1, 1 … … 391 623 # endif 392 624 ik = nmln(ji,jj) - 1 393 IF( jk <= ik ) THEN ; omlmask(ji,jj,jk) = 1. e0394 ELSE ; omlmask(ji,jj,jk) = 0. e0625 IF( jk <= ik ) THEN ; omlmask(ji,jj,jk) = 1._wp 626 ELSE ; omlmask(ji,jj,jk) = 0._wp 395 627 ENDIF 396 628 END DO … … 409 641 !----------------------------------------------------------------------- 410 642 ! 411 ! !* Slope at u points412 643 # if defined key_vectopt_loop 413 644 DO jj = 1, 1 … … 417 648 DO ji = 2, jpim1 418 649 # endif 419 ! horizontal and vertical density gradient at u-points 420 ik = MAX( 1, nmln(ji,jj) , nmln(ji+1,jj) ) 421 ik = MIN( ik, jpkm1 ) 422 zau = 1./ e1u(ji,jj) * zgru(ji,jj,ik) 423 zbu = 0.5*( zwy(ji,jj,ik) + zwy(ji+1,jj,ik) ) 424 ! bound the slopes: abs(zw.)<= 1/100 and zb..<0 425 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 426 zbu = MIN( zbu, -100.*ABS(zau), -7.e+3/fse3u(ji,jj,ik)*ABS(zau) ) 427 ! uslpml 428 uslpml (ji,jj) = zau / ( zbu - zeps ) * umask (ji,jj,ik) 429 END DO 430 END DO 431 CALL lbc_lnk( uslpml, 'U', -1. ) ! lateral boundary conditions (i-gradient => sign change) 432 433 ! !* Slope at v points 434 # if defined key_vectopt_loop 435 DO jj = 1, 1 436 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) 437 # else 438 DO jj = 2, jpjm1 439 DO ji = 2, jpim1 440 # endif 441 ! horizontal and vertical density gradient at v-points 442 ik = MAX( 1, nmln(ji,jj) , nmln(ji,jj+1) ) 443 ik = MIN( ik,jpkm1 ) 444 zav = 1./ e2v(ji,jj) * zgrv(ji,jj,ik) 445 zbv = 0.5*( zwy(ji,jj,ik) + zwy(ji,jj+1,ik) ) 446 ! bound the slopes: abs(zw.)<= 1/100 and zb..<0 447 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 448 zbv = MIN( zbv, -100.*ABS(zav), -7.e+3/fse3v(ji,jj,ik)*ABS( zav ) ) 449 ! vslpml 450 vslpml (ji,jj) = zav / ( zbv - zeps ) * vmask (ji,jj,ik) 451 END DO 452 END DO 453 CALL lbc_lnk( vslpml, 'V', -1. ) ! lateral boundary conditions (j-gradient => sign change) 454 455 ! !* Slopes at w points 456 # if defined key_vectopt_loop 457 DO jj = 1, 1 458 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) 459 # else 460 DO jj = 2, jpjm1 461 DO ji = 2, jpim1 462 # endif 463 ik = nmln(ji,jj) + 1 464 ik = MIN( ik, jpk ) 465 ikm1 = MAX ( 1, ik-1 ) 466 ! horizontal density i-gradient at w-points 467 zcoef1 = MAX( zeps, umask(ji-1,jj,ik )+umask(ji,jj,ik ) & 468 & +umask(ji-1,jj,ikm1)+umask(ji,jj,ikm1) ) 469 zcoef1 = 1. / ( zcoef1 * e1t (ji,jj) ) 470 zai = zcoef1 * ( zgru(ji ,jj,ik ) + zgru(ji ,jj,ikm1) & 471 & + zgru(ji-1,jj,ikm1) + zgru(ji-1,jj,ik ) ) * tmask (ji,jj,ik) 472 ! horizontal density j-gradient at w-points 473 zcoef2 = MAX( zeps, vmask(ji,jj-1,ik )+vmask(ji,jj,ikm1) & 474 & +vmask(ji,jj-1,ikm1)+vmask(ji,jj,ik ) ) 475 zcoef2 = 1.0 / ( zcoef2 * e2t (ji,jj) ) 476 zaj = zcoef2 * ( zgrv(ji,jj ,ik ) + zgrv(ji,jj ,ikm1) & 477 & + zgrv(ji,jj-1,ikm1) + zgrv(ji,jj-1,ik ) ) * tmask (ji,jj,ik) 478 ! vertical density gradient at w-point (from N^2) 479 zbw = zm05g * pn2 (ji,jj,ik) * ( prd (ji,jj,ik) + prd (ji,jj,ikm1) + 2. ) 480 ! bound the slopes: abs(zw.)<= 1/100 and zb..<0. 481 ! static instability: 482 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 483 zbi = MIN ( zbw,- 100.*ABS(zai), -7.e+3/fse3w(ji,jj,ik)*ABS(zai) ) 484 zbj = MIN ( zbw, -100.*ABS(zaj), -7.e+3/fse3w(ji,jj,ik)*ABS(zaj) ) 485 ! wslpiml and wslpjml 486 wslpiml (ji,jj) = zai / ( zbi - zeps) * tmask (ji,jj,ik) 487 wslpjml (ji,jj) = zaj / ( zbj - zeps) * tmask (ji,jj,ik) 488 END DO 489 END DO 490 CALL lbc_lnk( wslpiml, 'W', -1. ) ; CALL lbc_lnk( wslpjml, 'W', -1. ) ! lateral boundary conditions 650 ! !== Slope at u- & v-points just below the Mixed Layer ==! 651 ! 652 ! !- vertical density gradient for u- and v-slopes (from dzr at T-point) 653 iku = MIN( MAX( 1, nmln(ji,jj) , nmln(ji+1,jj) ) , jpkm1 ) ! ML (MAX of T-pts, bound by jpkm1) 654 ikv = MIN( MAX( 1, nmln(ji,jj) , nmln(ji,jj+1) ) , jpkm1 ) ! 655 zbu = 0.5_wp * ( p_dzr(ji,jj,iku) + p_dzr(ji+1,jj ,iku) ) 656 zbv = 0.5_wp * ( p_dzr(ji,jj,ikv) + p_dzr(ji ,jj+1,ikv) ) 657 ! !- horizontal density gradient at u- & v-points 658 zau = p_gru(ji,jj,iku) / e1u(ji,jj) 659 zav = p_grv(ji,jj,ikv) / e2v(ji,jj) 660 ! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0 661 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 662 zbu = MIN( zbu , -100._wp* ABS( zau ) , -7.e+3_wp/fse3u(ji,jj,ik)* ABS( zau ) ) 663 zbv = MIN( zbv , -100._wp* ABS( zav ) , -7.e+3_wp/fse3v(ji,jj,ik)* ABS( zav ) ) 664 ! !- Slope at u- & v-points (uslpml, vslpml) 665 uslpml(ji,jj) = zau / ( zbu - zeps ) * umask(ji,jj,ik) 666 vslpml(ji,jj) = zav / ( zbv - zeps ) * vmask(ji,jj,ik) 667 ! 668 ! !== i- & j-slopes at w-points just below the Mixed Layer ==! 669 ! 670 ik = MIN( nmln(ji,jj) + 1, jpk ) 671 ikm1 = MAX( 1, ik-1 ) 672 ! !- vertical density gradient for w-slope (from N^2) 673 zbw = zm1_2g * pn2 (ji,jj,ik) * ( prd (ji,jj,ik) + prd (ji,jj,ikm1) + 2. ) 674 ! !- horizontal density i- & j-gradient at w-points 675 zci = MAX( umask(ji-1,jj,ik ) + umask(ji,jj,ik ) & 676 & + umask(ji-1,jj,ikm1) + umask(ji,jj,ikm1) , zeps ) * e1t(ji,jj) 677 zcj = MAX( vmask(ji,jj-1,ik ) + vmask(ji,jj,ik ) & 678 & + vmask(ji,jj-1,ikm1) + vmask(ji,jj,ikm1) , zeps ) * e2t(ji,jj) 679 zai = ( p_gru(ji-1,jj,ik ) + p_gru(ji,jj,ik) & 680 & + p_gru(ji-1,jj,ikm1) + p_gru(ji,jj,ikm1 ) ) / zci * tmask(ji,jj,ik) 681 zaj = ( p_grv(ji,jj-1,ik ) + p_grv(ji,jj,ik ) & 682 & + p_grv(ji,jj-1,ikm1) + p_grv(ji,jj,ikm1) ) / zcj * tmask(ji,jj,ik) 683 ! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0. 684 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 685 zbi = MIN( zbw , -100._wp* ABS( zai ) , -7.e+3_wp/fse3w(ji,jj,ik)* ABS( zai ) ) 686 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/fse3w(ji,jj,ik)* ABS( zaj ) ) 687 ! !- i- & j-slope at w-points (wslpiml, wslpjml) 688 wslpiml(ji,jj) = zai / ( zbi - zeps ) * tmask (ji,jj,ik) 689 wslpjml(ji,jj) = zaj / ( zbj - zeps ) * tmask (ji,jj,ik) 690 END DO 691 END DO 692 !!gm this lbc_lnk should be useless.... 693 CALL lbc_lnk( uslpml , 'U', -1. ) ; CALL lbc_lnk( vslpml , 'V', -1. ) ! lateral boundary cond. (sign change) 694 CALL lbc_lnk( wslpiml, 'W', -1. ) ; CALL lbc_lnk( wslpjml, 'W', -1. ) ! lateral boundary conditions 491 695 ! 492 696 END SUBROUTINE ldf_slp_mxl … … 503 707 !!---------------------------------------------------------------------- 504 708 INTEGER :: ji, jj, jk ! dummy loop indices 709 INTEGER :: ierr ! local integer 505 710 !!---------------------------------------------------------------------- 506 711 507 712 IF(lwp) THEN 508 713 WRITE(numout,*) 509 WRITE(numout,*) 'ldf_slp : direction of lateral mixing'510 WRITE(numout,*) '~~~~~~~ '714 WRITE(numout,*) 'ldf_slp_init : direction of lateral mixing' 715 WRITE(numout,*) '~~~~~~~~~~~~' 511 716 ENDIF 512 513 ! Direction of lateral diffusion (tracers and/or momentum) 514 ! ------------------------------ 515 ! set the slope to zero (even in s-coordinates) 516 517 uslp (:,:,:) = 0.e0 518 vslp (:,:,:) = 0.e0 519 wslpi(:,:,:) = 0.e0 520 wslpj(:,:,:) = 0.e0 521 522 uslpml (:,:) = 0.e0 523 vslpml (:,:) = 0.e0 524 wslpiml(:,:) = 0.e0 525 wslpjml(:,:) = 0.e0 526 527 IF( (ln_traldf_hor .OR. ln_dynldf_hor) .AND. .NOT. (lk_vvl .AND. ln_rstart) ) THEN 528 IF(lwp) THEN 529 WRITE(numout,*) ' Horizontal mixing in s-coordinate: slope = slope of s-surfaces' 717 718 IF( ln_traldf_grif ) THEN ! Griffies operator : triad of slopes 719 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 ) 720 ALLOCATE( triadi (jpi,jpj,jpk,0:1,0:1) , triadj (jpi,jpj,jpk,0:1,0:1) , STAT=ierr ) 721 IF( ierr > 0 ) THEN 722 CALL ctl_stop( 'ldf_slp_init : unable to allocate Griffies operator slope ' ) ; RETURN 530 723 ENDIF 531 532 ! geopotential diffusion in s-coordinates on tracers and/or momentum 533 ! The slopes of s-surfaces are computed once (no call to ldfslp in step) 534 ! The slopes for momentum diffusion are i- or j- averaged of those on tracers 535 536 ! set the slope of diffusion to the slope of s-surfaces 537 ! ( c a u t i o n : minus sign as fsdep has positive value ) 538 DO jk = 1, jpk 539 DO jj = 2, jpjm1 540 DO ji = fs_2, fs_jpim1 ! vector opt. 541 uslp (ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept(ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * umask(ji,jj,jk) 542 vslp (ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept(ji,jj+1,jk) - fsdept(ji ,jj ,jk) ) * vmask(ji,jj,jk) 543 wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw(ji+1,jj,jk) - fsdepw(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5 544 wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw(ji,jj+1,jk) - fsdepw(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5 724 ! 725 IF( ln_dynldf_iso ) CALL ctl_stop( 'ldf_slp_init: Griffies operator on momentum not supported' ) 726 ! 727 IF( ( ln_traldf_hor .AND. ln_dynldf_hor ) .AND. ln_sco ) & 728 & CALL ctl_stop( 'ldf_slp_init: horizontal Griffies operator ', & 729 & 'in s-coordinate not supported' ) 730 ! 731 ELSE ! Madec operator : slopes at u-, v-, and w-points 732 ALLOCATE( uslp(jpi,jpj,jpk) , vslp(jpi,jpj,jpk) , wslpi(jpi,jpj,jpk) , wslpj(jpi,jpj,jpk) , & 733 & omlmask(jpi,jpj,jpk) , uslpml(jpi,jpj) , vslpml(jpi,jpj) , wslpiml(jpi,jpj) , wslpjml(jpi,jpj) , STAT=ierr ) 734 IF( ierr > 0 ) THEN 735 CALL ctl_stop( 'ldf_slp_init : unable to allocate Madec operator slope ' ) ; RETURN 736 ENDIF 737 738 ! Direction of lateral diffusion (tracers and/or momentum) 739 ! ------------------------------ 740 uslp (:,:,:) = 0._wp ; uslpml (:,:) = 0._wp ! set the slope to zero (even in s-coordinates) 741 vslp (:,:,:) = 0._wp ; vslpml (:,:) = 0._wp 742 wslpi(:,:,:) = 0._wp ; wslpiml(:,:) = 0._wp 743 wslpj(:,:,:) = 0._wp ; wslpjml(:,:) = 0._wp 744 745 !!gm I no longer understand this..... 746 IF( (ln_traldf_hor .OR. ln_dynldf_hor) .AND. .NOT. (lk_vvl .AND. ln_rstart) ) THEN 747 IF(lwp) THEN 748 WRITE(numout,*) ' Horizontal mixing in s-coordinate: slope = slope of s-surfaces' 749 ENDIF 750 751 ! geopotential diffusion in s-coordinates on tracers and/or momentum 752 ! The slopes of s-surfaces are computed once (no call to ldfslp in step) 753 ! The slopes for momentum diffusion are i- or j- averaged of those on tracers 754 755 ! set the slope of diffusion to the slope of s-surfaces 756 ! ( c a u t i o n : minus sign as fsdep has positive value ) 757 DO jk = 1, jpk 758 DO jj = 2, jpjm1 759 DO ji = fs_2, fs_jpim1 ! vector opt. 760 uslp (ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept(ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * umask(ji,jj,jk) 761 vslp (ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept(ji,jj+1,jk) - fsdept(ji ,jj ,jk) ) * vmask(ji,jj,jk) 762 wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw(ji+1,jj,jk) - fsdepw(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5 763 wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw(ji,jj+1,jk) - fsdepw(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5 764 END DO 545 765 END DO 546 766 END DO 547 END DO 548 ! Lateral boundary conditions on the slopes 549 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) 550 CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 551 ENDIF 552 ! 767 ! Lateral boundary conditions on the slopes 768 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) 769 CALL lbc_lnk( wslpi, 'W', -1. ) ; CALL lbc_lnk( wslpj, 'W', -1. ) 770 ENDIF 771 ENDIF ! 553 772 END SUBROUTINE ldf_slp_init 554 773 … … 564 783 WRITE(*,*) 'ldf_slp: You should not have seen this print! error?', kt, prd(1,1,1), pn2(1,1,1) 565 784 END SUBROUTINE ldf_slp 785 SUBROUTINE ldf_slp_init ! Dummy routine 786 END SUBROUTINE ldf_slp_init 566 787 #endif 567 788 -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldftra.F90
- Property svn:eol-style deleted
- Property svn:executable deleted
r1601 r2528 34 34 # include "vectopt_loop_substitute.h90" 35 35 !!---------------------------------------------------------------------- 36 !! NEMO/OPA 3. 2 , LOCEAN-IPSL (2009)36 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 37 37 !! $Id$ 38 !! Software governed by the CeCILL licence ( modipsl/doc/NEMO_CeCILL.txt)38 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 39 39 !!---------------------------------------------------------------------- 40 40 … … 67 67 NAMELIST/namtra_ldf/ ln_traldf_lap , ln_traldf_bilap, & 68 68 & ln_traldf_level, ln_traldf_hor , ln_traldf_iso, & 69 & rn_aht_0 , rn_ahtb_0 , rn_aeiv_0 69 & ln_traldf_grif , ln_traldf_gdia, & 70 & rn_aht_0 , rn_ahtb_0 , rn_aeiv_0, & 71 & rn_slpmax 70 72 !!---------------------------------------------------------------------- 71 73 … … 80 82 WRITE(numout,*) 'ldf_tra_init : lateral tracer physics' 81 83 WRITE(numout,*) '~~~~~~~~~~~~ ' 82 WRITE(numout,*) ' Namelist namtra_ldf : lateral mixing coefficients'84 WRITE(numout,*) ' Namelist namtra_ldf : lateral mixing parameters (type, direction, coefficients)' 83 85 WRITE(numout,*) ' laplacian operator ln_traldf_lap = ', ln_traldf_lap 84 86 WRITE(numout,*) ' bilaplacian operator ln_traldf_bilap = ', ln_traldf_bilap 87 WRITE(numout,*) ' iso-level ln_traldf_level = ', ln_traldf_level 88 WRITE(numout,*) ' horizontal (geopotential) ln_traldf_hor = ', ln_traldf_hor 89 WRITE(numout,*) ' iso-neutral ln_traldf_iso = ', ln_traldf_iso 90 WRITE(numout,*) ' iso-neutral (Griffies) ln_traldf_grif = ', ln_traldf_grif 91 WRITE(numout,*) ' Griffies strmfn diagnostics ln_traldf_gdia = ', ln_traldf_gdia 85 92 WRITE(numout,*) ' lateral eddy diffusivity rn_aht_0 = ', rn_aht_0 86 93 WRITE(numout,*) ' background hor. diffusivity rn_ahtb_0 = ', rn_ahtb_0 87 94 WRITE(numout,*) ' eddy induced velocity coef. rn_aeiv_0 = ', rn_aeiv_0 95 WRITE(numout,*) ' maximum isoppycnal slope rn_slpmax = ', rn_slpmax 96 WRITE(numout,*) ' + griffies operator internal controls not set via the namelist (experimental): ' 97 WRITE(numout,*) ' calculate triads twice l_triad_iso = ', l_triad_iso 98 WRITE(numout,*) ' no Shapiro filter l_no_smooth = ', l_no_smooth 88 99 WRITE(numout,*) 89 100 ENDIF -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldftra_c1d.h90
- Property svn:eol-style deleted
r1152 r2528 4 4 5 5 !!---------------------------------------------------------------------- 6 !! OPA 9.0 , LOCEAN-IPSL (2005)6 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 7 7 !! $Id$ 8 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt8 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 9 9 !!---------------------------------------------------------------------- 10 10 -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldftra_c2d.h90
- Property svn:eol-style deleted
r1152 r2528 4 4 5 5 !!---------------------------------------------------------------------- 6 !! OPA 9.0 , LOCEAN-IPSL (2005)6 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 7 7 !! $Id$ 8 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt8 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 9 9 !!---------------------------------------------------------------------- 10 10 -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldftra_c3d.h90
- Property svn:eol-style deleted
r1566 r2528 4 4 5 5 !!---------------------------------------------------------------------- 6 !! OPA 9.0 , LOCEAN-IPSL (2005)6 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 7 7 !! $Id$ 8 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt8 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 9 9 !!---------------------------------------------------------------------- 10 10 -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldftra_oce.F90
- Property svn:eol-style deleted
r1601 r2528 20 20 LOGICAL , PUBLIC :: ln_traldf_hor = .FALSE. !: horizontal (geopotential) direction 21 21 LOGICAL , PUBLIC :: ln_traldf_iso = .TRUE. !: iso-neutral direction 22 LOGICAL , PUBLIC :: ln_traldf_grif = .FALSE. !: griffies skew flux 23 LOGICAL , PUBLIC :: ln_traldf_gdia = .FALSE. !: griffies skew flux streamfunction diagnostics 22 24 REAL(wp), PUBLIC :: rn_aht_0 = 2000._wp !: lateral eddy diffusivity (m2/s) 23 25 REAL(wp), PUBLIC :: rn_ahtb_0 = 0._wp !: lateral background eddy diffusivity (m2/s) 24 26 REAL(wp), PUBLIC :: rn_aeiv_0 = 2000._wp !: eddy induced velocity coefficient (m2/s) 27 REAL(wp), PUBLIC :: rn_slpmax = 0.01_wp !: slope limit 25 28 26 29 REAL(wp), PUBLIC :: aht0, ahtb0, aeiv0 !!: OLD namelist names 30 LOGICAL , PUBLIC :: l_triad_iso = .FALSE. !: calculate triads twice 31 LOGICAL , PUBLIC :: l_no_smooth = .FALSE. !: no Shapiro smoothing 27 32 28 33 #if defined key_traldf_c3d … … 41 46 !! 'key_traldf_eiv' eddy induced velocity 42 47 !!---------------------------------------------------------------------- 43 LOGICAL, PUBLIC, PARAMETER :: lk_traldf_eiv = .TRUE. !: eddy induced velocity flag48 LOGICAL, PUBLIC, PARAMETER :: lk_traldf_eiv = .TRUE. !: eddy induced velocity flag 44 49 45 50 # if defined key_traldf_c3d … … 65 70 66 71 !!---------------------------------------------------------------------- 67 !! NEMO/OPA 3. 2 , LOCEAN-IPSL (2009)72 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 68 73 !! $Id$ 69 !! Software governed by the CeCILL licence ( modipsl/doc/NEMO_CeCILL.txt)74 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 70 75 !!===================================================================== 71 76 END MODULE ldftra_oce -
trunk/NEMOGCM/NEMO/OPA_SRC/LDF/ldftra_substitute.h90
- Property svn:eol-style deleted
r1152 r2528 6 6 !!---------------------------------------------------------------------- 7 7 !!---------------------------------------------------------------------- 8 !! OPA 9.0 , LOCEAN-IPSL (2005)8 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 9 9 !! $Id$ 10 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt10 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 11 11 !!---------------------------------------------------------------------- 12 12 #if defined key_traldf_c3d
Note: See TracChangeset
for help on using the changeset viewer.
