Changeset 2980
- Timestamp:
- 2011-10-24T11:29:46+02:00 (12 years ago)
- Location:
- branches/2011/dev_NOC_2011_MERGE
- Files:
-
- 12 edited
-
DOC/TexFiles/Namelist/namtra_ldf (modified) (1 diff)
-
NEMOGCM/CONFIG/GYRE/EXP00/namelist (modified) (2 diffs)
-
NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist (modified) (2 diffs)
-
NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist (modified) (2 diffs)
-
NEMOGCM/CONFIG/POMME/EXP00/namelist (modified) (2 diffs)
-
NEMOGCM/NEMO/OPA_SRC/LDF/ldfslp.F90 (modified) (35 diffs)
-
NEMOGCM/NEMO/OPA_SRC/LDF/ldftra.F90 (modified) (2 diffs)
-
NEMOGCM/NEMO/OPA_SRC/LDF/ldftra_oce.F90 (modified) (1 diff)
-
NEMOGCM/NEMO/OPA_SRC/TRA/eosbn2.F90 (modified) (21 diffs)
-
NEMOGCM/NEMO/OPA_SRC/TRA/traldf_iso_grif.F90 (modified) (20 diffs)
-
NEMOGCM/NEMO/TOP_SRC/TRP/trcadv.F90 (modified) (1 diff)
-
NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90 (modified) (12 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/dev_NOC_2011_MERGE/DOC/TexFiles/Namelist/namtra_ldf
r2540 r2980 9 9 ln_traldf_hor = .false. ! horizontal (geopotential) (require "key_ldfslp" when ln_sco=T) 10 10 ln_traldf_iso = .true. ! iso-neutral (require "key_ldfslp") 11 ln_traldf_grif = .false. ! griffies skew flux formulation (require "key_ldfslp") ! UNDER TEST, DO NOT USE 12 ln_traldf_gdia = .false. ! griffies operator strfn diagnostics (require "key_ldfslp") ! UNDER TEST, DO NOT USE 11 ln_traldf_grif = .false. ! griffies skew flux formulation (require "key_ldfslp") 12 ln_traldf_gdia = .false. ! griffies operator strfn diagnostics (require "key_ldfslp") 13 ln_triad_iso = .false. ! griffies operator calculates triads twice => pure lateral mixing in ML (require "key_ldfslp") 14 ln_botmix_grif = .false. ! griffies operator with lateral mixing on bottom (require "key_ldfslp") 13 15 ! ! Coefficient 14 16 rn_aht_0 = 2000. ! horizontal eddy diffusivity for tracers [m2/s] -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/CONFIG/GYRE/EXP00/namelist
r2715 r2980 463 463 ln_traadv_muscl2 = .false. ! MUSCL2 scheme + cen2 at boundaries 464 464 ln_traadv_ubs = .false. ! UBS scheme 465 ln_traadv_qck = .false. ! QU CIKEST scheme465 ln_traadv_qck = .false. ! QUICKEST scheme 466 466 / 467 467 !----------------------------------------------------------------------- … … 475 475 ln_traldf_hor = .false. ! horizontal (geopotential) (require "key_ldfslp" when ln_sco=T) 476 476 ln_traldf_iso = .true. ! iso-neutral (require "key_ldfslp") 477 ln_traldf_grif = .false. ! griffies skew flux formulation (require "key_ldfslp") ! UNDER TEST, DO NOT USE 478 ln_traldf_gdia = .false. ! griffies operator strfn diagnostics (require "key_ldfslp") ! UNDER TEST, DO NOT USE 477 ln_traldf_grif = .false. ! griffies skew flux formulation (require "key_ldfslp") 478 ln_traldf_gdia = .false. ! griffies operator strfn diagnostics (require "key_ldfslp") 479 ln_triad_iso = .false. ! griffies operator calculates triads twice => pure lateral mixing in ML (require "key_ldfslp") 480 ln_botmix_grif = .false. ! griffies operator with lateral mixing on bottom (require "key_ldfslp") 479 481 ! ! Coefficient 480 482 rn_aht_0 = 1000. ! horizontal eddy diffusivity for tracers [m2/s] -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist
r2715 r2980 463 463 ln_traadv_muscl2 = .false. ! MUSCL2 scheme + cen2 at boundaries 464 464 ln_traadv_ubs = .false. ! UBS scheme 465 ln_traadv_qck = .false. ! QU CIKEST scheme465 ln_traadv_qck = .false. ! QUICKEST scheme 466 466 / 467 467 !----------------------------------------------------------------------- … … 475 475 ln_traldf_hor = .false. ! horizontal (geopotential) (require "key_ldfslp" when ln_sco=T) 476 476 ln_traldf_iso = .true. ! iso-neutral (require "key_ldfslp") 477 ln_traldf_grif = .false. ! griffies skew flux formulation (require "key_ldfslp") ! UNDER TEST, DO NOT USE 478 ln_traldf_gdia = .false. ! griffies operator strfn diagnostics (require "key_ldfslp") ! UNDER TEST, DO NOT USE 479 ! ! Coefficient 477 ln_traldf_grif = .false. ! griffies skew flux formulation (require "key_ldfslp") 478 ln_traldf_gdia = .false. ! griffies operator strfn diagnostics (require "key_ldfslp") 479 ln_triad_iso = .false. ! griffies operator calculates triads twice => pure lateral mixing in ML (require "key_ldfslp") 480 ln_botmix_grif = .false. ! griffies operator with lateral mixing on bottom (require "key_ldfslp") 481 ! Coefficient 480 482 rn_aht_0 = 2000. ! horizontal eddy diffusivity for tracers [m2/s] 481 483 rn_ahtb_0 = 0. ! background eddy diffusivity for ldf_iso [m2/s] -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist
r2715 r2980 463 463 ln_traadv_muscl2 = .false. ! MUSCL2 scheme + cen2 at boundaries 464 464 ln_traadv_ubs = .false. ! UBS scheme 465 ln_traadv_qck = .false. ! QU CIKEST scheme465 ln_traadv_qck = .false. ! QUICKEST scheme 466 466 / 467 467 !----------------------------------------------------------------------- … … 475 475 ln_traldf_hor = .false. ! horizontal (geopotential) (require "key_ldfslp" when ln_sco=T) 476 476 ln_traldf_iso = .true. ! iso-neutral (require "key_ldfslp") 477 ln_traldf_grif = .false. ! griffies skew flux formulation (require "key_ldfslp") ! UNDER TEST, DO NOT USE 478 ln_traldf_gdia = .false. ! griffies operator strfn diagnostics (require "key_ldfslp") ! UNDER TEST, DO NOT USE 479 ! ! Coefficient 477 ln_traldf_grif = .false. ! griffies skew flux formulation (require "key_ldfslp") 478 ln_traldf_gdia = .false. ! griffies operator strfn diagnostics (require "key_ldfslp") 479 ln_triad_iso = .false. ! griffies operator calculates triads twice => pure lateral mixing in ML (require "key_ldfslp") 480 ln_botmix_grif = .false. ! griffies operator with lateral mixing on bottom (require "key_ldfslp") 481 ! Coefficient 480 482 rn_aht_0 = 2000. ! horizontal eddy diffusivity for tracers [m2/s] 481 483 rn_ahtb_0 = 0. ! background eddy diffusivity for ldf_iso [m2/s] -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/CONFIG/POMME/EXP00/namelist
r2650 r2980 463 463 ln_traadv_muscl2 = .false. ! MUSCL2 scheme + cen2 at boundaries 464 464 ln_traadv_ubs = .false. ! UBS scheme 465 ln_traadv_qck = .false. ! QU CIKEST scheme465 ln_traadv_qck = .false. ! QUICKEST scheme 466 466 / 467 467 !----------------------------------------------------------------------- … … 475 475 ln_traldf_hor = .false. ! horizontal (geopotential) (require "key_ldfslp" when ln_sco=T) 476 476 ln_traldf_iso = .true. ! iso-neutral (require "key_ldfslp") 477 ln_traldf_grif = .false. ! griffies skew flux formulation (require "key_ldfslp") ! UNDER TEST, DO NOT USE 478 ln_traldf_gdia = .false. ! griffies operator strfn diagnostics (require "key_ldfslp") ! UNDER TEST, DO NOT USE 477 ln_traldf_grif = .false. ! griffies skew flux formulation (require "key_ldfslp") 478 ln_traldf_gdia = .false. ! griffies operator strfn diagnostics (require "key_ldfslp") 479 ln_triad_iso = .false. ! griffies operator calculates triads twice => pure lateral mixing in ML (require "key_ldfslp") 480 ln_botmix_grif = .false. ! griffies operator with lateral mixing on bottom (require "key_ldfslp") 479 481 ! ! Coefficient 480 482 rn_aht_0 = 300. ! horizontal eddy diffusivity for tracers [m2/s] -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/NEMO/OPA_SRC/LDF/ldfslp.F90
r2772 r2980 46 46 ! !! Griffies operator 47 47 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: wslp2 !: wslp**2 from Griffies quarter cells 48 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:,:) :: triadi_g, triadj_g !: skew flux slopes relative to geopotentials 48 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:,:) :: triadi_g, triadj_g !: skew flux slopes relative to geopotentials 49 49 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:,:) :: triadi , triadj !: isoneutral slopes relative to model-coordinate 50 50 … … 58 58 59 59 ! Workspace arrays for ldf_slp_grif. These could be replaced by several 3D and 2D workspace 60 ! arrays from the wrk_nemo module with a bit of code re-writing. The 4D workspace 60 ! arrays from the wrk_nemo module with a bit of code re-writing. The 4D workspace 61 61 ! arrays can't be used here because of the zero-indexing of some of the ranks. ARPDBG. 62 62 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: zdzrho , zdyrho, zdxrho ! Horizontal and vertical density gradients … … 93 93 !!---------------------------------------------------------------------- 94 94 !! *** ROUTINE ldf_slp *** 95 !! 95 !! 96 96 !! ** Purpose : Compute the slopes of neutral surface (slope of isopycnal 97 97 !! surfaces referenced locally) (ln_traldf_iso=T). 98 !! 99 !! ** Method : The slope in the i-direction is computed at U- and 100 !! W-points (uslp, wslpi) and the slope in the j-direction is 98 !! 99 !! ** Method : The slope in the i-direction is computed at U- and 100 !! W-points (uslp, wslpi) and the slope in the j-direction is 101 101 !! computed at V- and W-points (vslp, wslpj). 102 102 !! They are bounded by 1/100 over the whole ocean, and within the … … 112 112 !! bottom slope (ln_sco=T) at level jpk in inildf] 113 113 !! 114 !! ** Action : - uslp, wslpi, and vslp, wslpj, the i- and j-slopes 114 !! ** Action : - uslp, wslpi, and vslp, wslpj, the i- and j-slopes 115 115 !! of now neutral surfaces at u-, w- and v- w-points, resp. 116 116 !!---------------------------------------------------------------------- … … 127 127 INTEGER :: ii0, ii1, iku ! temporary integer 128 128 INTEGER :: ij0, ij1, ikv ! temporary integer 129 REAL(wp) :: zeps, zm1_g, zm1_2g, z1_16 129 REAL(wp) :: zeps, zm1_g, zm1_2g, z1_16, zcofw ! local scalars 130 130 REAL(wp) :: zci, zfi, zau, zbu, zai, zbi ! - - 131 131 REAL(wp) :: zcj, zfj, zav, zbv, zaj, zbj ! - - … … 148 148 DO jj = 1, jpjm1 149 149 DO ji = 1, fs_jpim1 ! vector opt. 150 zgru(ji,jj,jk) = umask(ji,jj,jk) * ( prd(ji+1,jj ,jk) - prd(ji,jj,jk) ) 151 zgrv(ji,jj,jk) = vmask(ji,jj,jk) * ( prd(ji ,jj+1,jk) - prd(ji,jj,jk) ) 150 zgru(ji,jj,jk) = umask(ji,jj,jk) * ( prd(ji+1,jj ,jk) - prd(ji,jj,jk) ) 151 zgrv(ji,jj,jk) = vmask(ji,jj,jk) * ( prd(ji ,jj+1,jk) - prd(ji,jj,jk) ) 152 152 END DO 153 153 END DO 154 154 END DO 155 155 IF( ln_zps ) THEN ! partial steps correction at the bottom ocean level 156 # if defined key_vectopt_loop 156 # if defined key_vectopt_loop 157 157 DO jj = 1, 1 158 158 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) … … 161 161 DO ji = 1, jpim1 162 162 # endif 163 zgru(ji,jj,mbku(ji,jj)) = gru(ji,jj) 164 zgrv(ji,jj,mbkv(ji,jj)) = grv(ji,jj) 163 zgru(ji,jj,mbku(ji,jj)) = gru(ji,jj) 164 zgrv(ji,jj,mbkv(ji,jj)) = grv(ji,jj) 165 165 END DO 166 166 END DO … … 181 181 CALL ldf_slp_mxl( prd, pn2, zgru, zgrv, zdzr ) ! output: uslpml, vslpml, wslpiml, wslpjml 182 182 183 183 184 184 ! I. slopes at u and v point | uslp = d/di( prd ) / d/dz( prd ) 185 185 ! =========================== | vslp = d/dj( prd ) / d/dz( prd ) 186 ! 186 ! 187 187 DO jk = 2, jpkm1 !* Slopes at u and v points 188 188 DO jj = 2, jpjm1 … … 225 225 DO jk = 2, jpkm1 226 226 DO jj = 2, jpjm1, MAX(1, jpj-3) ! rows jj=2 and =jpjm1 only 227 DO ji = 2, jpim1 227 DO ji = 2, jpim1 228 228 uslp(ji,jj,jk) = z1_16 * ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 229 229 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & … … 266 266 ! II. slopes at w point | wslpi = mij( d/di( prd ) / d/dz( prd ) 267 267 ! =========================== | wslpj = mij( d/dj( prd ) / d/dz( prd ) 268 ! 268 ! 269 269 DO jk = 2, jpkm1 270 270 DO jj = 2, jpjm1 … … 308 308 DO jj = 2, jpjm1, MAX(1, jpj-3) ! rows jj=2 and =jpjm1 only 309 309 DO ji = 2, jpim1 310 zcofw = tmask(ji,jj,jk) * z1_16 310 311 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 311 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) &312 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) &313 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) &314 & + 4.* zwz(ji ,jj ,jk) ) * z1_16 * tmask(ji,jj,jk)312 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 313 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 314 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 315 & + 4.* zwz(ji ,jj ,jk) ) * zcofw 315 316 316 317 wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 317 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) &318 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) &319 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) &320 & + 4.* zww(ji ,jj ,jk) ) * z1_16 * tmask(ji,jj,jk)321 END DO 322 END DO 318 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 319 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 320 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 321 & + 4.* zww(ji ,jj ,jk) ) * zcofw 322 END DO 323 END DO 323 324 DO jj = 3, jpj-2 ! other rows 324 325 DO ji = fs_2, fs_jpim1 ! vector opt. 326 zcofw = tmask(ji,jj,jk) * z1_16 325 327 wslpi(ji,jj,jk) = ( zwz(ji-1,jj-1,jk) + zwz(ji+1,jj-1,jk) & 326 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) &327 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) &328 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) &329 & + 4.* zwz(ji ,jj ,jk) ) * z1_16 * tmask(ji,jj,jk)328 & + zwz(ji-1,jj+1,jk) + zwz(ji+1,jj+1,jk) & 329 & + 2.*( zwz(ji ,jj-1,jk) + zwz(ji-1,jj ,jk) & 330 & + zwz(ji+1,jj ,jk) + zwz(ji ,jj+1,jk) ) & 331 & + 4.* zwz(ji ,jj ,jk) ) * zcofw 330 332 331 333 wslpj(ji,jj,jk) = ( zww(ji-1,jj-1,jk) + zww(ji+1,jj-1,jk) & 332 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) &333 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) &334 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) &335 & + 4.* zww(ji ,jj ,jk) ) * z1_16 * tmask(ji,jj,jk)334 & + zww(ji-1,jj+1,jk) + zww(ji+1,jj+1,jk) & 335 & + 2.*( zww(ji ,jj-1,jk) + zww(ji-1,jj ,jk) & 336 & + zww(ji+1,jj ,jk) + zww(ji ,jj+1,jk) ) & 337 & + 4.* zww(ji ,jj ,jk) ) * zcofw 336 338 END DO 337 339 END DO … … 346 348 END DO 347 349 END DO 348 349 ! III. Specific grid points 350 ! =========================== 351 ! 350 351 ! III. Specific grid points 352 ! =========================== 353 ! 352 354 IF( cp_cfg == "orca" .AND. jp_cfg == 4 ) THEN ! ORCA_R4 configuration: horizontal diffusion in specific area 353 355 ! ! Gibraltar Strait … … 368 370 ENDIF 369 371 370 ! IV. Lateral boundary conditions 372 373 ! IV. Lateral boundary conditions 371 374 ! =============================== 372 375 CALL lbc_lnk( uslp , 'U', -1. ) ; CALL lbc_lnk( vslp , 'V', -1. ) … … 382 385 ! 383 386 END SUBROUTINE ldf_slp 384 387 385 388 386 389 SUBROUTINE ldf_slp_grif ( kt ) … … 390 393 !! ** Purpose : Compute the squared slopes of neutral surfaces (slope 391 394 !! of iso-pycnal surfaces referenced locally) (ln_traldf_grif=T) 392 !! at W-points using the Griffies quarter-cells. 393 !! 394 !! ** Method : calculates alpha and beta at T-points 395 !! at W-points using the Griffies quarter-cells. 396 !! 397 !! ** Method : calculates alpha and beta at T-points 395 398 !! 396 399 !! ** Action : - triadi_g, triadj_g T-pts i- and j-slope triads relative to geopot. (used for eiv) … … 399 402 !!---------------------------------------------------------------------- 400 403 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 401 USE oce , ONLY: zdit => ua , zdis => va ! (ua,va) used as workspace 402 USE oce , ONLY: zdjt => ta , zdjs => sa ! (ta,sa) used as workspace 403 USE wrk_nemo, ONLY: zdkt => wrk_3d_2 , zdks => wrk_3d_3 ! 3D workspace 404 USE wrk_nemo, ONLY: zalpha => wrk_3d_4 , zbeta => wrk_3d_5 ! alpha, beta at T points, at depth fsgdept 404 USE oce , ONLY: zalbet => ua ! use ua as workspace 405 405 USE wrk_nemo, ONLY: z1_mlbw => wrk_2d_1 406 ! 407 INTEGER, INTENT( in ) :: kt ! ocean time-step index408 ! 406 !! 407 INTEGER, INTENT( in ) :: kt ! ocean time-step index 408 !! 409 409 INTEGER :: ji, jj, jk, jl, ip, jp, kp ! dummy loop indices 410 INTEGER :: iku, ikv 411 REAL(wp) :: zfacti, zfactj , zatempw,zatempu,zatempv! local scalars412 REAL(wp) :: z bu, zbv, zbti, zbtj ! - -410 INTEGER :: iku, ikv ! local integer 411 REAL(wp) :: zfacti, zfactj ! local scalars 412 REAL(wp) :: zdit, zdis, zdjt, zdjs, zdkt, zdks, zbu, zbv, zbti, zbtj 413 413 REAL(wp) :: zdxrho_raw, zti_coord, zti_raw, zti_lim, zti_lim2, zti_g_raw, zti_g_lim 414 414 REAL(wp) :: zdyrho_raw, ztj_coord, ztj_raw, ztj_lim, ztj_lim2, ztj_g_raw, ztj_g_lim 415 415 REAL(wp) :: zdzrho_raw 416 REAL(wp) :: zbeta0 416 417 !!---------------------------------------------------------------------- 417 418 … … 424 425 !--------------------------------! 425 426 ! 426 CALL eos_alpbet( tsb, zalpha, zbeta ) !== before thermal and haline expension coeff. at T-points ==! 427 ! 428 DO jk = 1, jpkm1 !== before lateral T & S gradients at T-level jk ==! 429 DO jj = 1, jpjm1 430 DO ji = 1, fs_jpim1 ! vector opt. 431 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 432 zdis(ji,jj,jk) = ( sb(ji+1,jj,jk) - sb(ji,jj,jk) ) * umask(ji,jj,jk) 433 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 434 zdjs(ji,jj,jk) = ( sb(ji,jj+1,jk) - sb(ji,jj,jk) ) * vmask(ji,jj,jk) 435 END DO 436 END DO 437 END DO 438 IF( ln_zps ) THEN ! partial steps: correction at the last level 427 CALL eos_alpbet( tsb, zalbet, zbeta0 ) !== before local thermal/haline expension ratio at T-points ==! 428 ! 429 DO jl = 0, 1 !== unmasked before density i- j-, k-gradients ==! 430 ! 431 ip = jl ; jp = jl ! guaranteed nonzero gradients ( absolute value larger than repsln) 432 DO jk = 1, jpkm1 ! done each pair of triad 433 DO jj = 1, jpjm1 ! NB: not masked ==> a minimum value is set 434 DO ji = 1, fs_jpim1 ! vector opt. 435 zdit = ( tsb(ji+1,jj,jk,jp_tem) - tsb(ji,jj,jk,jp_tem) ) ! i-gradient of T & S at u-point 436 zdis = ( tsb(ji+1,jj,jk,jp_sal) - tsb(ji,jj,jk,jp_sal) ) 437 zdjt = ( tsb(ji,jj+1,jk,jp_tem) - tsb(ji,jj,jk,jp_tem) ) ! j-gradient of T & S at v-point 438 zdjs = ( tsb(ji,jj+1,jk,jp_sal) - tsb(ji,jj,jk,jp_sal) ) 439 zdxrho_raw = ( - zalbet(ji+ip,jj ,jk) * zdit + zbeta0*zdis ) / e1u(ji,jj) 440 zdyrho_raw = ( - zalbet(ji ,jj+jp,jk) * zdjt + zbeta0*zdjs ) / 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 ! 447 IF( ln_zps.and.l_grad_zps ) THEN ! partial steps: correction of i- & j-grad on bottom 439 448 # if defined key_vectopt_loop 440 DO jj = 1, 1441 DO ji = 1, jpij-jpi! vector opt. (forced unrolling)449 DO jj = 1, 1 450 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 442 451 # else 443 DO jj = 1, jpjm1444 DO ji = 1, jpim1452 DO jj = 1, jpjm1 453 DO ji = 1, jpim1 445 454 # endif 446 zdit(ji,jj,mbku(ji,jj)) = gtsu(ji,jj,jp_tem) ! i-gradient of T and S 447 zdis(ji,jj,mbku(ji,jj)) = gtsu(ji,jj,jp_sal) 448 zdjt(ji,jj,mbkv(ji,jj)) = gtsv(ji,jj,jp_tem) ! j-gradient of T and S 449 zdjs(ji,jj,mbkv(ji,jj)) = gtsv(ji,jj,jp_sal) 450 END DO 451 END DO 452 ENDIF 453 ! 454 zdkt(:,:,1) = 0._wp !== before vertical T & S gradient at w-level ==! 455 zdks(:,:,1) = 0._wp 456 DO jk = 2, jpk 457 zdkt(:,:,jk) = ( tb(:,:,jk-1) - tb(:,:,jk) ) * tmask(:,:,jk) 458 zdks(:,:,jk) = ( sb(:,:,jk-1) - sb(:,:,jk) ) * tmask(:,:,jk) 459 END DO 460 ! 461 ! 462 DO jl = 0, 1 !== density i-, j-, and k-gradients ==! 463 ip = jl ; jp = jl ! guaranteed nonzero gradients ( absolute value larger than repsln) 464 DO jk = 1, jpkm1 ! done each pair of triad 465 DO jj = 1, jpjm1 ! NB: not masked due to the minimum value set 466 DO ji = 1, fs_jpim1 ! vector opt. 467 zdxrho_raw = ( zalpha(ji+ip,jj ,jk) * zdit(ji,jj,jk) + zbeta(ji+ip,jj ,jk) * zdis(ji,jj,jk) ) / e1u(ji,jj) 468 zdyrho_raw = ( zalpha(ji ,jj+jp,jk) * zdjt(ji,jj,jk) + zbeta(ji ,jj+jp,jk) * zdjs(ji,jj,jk) ) / e2v(ji,jj) 469 zdxrho(ji+ip,jj ,jk,1-ip) = SIGN( MAX( repsln, ABS( zdxrho_raw ) ), zdxrho_raw ) ! keep the sign 470 zdyrho(ji ,jj+jp,jk,1-jp) = SIGN( MAX( repsln, ABS( zdyrho_raw ) ), zdyrho_raw ) 455 iku = mbku(ji,jj) ; ikv = mbkv(ji,jj) ! last ocean level (u- & v-points) 456 zdit = gtsu(ji,jj,jp_tem) ; zdjt = gtsv(ji,jj,jp_tem) ! i- & j-gradient of Temperature 457 zdis = gtsu(ji,jj,jp_sal) ; zdjs = gtsv(ji,jj,jp_sal) ! i- & j-gradient of Salinity 458 zdxrho_raw = ( - zalbet(ji+ip,jj ,iku) * zdit + zbeta0*zdis ) / e1u(ji,jj) 459 zdyrho_raw = ( - zalbet(ji ,jj+jp,ikv) * zdjt + zbeta0*zdjs ) / e2v(ji,jj) 460 zdxrho(ji+ip,jj ,iku,1-ip) = SIGN( MAX( repsln, ABS( zdxrho_raw ) ), zdxrho_raw ) ! keep the sign 461 zdyrho(ji ,jj+jp,ikv,1-jp) = SIGN( MAX( repsln, ABS( zdyrho_raw ) ), zdyrho_raw ) 471 462 END DO 472 463 END DO 473 END DO 474 END DO 475 DO kp = 0, 1 !== density i-, j-, and k-gradients ==! 476 DO jk = 1, jpkm1 ! done each pair of triad 477 DO jj = 1, jpj ! NB: not masked due to the minimum value set 478 DO ji = 1, jpi ! vector opt. 479 zdzrho_raw = ( zalpha(ji,jj,jk) * zdkt(ji,jj,jk+kp) + zbeta(ji,jj,jk) * zdks(ji,jj,jk+kp) ) & 480 & / fse3w(ji,jj,jk+kp) 481 zdzrho(ji ,jj ,jk, kp) = - MIN( - repsln, zdzrho_raw ) ! force zdzrho >= repsln 482 END DO 483 END DO 484 END DO 485 END DO 486 ! 487 DO jj = 1, jpj !== Reciprocal depth of the w-point below ML base ==! 464 ENDIF 465 ! 466 END DO 467 468 DO kp = 0, 1 !== unmasked before density i- j-, k-gradients ==! 469 DO jk = 1, jpkm1 ! done each pair of triad 470 DO jj = 1, jpj ! NB: not masked ==> a minimum value is set 471 DO ji = 1, jpi ! vector opt. 472 IF( jk+kp > 1 ) THEN ! k-gradient of T & S a jk+kp 473 zdkt = ( tsb(ji,jj,jk+kp-1,jp_tem) - tsb(ji,jj,jk+kp,jp_tem) ) 474 zdks = ( tsb(ji,jj,jk+kp-1,jp_sal) - tsb(ji,jj,jk+kp,jp_sal) ) 475 ELSE 476 zdkt = 0._wp ! 1st level gradient set to zero 477 zdks = 0._wp 478 ENDIF 479 zdzrho_raw = ( - zalbet(ji ,jj ,jk) * zdkt + zbeta0*zdks ) / fse3w(ji,jj,jk+kp) 480 zdzrho(ji ,jj ,jk, kp) = - MIN( - repsln, zdzrho_raw ) ! force zdzrho >= repsln 481 END DO 482 END DO 483 END DO 484 END DO 485 ! 486 DO jj = 1, jpj !== Reciprocal depth of the w-point below ML base ==! 488 487 DO ji = 1, jpi 489 488 jk = MIN( nmln(ji,jj), mbkt(ji,jj) ) + 1 ! MIN in case ML depth is the ocean depth … … 492 491 END DO 493 492 ! 494 ! !== intialisations to zero ==!495 ! 496 wslp2 (:,:,:) = 0._wp 497 triadi_g(:,:,1,:,:) = 0._wp ; triadi_g(:,:,jpk,:,:) = 0._wp 493 ! !== intialisations to zero ==! 494 ! 495 wslp2 (:,:,:) = 0._wp ! wslp2 will be cumulated 3D field set to zero 496 triadi_g(:,:,1,:,:) = 0._wp ; triadi_g(:,:,jpk,:,:) = 0._wp ! set surface and bottom slope to zero 498 497 triadj_g(:,:,1,:,:) = 0._wp ; triadj_g(:,:,jpk,:,:) = 0._wp 499 !!gm _iso set to zero missing500 triadi (:,:,1,:,:) = 0._wp ; triadj (:,:,jpk,:,:) = 0._wp! set surface and bottom slope to zero501 triadj (:,:,1,:,:) = 0._wp ; triadj(:,:,jpk,:,:) = 0._wp502 498 !!gm _iso set to zero missing 499 triadi (:,:,1,:,:) = 0._wp ; triadj (:,:,jpk,:,:) = 0._wp ! set surface and bottom slope to zero 500 triadj (:,:,1,:,:) = 0._wp ; triadj (:,:,jpk,:,:) = 0._wp 501 503 502 !-------------------------------------! 504 503 ! Triads just below the Mixed Layer ! 505 504 !-------------------------------------! 506 505 ! 507 DO jl = 0, 1 ! calculate slope of the 4 triads immediately ONE level below mixed-layer base508 DO kp = 0, 1 ! with only the slope-max limit and MASKED506 DO jl = 0, 1 ! calculate slope of the 4 triads immediately ONE level below mixed-layer base 507 DO kp = 0, 1 ! with only the slope-max limit and MASKED 509 508 DO jj = 1, jpjm1 510 509 DO ji = 1, fs_jpim1 511 510 ip = jl ; jp = jl 512 511 jk = MIN( nmln(ji+ip,jj) , mbkt(ji+ip,jj) ) + 1 ! ML level+1 (MIN in case ML depth is the ocean depth) 512 ! Add s-coordinate slope at t-points (do this by *subtracting* gradient of depth) 513 513 zti_g_raw = ( zdxrho(ji+ip,jj,jk-kp,1-ip) / zdzrho(ji+ip,jj,jk-kp,kp) & 514 & +( fsdept(ji+1,jj,jk-kp) - fsdept(ji,jj,jk-kp) ) / e1u(ji,jj) ) * umask(ji,jj,jk)514 & - ( fsdept(ji+1,jj,jk-kp) - fsdept(ji,jj,jk-kp) ) / e1u(ji,jj) ) * umask(ji,jj,jk) 515 515 jk = MIN( nmln(ji,jj+jp) , mbkt(ji,jj+jp) ) + 1 516 516 ztj_g_raw = ( zdyrho(ji,jj+jp,jk-kp,1-jp) / zdzrho(ji,jj+jp,jk-kp,kp) & 517 & +( fsdept(ji,jj+1,jk-kp) - fsdept(ji,jj,jk-kp) ) / e2v(ji,jj) ) * vmask(ji,jj,jk)517 & - ( fsdept(ji,jj+1,jk-kp) - fsdept(ji,jj,jk-kp) ) / e2v(ji,jj) ) * vmask(ji,jj,jk) 518 518 zti_mlb(ji+ip,jj ,1-ip,kp) = SIGN( MIN( rn_slpmax, ABS( zti_g_raw ) ), zti_g_raw ) 519 519 ztj_mlb(ji ,jj+jp,1-jp,kp) = SIGN( MIN( rn_slpmax, ABS( ztj_g_raw ) ), ztj_g_raw ) … … 527 527 !-------------------------------------! 528 528 ! 529 DO kp = 0, 1 ! k-index of triads529 DO kp = 0, 1 ! k-index of triads 530 530 DO jl = 0, 1 531 ip = jl ; jp = jl ! i- and j-indices of triads (i-k and j-k planes)531 ip = jl ; jp = jl ! i- and j-indices of triads (i-k and j-k planes) 532 532 DO jk = 1, jpkm1 533 533 DO jj = 1, jpjm1 534 DO ji = 1, fs_jpim1 ! vector opt.534 DO ji = 1, fs_jpim1 ! vector opt. 535 535 ! 536 536 ! Calculate slope relative to geopotentials used for GM skew fluxes 537 ! For s-coordinate, subtract slope at t-points (equivalent to *adding* gradient of depth)537 ! Add s-coordinate slope at t-points (do this by *subtracting* gradient of depth) 538 538 ! Limit by slope *relative to geopotentials* by rn_slpmax, and mask by psi-point 539 539 ! masked by umask taken at the level of dz(rho) … … 544 544 ztj_raw = zdyrho(ji ,jj+jp,jk,1-jp) / zdzrho(ji ,jj+jp,jk,kp) 545 545 zti_coord = ( fsdept(ji+1,jj ,jk) - fsdept(ji,jj,jk) ) / e1u(ji,jj) 546 ztj_coord = ( fsdept(ji ,jj+1,jk) - fsdept(ji,jj,jk) ) / e2v(ji,jj) 547 ! unmasked 548 zti_g_raw = zti_raw + zti_coord ! ref to geopot surfaces 549 ztj_g_raw = ztj_raw + ztj_coord 546 ztj_coord = ( fsdept(ji ,jj+1,jk) - fsdept(ji,jj,jk) ) / e2v(ji,jj) ! unmasked 547 zti_g_raw = zti_raw - zti_coord ! ref to geopot surfaces 548 ztj_g_raw = ztj_raw - ztj_coord 550 549 zti_g_lim = SIGN( MIN( rn_slpmax, ABS( zti_g_raw ) ), zti_g_raw ) 551 550 ztj_g_lim = SIGN( MIN( rn_slpmax, ABS( ztj_g_raw ) ), ztj_g_raw ) … … 562 561 ztj_g_lim = zfactj * ztj_g_lim & 563 562 & + ( 1._wp - zfactj ) * ztj_mlb(ji,jj+jp,1-jp,kp) & 564 & * fsdepw(ji,jj+jp,jk+kp) * z1_mlbw(ji,jj+jp) ! masked565 ! 566 triadi_g(ji+ip,jj ,jk,1-ip,kp) = zti_g_lim * umask(ji,jj,jk+kp) 563 & * fsdepw(ji,jj+jp,jk+kp) * z1_mlbw(ji,jj+jp) 564 ! 565 triadi_g(ji+ip,jj ,jk,1-ip,kp) = zti_g_lim * umask(ji,jj,jk+kp) ! masked 567 566 triadj_g(ji ,jj+jp,jk,1-jp,kp) = ztj_g_lim * vmask(ji,jj,jk+kp) 568 567 ! … … 574 573 ! Equivalent to tapering A_iso = sx_limited**2/(real slope)**2 575 574 ! 576 zti_lim = zti_g_lim - zti_coord ! remove the coordinate slope ==> relative to coordinate surfaces 577 ztj_lim = ztj_g_lim - ztj_coord 578 zti_lim2 = zti_lim * zti_lim * umask(ji,jj,jk+kp) ! square of limited slopes ! masked <<== 579 ztj_lim2 = ztj_lim * ztj_lim * vmask(ji,jj,jk+kp) 580 ! 575 zti_lim = ( zti_g_lim + zti_coord ) * umask(ji,jj,jk+kp) ! remove coordinate slope => relative to coordinate surfaces 576 ztj_lim = ( ztj_g_lim + ztj_coord ) * vmask(ji,jj,jk+kp) 577 zti_lim2 = zti_lim * zti_lim ! square of limited slopes ! masked <<== 578 ztj_lim2 = ztj_lim * ztj_lim 579 ! 580 IF( ln_triad_iso ) THEN 581 zti_raw = zti_lim2 / zti_raw 582 ztj_raw = ztj_lim2 / ztj_raw 583 zti_raw = SIGN( MIN( ABS(zti_lim), ABS( zti_raw ) ), zti_raw ) 584 ztj_raw = SIGN( MIN( ABS(ztj_lim), ABS( ztj_raw ) ), ztj_raw ) 585 zti_lim = zfacti * zti_lim & 586 & + ( 1._wp - zfacti ) * zti_raw 587 ztj_lim = zfactj * ztj_lim & 588 & + ( 1._wp - zfactj ) * ztj_raw 589 ENDIF 590 triadi(ji+ip,jj ,jk,1-ip,kp) = zti_lim 591 triadj(ji ,jj+jp,jk,1-jp,kp) = ztj_lim 592 ! 581 593 zbu = e1u(ji ,jj) * e2u(ji ,jj) * fse3u(ji ,jj,jk ) 582 594 zbv = e1v(ji ,jj) * e2v(ji ,jj) * fse3v(ji ,jj,jk ) … … 584 596 zbtj = e1t(ji,jj+jp) * e2t(ji,jj+jp) * fse3w(ji,jj+jp,jk+kp) 585 597 ! 586 triadi(ji+ip,jj ,jk,1-ip,kp) = zti_lim2 / zti_raw ! masked 587 triadj(ji ,jj+jp,jk,1-jp,kp) = ztj_lim2 / ztj_raw 588 ! 589 !!gm this may inhibit vectorization on Vect Computers, and even on scalar computers.... ==> to be checked 590 wslp2 (ji+ip,jj,jk+kp) = wslp2(ji+ip,jj,jk+kp) + 0.25_wp * zbu / zbti * zti_lim2 ! masked 598 !!gm this may inhibit vectorization on Vect Computers, and even on scalar computers.... ==> to be checked 599 ! agn may need to change to using ztj_g_lim**2, as wslp2 just used for eddy growth rate, needs *real* slope 600 wslp2 (ji+ip,jj,jk+kp) = wslp2(ji+ip,jj,jk+kp) + 0.25_wp * zbu / zbti * zti_lim2 ! masked 591 601 wslp2 (ji,jj+jp,jk+kp) = wslp2(ji,jj+jp,jk+kp) + 0.25_wp * zbv / zbtj * ztj_lim2 592 602 END DO … … 597 607 ! 598 608 wslp2(:,:,1) = 0._wp ! force the surface wslp to zero 599 609 600 610 CALL lbc_lnk( wslp2, 'W', 1. ) ! lateral boundary confition on wslp2 only ==>>> gm : necessary ? to be checked 601 611 ! … … 610 620 !! *** ROUTINE ldf_slp_mxl *** 611 621 !! 612 !! ** Purpose : Compute the slopes of iso-neutral surface just below 622 !! ** Purpose : Compute the slopes of iso-neutral surface just below 613 623 !! the mixed layer. 614 624 !! … … 619 629 !! 620 630 !! ** Action : uslpml, wslpiml : i- & j-slopes of neutral surfaces 621 !! vslpml, wslpjml just below the mixed layer 631 !! vslpml, wslpjml just below the mixed layer 622 632 !! omlmask : mixed layer mask 623 633 !!---------------------------------------------------------------------- … … 627 637 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: p_dzr ! z-gradient of density (T-point) 628 638 !! 629 INTEGER :: ji , jj , jk ! dummy loop indices630 INTEGER :: iku, ikv, ik, ikm1 ! local integers639 INTEGER :: ji , jj , jk ! dummy loop indices 640 INTEGER :: iku, ikv, ik, ikm1 ! local integers 631 641 REAL(wp) :: zeps, zm1_g, zm1_2g ! local scalars 632 642 REAL(wp) :: zci, zfi, zau, zbu, zai, zbi ! - - … … 644 654 wslpjml(1,:) = 0._wp ; wslpjml(jpi,:) = 0._wp 645 655 ! 646 ! !== surface mixed layer mask !647 DO jk = 1, jpk ! =1 inside the mixed layer, =0 otherwise656 ! !== surface mixed layer mask ! 657 DO jk = 1, jpk ! =1 inside the mixed layer, =0 otherwise 648 658 # if defined key_vectopt_loop 649 659 DO jj = 1, 1 650 DO ji = 1, jpij ! vector opt. (forced unrolling)660 DO ji = 1, jpij ! vector opt. (forced unrolling) 651 661 # else 652 662 DO jj = 1, jpj … … 679 689 DO ji = 2, jpim1 680 690 # endif 681 ! !== Slope at u- & v-points just below the Mixed Layer ==!691 ! !== Slope at u- & v-points just below the Mixed Layer ==! 682 692 ! 683 ! 693 ! !- vertical density gradient for u- and v-slopes (from dzr at T-point) 684 694 iku = MIN( MAX( 1, nmln(ji,jj) , nmln(ji+1,jj) ) , jpkm1 ) ! ML (MAX of T-pts, bound by jpkm1) 685 ikv = MIN( MAX( 1, nmln(ji,jj) , nmln(ji,jj+1) ) , jpkm1 ) ! 695 ikv = MIN( MAX( 1, nmln(ji,jj) , nmln(ji,jj+1) ) , jpkm1 ) ! 686 696 zbu = 0.5_wp * ( p_dzr(ji,jj,iku) + p_dzr(ji+1,jj ,iku) ) 687 697 zbv = 0.5_wp * ( p_dzr(ji,jj,ikv) + p_dzr(ji ,jj+1,ikv) ) 688 ! 698 ! !- horizontal density gradient at u- & v-points 689 699 zau = p_gru(ji,jj,iku) / e1u(ji,jj) 690 700 zav = p_grv(ji,jj,ikv) / e2v(ji,jj) 691 ! 692 ! 701 ! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0 702 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 693 703 zbu = MIN( zbu , -100._wp* ABS( zau ) , -7.e+3_wp/fse3u(ji,jj,iku)* ABS( zau ) ) 694 704 zbv = MIN( zbv , -100._wp* ABS( zav ) , -7.e+3_wp/fse3v(ji,jj,ikv)* ABS( zav ) ) 695 ! 705 ! !- Slope at u- & v-points (uslpml, vslpml) 696 706 uslpml(ji,jj) = zau / ( zbu - zeps ) * umask(ji,jj,iku) 697 707 vslpml(ji,jj) = zav / ( zbv - zeps ) * vmask(ji,jj,ikv) 698 708 ! 699 ! !== i- & j-slopes at w-points just below the Mixed Layer ==!709 ! !== i- & j-slopes at w-points just below the Mixed Layer ==! 700 710 ! 701 711 ik = MIN( nmln(ji,jj) + 1, jpk ) 702 712 ikm1 = MAX( 1, ik-1 ) 703 ! 713 ! !- vertical density gradient for w-slope (from N^2) 704 714 zbw = zm1_2g * pn2 (ji,jj,ik) * ( prd (ji,jj,ik) + prd (ji,jj,ikm1) + 2. ) 705 ! 715 ! !- horizontal density i- & j-gradient at w-points 706 716 zci = MAX( umask(ji-1,jj,ik ) + umask(ji,jj,ik ) & 707 & + umask(ji-1,jj,ikm1) + umask(ji,jj,ikm1) , zeps ) * e1t(ji,jj) 717 & + umask(ji-1,jj,ikm1) + umask(ji,jj,ikm1) , zeps ) * e1t(ji,jj) 708 718 zcj = MAX( vmask(ji,jj-1,ik ) + vmask(ji,jj,ik ) & 709 719 & + vmask(ji,jj-1,ikm1) + vmask(ji,jj,ikm1) , zeps ) * e2t(ji,jj) … … 712 722 zaj = ( p_grv(ji,jj-1,ik ) + p_grv(ji,jj,ik ) & 713 723 & + p_grv(ji,jj-1,ikm1) + p_grv(ji,jj,ikm1) ) / zcj * tmask(ji,jj,ik) 714 ! 715 ! 724 ! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0. 725 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 716 726 zbi = MIN( zbw , -100._wp* ABS( zai ) , -7.e+3_wp/fse3w(ji,jj,ik)* ABS( zai ) ) 717 727 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/fse3w(ji,jj,ik)* ABS( zaj ) ) 718 ! 728 ! !- i- & j-slope at w-points (wslpiml, wslpjml) 719 729 wslpiml(ji,jj) = zai / ( zbi - zeps ) * tmask (ji,jj,ik) 720 730 wslpjml(ji,jj) = zaj / ( zbj - zeps ) * tmask (ji,jj,ik) 721 731 END DO 722 732 END DO 723 !!gm this lbc_lnk should be useless....733 !!gm this lbc_lnk should be useless.... 724 734 CALL lbc_lnk( uslpml , 'U', -1. ) ; CALL lbc_lnk( vslpml , 'V', -1. ) ! lateral boundary cond. (sign change) 725 735 CALL lbc_lnk( wslpiml, 'W', -1. ) ; CALL lbc_lnk( wslpjml, 'W', -1. ) ! lateral boundary conditions … … 734 744 !! ** Purpose : Initialization for the isopycnal slopes computation 735 745 !! 736 !! ** Method : read the nammbf namelist and check the parameter 737 !! 746 !! ** Method : read the nammbf namelist and check the parameter 747 !! values called by tra_dmp at the first timestep (nit000) 738 748 !!---------------------------------------------------------------------- 739 749 INTEGER :: ji, jj, jk ! dummy loop indices 740 750 INTEGER :: ierr ! local integer 741 751 !!---------------------------------------------------------------------- 742 743 IF(lwp) THEN 752 753 IF(lwp) THEN 744 754 WRITE(numout,*) 745 755 WRITE(numout,*) 'ldf_slp_init : direction of lateral mixing' 746 756 WRITE(numout,*) '~~~~~~~~~~~~' 747 757 ENDIF 748 758 749 759 IF( ln_traldf_grif ) THEN ! Griffies operator : triad of slopes 750 760 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 ) … … 755 765 IF( ln_dynldf_iso ) CALL ctl_stop( 'ldf_slp_init: Griffies operator on momentum not supported' ) 756 766 ! 757 IF( ( ln_traldf_hor .OR. ln_dynldf_hor ) .AND. ln_sco ) &758 CALL ctl_stop( 'ldf_slp_init: horizontal Griffies operator in s-coordinate not supported' )759 !760 767 ELSE ! Madec operator : slopes at u-, v-, and w-points 761 768 ALLOCATE( uslp(jpi,jpj,jpk) , vslp(jpi,jpj,jpk) , wslpi(jpi,jpj,jpk) , wslpj(jpi,jpj,jpk) , & … … 770 777 wslpj(:,:,:) = 0._wp ; wslpjml(:,:) = 0._wp 771 778 772 !!gm I no longer understand this.....779 !!gm I no longer understand this..... 773 780 IF( (ln_traldf_hor .OR. ln_dynldf_hor) .AND. .NOT. (lk_vvl .AND. ln_rstart) ) THEN 774 781 IF(lwp) WRITE(numout,*) ' Horizontal mixing in s-coordinate: slope = slope of s-surfaces' … … 803 810 LOGICAL, PUBLIC, PARAMETER :: lk_ldfslp = .FALSE. !: slopes flag 804 811 CONTAINS 805 SUBROUTINE ldf_slp( kt, prd, pn2 ) 806 INTEGER, INTENT(in) :: kt 812 SUBROUTINE ldf_slp( kt, prd, pn2 ) ! Dummy routine 813 INTEGER, INTENT(in) :: kt 807 814 REAL, DIMENSION(:,:,:), INTENT(in) :: prd, pn2 808 815 WRITE(*,*) 'ldf_slp: You should not have seen this print! error?', kt, prd(1,1,1), pn2(1,1,1) … … 812 819 WRITE(*,*) 'ldf_slp_grif: You should not have seen this print! error?', kt 813 820 END SUBROUTINE ldf_slp_grif 814 SUBROUTINE ldf_slp_init ! Dummy routine821 SUBROUTINE ldf_slp_init ! Dummy routine 815 822 END SUBROUTINE ldf_slp_init 816 823 #endif -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/NEMO/OPA_SRC/LDF/ldftra.F90
r2715 r2980 67 67 & ln_traldf_level, ln_traldf_hor , ln_traldf_iso, & 68 68 & ln_traldf_grif , ln_traldf_gdia, & 69 & ln_triad_iso , ln_botmix_grif, & 69 70 & rn_aht_0 , rn_ahtb_0 , rn_aeiv_0, & 70 71 & rn_slpmax … … 94 95 WRITE(numout,*) ' maximum isoppycnal slope rn_slpmax = ', rn_slpmax 95 96 WRITE(numout,*) ' + griffies operator internal controls not set via the namelist (experimental): ' 96 WRITE(numout,*) ' calculate triads twice l _triad_iso = ', l_triad_iso97 WRITE(numout,*) ' no Shapiro filter l_no_smooth = ', l_no_smooth97 WRITE(numout,*) ' calculate triads twice ln_triad_iso = ', ln_triad_iso 98 WRITE(numout,*) ' GM -->lat mixing on bottom ln_botmix_grif = ', ln_botmix_grif 98 99 WRITE(numout,*) 99 100 ENDIF -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/NEMO/OPA_SRC/LDF/ldftra_oce.F90
r2715 r2980 32 32 33 33 REAL(wp), PUBLIC :: aht0, ahtb0, aeiv0 !!: OLD namelist names 34 LOGICAL , PUBLIC :: l_triad_iso = .FALSE. !: calculate triads twice 35 LOGICAL , PUBLIC :: l_no_smooth = .FALSE. !: no Shapiro smoothing 34 LOGICAL , PUBLIC :: ln_triad_iso = .FALSE. !: calculate triads twice 35 LOGICAL , PUBLIC :: ln_botmix_grif = .FALSE. !: mixing on bottom 36 LOGICAL , PUBLIC :: l_grad_zps = .FALSE. !: special treatment for Horz Tgradients w partial steps 36 37 37 38 #if defined key_traldf_c3d -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/NEMO/OPA_SRC/TRA/eosbn2.F90
r2715 r2980 3 3 !! *** MODULE eosbn2 *** 4 4 !! Ocean diagnostic variable : equation of state - in situ and potential density 5 !! - Brunt-Vaisala frequency 5 !! - Brunt-Vaisala frequency 6 6 !!============================================================================== 7 7 !! History : OPA ! 1989-03 (O. Marti) Original code … … 27 27 !! eos_insitu_2d : Compute the in situ density for 2d fields 28 28 !! eos_bn2 : Compute the Brunt-Vaisala frequency 29 !! eos_alpbet : calculates the in situ thermal and haline expansion coeff.29 !! eos_alpbet : calculates the in situ thermal/haline expansion ratio 30 30 !! tfreez : Compute the surface freezing temperature 31 31 !! eos_init : set eos parameters (namelist) … … 41 41 PRIVATE 42 42 43 ! !! * Interface 43 ! !! * Interface 44 44 INTERFACE eos 45 45 MODULE PROCEDURE eos_insitu, eos_insitu_pot, eos_insitu_2d 46 END INTERFACE 46 END INTERFACE 47 47 INTERFACE bn2 48 48 MODULE PROCEDURE eos_bn2 49 END INTERFACE 49 END INTERFACE 50 50 51 51 PUBLIC eos ! called by step, istate, tranpc and zpsgrd modules … … 61 61 62 62 REAL(wp), PUBLIC :: ralpbet !: alpha / beta ratio 63 63 64 64 !! * Substitutions 65 65 # include "domzgr_substitute.h90" … … 75 75 !!---------------------------------------------------------------------- 76 76 !! *** ROUTINE eos_insitu *** 77 !! 78 !! ** Purpose : Compute the in situ density (ratio rho/rau0) from 77 !! 78 !! ** Purpose : Compute the in situ density (ratio rho/rau0) from 79 79 !! potential temperature and salinity using an equation of state 80 80 !! defined through the namelist parameter nn_eos. … … 134 134 !CDIR NOVERRCHK 135 135 zws(:,:,:) = SQRT( ABS( pts(:,:,:,jp_sal) ) ) 136 ! 136 ! 137 137 DO jk = 1, jpkm1 138 138 DO jj = 1, jpj … … 199 199 !!---------------------------------------------------------------------- 200 200 !! *** ROUTINE eos_insitu_pot *** 201 !! 201 !! 202 202 !! ** Purpose : Compute the in situ density (ratio rho/rau0) and the 203 203 !! potential volumic mass (Kg/m3) from potential temperature and 204 !! salinity fields using an equation of state defined through the 204 !! salinity fields using an equation of state defined through the 205 205 !! namelist parameter nn_eos. 206 206 !! … … 230 230 !! nn_eos = 2 : linear equation of state function of temperature and 231 231 !! salinity 232 !! prd(t,s) = ( rho(t,s) - rau0 ) / rau0 232 !! prd(t,s) = ( rho(t,s) - rau0 ) / rau0 233 233 !! = rn_beta * s - rn_alpha * tn - 1. 234 234 !! rhop(t,s) = rho(t,s) … … 265 265 !CDIR NOVERRCHK 266 266 zws(:,:,:) = SQRT( ABS( pts(:,:,:,jp_sal) ) ) 267 ! 267 ! 268 268 DO jk = 1, jpkm1 269 269 DO jj = 1, jpj … … 336 336 !! *** ROUTINE eos_insitu_2d *** 337 337 !! 338 !! ** Purpose : Compute the in situ density (ratio rho/rau0) from 338 !! ** Purpose : Compute the in situ density (ratio rho/rau0) from 339 339 !! potential temperature and salinity using an equation of state 340 340 !! defined through the namelist parameter nn_eos. * 2D field case … … 374 374 ! ! 2 : salinity [psu] 375 375 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ) :: pdep ! depth [m] 376 REAL(wp), DIMENSION(jpi,jpj) , INTENT( out) :: prd ! in situ density 376 REAL(wp), DIMENSION(jpi,jpj) , INTENT( out) :: prd ! in situ density 377 377 !! 378 378 INTEGER :: ji, jj ! dummy loop indices … … 449 449 DO jj = 1, jpjm1 450 450 DO ji = 1, fs_jpim1 ! vector opt. 451 prd(ji,jj) = ( rn_beta * pts(ji,jj,jp_sal) - rn_alpha * pts(ji,jj,jp_tem) ) * tmask(ji,jj,1) 451 prd(ji,jj) = ( rn_beta * pts(ji,jj,jp_sal) - rn_alpha * pts(ji,jj,jp_tem) ) * tmask(ji,jj,1) 452 452 END DO 453 453 END DO … … 468 468 !! ** Purpose : Compute the local Brunt-Vaisala frequency at the time- 469 469 !! step of the input arguments 470 !! 470 !! 471 471 !! ** Method : 472 472 !! * nn_eos = 0 : UNESCO sea water properties … … 482 482 !! N^2 = grav * (rn_alpha * dk[ t ] - rn_beta * dk[ s ] ) / e3w 483 483 !! The use of potential density to compute N^2 introduces e r r o r 484 !! in the sign of N^2 at great depths. We recommand the use of 484 !! in the sign of N^2 at great depths. We recommand the use of 485 485 !! nn_eos = 0, except for academical studies. 486 486 !! Macro-tasked on horizontal slab (jk-loop) … … 497 497 !! 498 498 INTEGER :: ji, jj, jk ! dummy loop indices 499 REAL(wp) :: zgde3w, zt, zs, zh, zalbet, zbeta ! local scalars 499 REAL(wp) :: zgde3w, zt, zs, zh, zalbet, zbeta ! local scalars 500 500 #if defined key_zdfddm 501 501 REAL(wp) :: zds ! local scalars … … 504 504 505 505 ! pn2 : interior points only (2=< jk =< jpkm1 ) 506 ! -------------------------- 506 ! -------------------------- 507 507 ! 508 508 SELECT CASE( nn_eos ) … … 542 542 & - 0.121555e-07_wp ) * zh 543 543 ! 544 pn2(ji,jj,jk) = zgde3w * zbeta * tmask(ji,jj,jk) & ! N^2 544 pn2(ji,jj,jk) = zgde3w * zbeta * tmask(ji,jj,jk) & ! N^2 545 545 & * ( zalbet * ( pts(ji,jj,jk-1,jp_tem) - pts(ji,jj,jk,jp_tem) ) & 546 546 & - ( pts(ji,jj,jk-1,jp_sal) - pts(ji,jj,jk,jp_sal) ) ) … … 565 565 & - rn_beta * ( pts(:,:,jk-1,jp_sal) - pts(:,:,jk,jp_sal) ) ) & 566 566 & / fse3w(:,:,jk) * tmask(:,:,jk) 567 END DO 567 END DO 568 568 #if defined key_zdfddm 569 569 DO jk = 2, jpkm1 ! Rrau = (alpha / beta) (dk[t] / dk[s]) 570 570 DO jj = 1, jpj 571 571 DO ji = 1, jpi 572 zds = ( pts(ji,jj,jk-1,jp_sal) - pts(ji,jj,jk,jp_sal) ) 572 zds = ( pts(ji,jj,jk-1,jp_sal) - pts(ji,jj,jk,jp_sal) ) 573 573 IF ( ABS( zds ) <= 1.e-20_wp ) zds = 1.e-20_wp 574 574 rrau(ji,jj,jk) = ralpbet * ( pts(ji,jj,jk-1,jp_tem) - pts(ji,jj,jk,jp_tem) ) / zds … … 587 587 588 588 589 SUBROUTINE eos_alpbet( pts, palp h, pbeta)590 !!---------------------------------------------------------------------- 591 !! *** ROUTINE ldf_slp_grif***592 !! 593 !! ** Purpose : Calculates the thermal and haline expansion coefficients at T-points.594 !! 595 !! ** Method : calculates alpha and beta at T-points589 SUBROUTINE eos_alpbet( pts, palpbet, beta0 ) 590 !!---------------------------------------------------------------------- 591 !! *** ROUTINE eos_alpbet *** 592 !! 593 !! ** Purpose : Calculates the in situ thermal/haline expansion ratio at T-points 594 !! 595 !! ** Method : calculates alpha / beta ratio at T-points 596 596 !! * nn_eos = 0 : UNESCO sea water properties 597 !! The brunt-vaisala frequency is computed using the polynomial 598 !! polynomial expression of McDougall (1987): 599 !! N^2 = grav * beta * ( alpha/beta*dk[ t ] - dk[ s ] )/e3w 600 !! If lk_zdfddm=T, the heat/salt buoyancy flux ratio Rrau is 601 !! computed and used in zdfddm module : 602 !! Rrau = alpha/beta * ( dk[ t ] / dk[ s ] ) 597 !! The alpha/beta ratio is returned as 3-D array palpbet using the polynomial 598 !! polynomial expression of McDougall (1987). 599 !! Scalar beta0 is returned = 1. 603 600 !! * nn_eos = 1 : linear equation of state (temperature only) 604 !! N^2 = grav * rn_alpha * dk[ t ]/e3w 601 !! The ratio is undefined, so we return alpha as palpbet 602 !! Scalar beta0 is returned = 0. 605 603 !! * nn_eos = 2 : linear equation of state (temperature & salinity) 606 !! N^2 = grav * (rn_alpha * dk[ t ] - rn_beta * dk[ s ] ) / e3w 607 !! * nn_eos = 3 : Jackett JAOT 2003 ??? 608 !! 609 !! ** Action : - palph, pbeta : thermal and haline expansion coeff. at T-point 610 !!---------------------------------------------------------------------- 611 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! pot. temperature & salinity 612 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT( out) :: palph, pbeta ! thermal & haline expansion coeff. 613 ! 604 !! The alpha/beta ratio is returned as ralpbet 605 !! Scalar beta0 is returned = 1. 606 !! 607 !! ** Action : - palpbet : thermal/haline expansion ratio at T-points 608 !! : beta0 : 1. or 0. 609 !!---------------------------------------------------------------------- 610 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! pot. temperature & salinity 611 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT( out) :: palpbet ! thermal/haline expansion ratio 612 REAL(wp), INTENT( out) :: beta0 ! set = 1 except with case 1 eos, rho=rho(T) 613 !! 614 614 INTEGER :: ji, jj, jk ! dummy loop indices 615 REAL(wp) :: zt, zs, zh ! local scalars 615 REAL(wp) :: zt, zs, zh ! local scalars 616 616 !!---------------------------------------------------------------------- 617 617 ! … … 624 624 zt = pts(ji,jj,jk,jp_tem) ! potential temperature 625 625 zs = pts(ji,jj,jk,jp_sal) - 35._wp ! salinity anomaly (s-35) 626 zh = fsdept(ji,jj,jk) ! depth in meters 627 ! 628 pbeta(ji,jj,jk) = ( ( -0.415613e-09_wp * zt + 0.555579e-07_wp ) * zt & 629 & - 0.301985e-05_wp ) * zt & 630 & + 0.785567e-03_wp & 631 & + ( 0.515032e-08_wp * zs & 632 & + 0.788212e-08_wp * zt - 0.356603e-06_wp ) * zs & 633 & + ( ( 0.121551e-17_wp * zh & 634 & - 0.602281e-15_wp * zs & 635 & - 0.175379e-14_wp * zt + 0.176621e-12_wp ) * zh & 636 & + 0.408195e-10_wp * zs & 637 & + ( - 0.213127e-11_wp * zt + 0.192867e-09_wp ) * zt & 638 & - 0.121555e-07_wp ) * zh 639 ! 640 palph(ji,jj,jk) = - pbeta(ji,jj,jk) * & 641 & ((( ( - 0.255019e-07_wp * zt + 0.298357e-05_wp ) * zt & 642 & - 0.203814e-03_wp ) * zt & 643 & + 0.170907e-01_wp ) * zt & 644 & + 0.665157e-01_wp & 645 & + ( - 0.678662e-05_wp * zs & 646 & - 0.846960e-04_wp * zt + 0.378110e-02_wp ) * zs & 647 & + ( ( - 0.302285e-13_wp * zh & 648 & - 0.251520e-11_wp * zs & 649 & + 0.512857e-12_wp * zt * zt ) * zh & 650 & - 0.164759e-06_wp * zs & 651 & +( 0.791325e-08_wp * zt - 0.933746e-06_wp ) * zt & 652 & + 0.380374e-04_wp ) * zh) 626 zh = fsdept(ji,jj,jk) ! depth in meters 627 ! 628 palpbet(ji,jj,jk) = & 629 & ( ( ( - 0.255019e-07_wp * zt + 0.298357e-05_wp ) * zt & 630 & - 0.203814e-03_wp ) * zt & 631 & + 0.170907e-01_wp ) * zt & 632 & + 0.665157e-01_wp & 633 & + ( - 0.678662e-05_wp * zs & 634 & - 0.846960e-04_wp * zt + 0.378110e-02_wp ) * zs & 635 & + ( ( - 0.302285e-13_wp * zh & 636 & - 0.251520e-11_wp * zs & 637 & + 0.512857e-12_wp * zt * zt ) * zh & 638 & - 0.164759e-06_wp * zs & 639 & +( 0.791325e-08_wp * zt - 0.933746e-06_wp ) * zt & 640 & + 0.380374e-04_wp ) * zh 653 641 END DO 654 642 END DO 655 643 END DO 656 ! 657 CASE ( 1 ) 658 palph(:,:,:) = - rn_alpha 659 pbeta(:,:,:) = 0._wp 660 ! 661 CASE ( 2 ) 662 palph(:,:,:) = - rn_alpha 663 pbeta(:,:,:) = rn_beta 644 beta0 = 1._wp 645 ! 646 CASE ( 1 ) !== Linear formulation = F( temperature ) ==! 647 palpbet(:,:,:) = rn_alpha 648 beta0 = 0._wp 649 ! 650 CASE ( 2 ) !== Linear formulation = F( temperature , salinity ) ==! 651 palpbet(:,:,:) = ralpbet 652 beta0 = 1._wp 664 653 ! 665 654 CASE DEFAULT … … 747 736 748 737 !!====================================================================== 749 END MODULE eosbn2 738 END MODULE eosbn2 -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/NEMO/OPA_SRC/TRA/traldf_iso_grif.F90
r2715 r2980 4 4 !! Ocean tracers: horizontal component of the lateral tracer mixing trend 5 5 !!====================================================================== 6 !! History : 3.3 ! 2010-10 (G. Nurser, C. Harris, G. Madec) 6 !! History : 3.3 ! 2010-10 (G. Nurser, C. Harris, G. Madec) 7 7 !! ! Griffies operator version adapted from traldf_iso.F90 8 8 !!---------------------------------------------------------------------- … … 11 11 !! 'key_ldfslp' slope of the lateral diffusive direction 12 12 !!---------------------------------------------------------------------- 13 !! tra_ldf_iso_grif : update the tracer trend with the horizontal component 14 !! of the Griffies iso-neutral laplacian operator 13 !! tra_ldf_iso_grif : update the tracer trend with the horizontal component 14 !! of the Griffies iso-neutral laplacian operator 15 15 !!---------------------------------------------------------------------- 16 16 USE oce ! ocean dynamics and active tracers … … 34 34 REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE, SAVE :: psix_eiv, psiy_eiv !: eiv stream function (diag only) 35 35 REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE, SAVE :: ah_wslp2 !: aeiv*w-slope^2 36 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, SAVE :: zdkt ! atypic workspace36 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, SAVE :: zdkt3d !: vertical tracer gradient at 2 levels 37 37 38 38 !! * Substitutions … … 53 53 !! *** ROUTINE tra_ldf_iso_grif *** 54 54 !! 55 !! ** Purpose : Compute the before horizontal tracer (t & s) diffusive 56 !! trend for a laplacian tensor (ezxcept the dz[ dz[.] ] term) and 55 !! ** Purpose : Compute the before horizontal tracer (t & s) diffusive 56 !! trend for a laplacian tensor (ezxcept the dz[ dz[.] ] term) and 57 57 !! add it to the general trend of tracer equation. 58 58 !! 59 !! ** Method : The horizontal component of the lateral diffusive trends 59 !! ** Method : The horizontal component of the lateral diffusive trends 60 60 !! is provided by a 2nd order operator rotated along neural or geopo- 61 61 !! tential surfaces to which an eddy induced advection can be added … … 67 67 !! 68 68 !! 2nd part : horizontal fluxes of the lateral mixing operator 69 !! ======== 69 !! ======== 70 70 !! zftu = (aht+ahtb0) e2u*e3u/e1u di[ tb ] 71 71 !! - aht e2u*uslp dk[ mi(mk(tb)) ] … … 99 99 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 100 100 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields 101 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 101 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 102 102 REAL(wp) , INTENT(in ) :: pahtb0 ! background diffusion coef 103 103 ! … … 108 108 REAL(wp) :: zmskv, zabe2, zcof2, zcoef4 ! - - 109 109 REAL(wp) :: zcoef0, zbtr ! - - 110 !REAL(wp), POINTER, DIMENSION(:,:,:) :: zdkt ! 2D+1 workspace111 110 ! 112 111 REAL(wp) :: zslope_skew, zslope_iso, zslope2, zbu, zbv … … 121 120 CALL ctl_stop('tra_ldf_iso_grif: requested workspace arrays unavailable.') ; RETURN 122 121 ENDIF 123 ! ARP - line below uses 'bounds re-mapping' which is only defined in 124 ! Fortran 2003 and up. We would be OK if code was written to use 125 ! zdkt(:,:,1:2) instead as then wouldn't need to re-map bounds. 126 ! As it is, we make zdkt a module array and allocate it in _alloc(). 127 !zdkt(1:jpi,1:jpj,0:1) => wrk_3d_9(:,:,1:2) 128 129 IF( kt == nit000 ) THEN 122 123 IF( kt == nit000.AND..NOT.ALLOCATED(ah_wslp2) ) THEN 130 124 IF(lwp) WRITE(numout,*) 131 125 IF(lwp) WRITE(numout,*) 'tra_ldf_iso_grif : rotated laplacian diffusion operator on ', cdtype 132 IF(lwp) WRITE(numout,*) ' WARNING: STILL UNDER TEST, NOT RECOMMENDED. USE AT YOUR OWN PERIL'133 126 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 134 ALLOCATE( ah_wslp2(jpi,jpj,jpk) , zdkt (jpi,jpj,0:1), STAT=ierr )127 ALLOCATE( ah_wslp2(jpi,jpj,jpk) , zdkt3d(jpi,jpj,0:1), STAT=ierr ) 135 128 IF( lk_mpp ) CALL mpp_sum ( ierr ) 136 129 IF( ierr > 0 ) CALL ctl_stop('STOP', 'tra_ldf_iso_grif: unable to allocate arrays') … … 143 136 144 137 !!---------------------------------------------------------------------- 145 !! 0 - calculate ah_wslp2, psix_eiv, psiy_eiv 138 !! 0 - calculate ah_wslp2, psix_eiv, psiy_eiv 146 139 !!---------------------------------------------------------------------- 147 140 148 !!gm Future development: consider using Ah defined at T-points and attached to the 4 t-point triads141 !!gm Future development: consider using Ah defined at T-points and attached to the 4 t-point triads 149 142 150 143 ah_wslp2(:,:,:) = 0._wp … … 159 152 DO jj = 1, jpjm1 160 153 DO ji = 1, fs_jpim1 154 ze1ur = 1._wp / e1u(ji,jj) 161 155 ze3wr = 1._wp / fse3w(ji+ip,jj,jk+kp) 162 156 zbu = 0.25_wp * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) 163 zah = fsahtu(ji,jj,jk) ! 157 zah = fsahtu(ji,jj,jk) ! fsaht(ji+ip,jj,jk) 164 158 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 165 zslope2 = zslope_skew - ( fsdept(ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * ze1ur * umask(ji,jj,jk+kp) 159 ! Subtract s-coordinate slope at t-points to give slope rel to s surfaces 160 ! (do this by *adding* gradient of depth) 161 zslope2 = zslope_skew + ( fsdept(ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * ze1ur * umask(ji,jj,jk+kp) 166 162 zslope2 = zslope2 *zslope2 167 163 ah_wslp2(ji+ip,jj,jk+kp) = ah_wslp2(ji+ip,jj,jk+kp) & 168 164 & + zah * ( zbu * ze3wr / ( e1t(ji+ip,jj) * e2t(ji+ip,jj) ) ) * zslope2 169 165 IF( ln_traldf_gdia ) THEN 170 zaei_slp = fsaeiw(ji+ip,jj,jk) * zslope_skew !fsaeit(ji+ip,jj,jk)*zslope_skew166 zaei_slp = fsaeiw(ji+ip,jj,jk) * zslope_skew ! fsaeit(ji+ip,jj,jk)*zslope_skew 171 167 psix_eiv(ji,jj,jk+kp) = psix_eiv(ji,jj,jk+kp) + 0.25_wp * zaei_slp 172 168 ENDIF … … 182 178 DO jj = 1, jpjm1 183 179 DO ji=1,fs_jpim1 180 ze2vr = 1._wp / e2v(ji,jj) 184 181 ze3wr = 1.0_wp / fse3w(ji,jj+jp,jk+kp) 185 182 zbv = 0.25_wp * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) 186 zah = fsaht u(ji,jj,jk) !fsaht(ji,jj+jp,jk)183 zah = fsahtv(ji,jj,jk) ! fsaht(ji,jj+jp,jk) 187 184 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 188 zslope2 = zslope_skew - ( fsdept(ji,jj+1,jk) - fsdept(ji,jj,jk) ) * ze2vr * vmask(ji,jj,jk+kp) 185 ! Subtract s-coordinate slope at t-points to give slope rel to s surfaces 186 ! (do this by *adding* gradient of depth) 187 zslope2 = zslope_skew + ( fsdept(ji,jj+1,jk) - fsdept(ji,jj,jk) ) * ze2vr * vmask(ji,jj,jk+kp) 189 188 zslope2 = zslope2 * zslope2 190 189 ah_wslp2(ji,jj+jp,jk+kp) = ah_wslp2(ji,jj+jp,jk+kp) & 191 190 & + zah * ( zbv * ze3wr / ( e1t(ji,jj+jp) * e2t(ji,jj+jp) ) ) * zslope2 192 191 IF( ln_traldf_gdia ) THEN 193 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew !fsaeit(ji,jj+jp,jk)*zslope_skew192 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew ! fsaeit(ji,jj+jp,jk)*zslope_skew 194 193 psiy_eiv(ji,jj,jk+kp) = psiy_eiv(ji,jj,jk+kp) + 0.25_wp * zaei_slp 195 194 ENDIF … … 207 206 zftu(:,:,:) = 0._wp 208 207 zftv(:,:,:) = 0._wp 209 ! 208 ! 210 209 DO jk = 1, jpkm1 !== before lateral T & S gradients at T-level jk ==! 211 210 DO jj = 1, jpjm1 … … 216 215 END DO 217 216 END DO 218 IF( ln_zps ) THEN! partial steps: correction at the last level217 IF( ln_zps.and.l_grad_zps ) THEN ! partial steps: correction at the last level 219 218 # if defined key_vectopt_loop 220 219 DO jj = 1, 1 … … 224 223 DO ji = 1, jpim1 225 224 # endif 226 zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn) 227 zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 225 zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn) 226 zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 228 227 END DO 229 228 END DO … … 237 236 ! 238 237 ! !== Vertical tracer gradient at level jk and jk+1 239 zdkt (:,:,1) = ( ptb(:,:,jk,jn) - ptb(:,:,jk+1,jn) ) * tmask(:,:,jk+1)238 zdkt3d(:,:,1) = ( ptb(:,:,jk,jn) - ptb(:,:,jk+1,jn) ) * tmask(:,:,jk+1) 240 239 ! 241 ! ! surface boundary condition: zdkt(jk=1)=zdkt(jk=2)242 IF( jk == 1 ) THEN ; zdkt (:,:,0) = zdkt(:,:,1)243 ELSE ; zdkt (:,:,0) = ( ptb(:,:,jk-1,jn) - ptb(:,:,jk,jn) ) * tmask(:,:,jk)240 ! ! surface boundary condition: zdkt3d(jk=0)=zdkt3d(jk=1) 241 IF( jk == 1 ) THEN ; zdkt3d(:,:,0) = zdkt3d(:,:,1) 242 ELSE ; zdkt3d(:,:,0) = ( ptb(:,:,jk-1,jn) - ptb(:,:,jk,jn) ) * tmask(:,:,jk) 244 243 ENDIF 245 244 246 IF( .NOT. l_triad_iso ) THEN 247 triadi = triadi_g 248 triadj = triadj_g 249 ENDIF 250 251 DO ip = 0, 1 !== Horizontal & vertical fluxes 252 DO kp = 0, 1 253 DO jj = 1, jpjm1 254 DO ji = 1, fs_jpim1 255 ze1ur = 1._wp / e1u(ji,jj) 256 zdxt = zdit(ji,jj,jk) * ze1ur 257 ze3wr = 1._wp / fse3w(ji+ip,jj,jk+kp) 258 zdzt = zdkt(ji+ip,jj,kp) * ze3wr 259 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 260 zslope_iso = triadi(ji+ip,jj,jk,1-ip,kp) 261 262 zbu = 0.25_wp * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) 263 zah = fsahtu(ji,jj,jk) !*umask(ji,jj,jk+kp) !fsaht(ji+ip,jj,jk) ===>> ???? 264 zah_slp = zah * zslope_iso 265 zaei_slp = fsaeiw(ji+ip,jj,jk) * zslope_skew !fsaeit(ji+ip,jj,jk)*zslope_skew 266 zftu(ji,jj,jk) = zftu(ji,jj,jk) - ( zah * zdxt + (zah_slp - zaei_slp) * zdzt ) * zbu * ze1ur 267 ztfw(ji+ip,jj,jk+kp) = ztfw(ji+ip,jj,jk+kp) - (zah_slp + zaei_slp) * zdxt * zbu * ze3wr 245 246 IF (ln_botmix_grif) THEN 247 DO ip = 0, 1 !== Horizontal & vertical fluxes 248 DO kp = 0, 1 249 DO jj = 1, jpjm1 250 DO ji = 1, fs_jpim1 251 ze1ur = 1._wp / e1u(ji,jj) 252 zdxt = zdit(ji,jj,jk) * ze1ur 253 ze3wr = 1._wp / fse3w(ji+ip,jj,jk+kp) 254 zdzt = zdkt3d(ji+ip,jj,kp) * ze3wr 255 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 256 zslope_iso = triadi(ji+ip,jj,jk,1-ip,kp) 257 258 zbu = 0.25_wp * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) 259 ! ln_botmix_grif is .T. don't mask zah for bottom half cells 260 zah = fsahtu(ji,jj,jk) !*umask(ji,jj,jk+kp) !fsaht(ji+ip,jj,jk) ===>> ???? 261 zah_slp = zah * zslope_iso 262 zaei_slp = fsaeiw(ji+ip,jj,jk) * zslope_skew !fsaeit(ji+ip,jj,jk)*zslope_skew 263 zftu(ji,jj,jk) = zftu(ji,jj,jk) - ( zah * zdxt + (zah_slp - zaei_slp) * zdzt ) * zbu * ze1ur 264 ztfw(ji+ip,jj,jk+kp) = ztfw(ji+ip,jj,jk+kp) - (zah_slp + zaei_slp) * zdxt * zbu * ze3wr 265 END DO 268 266 END DO 269 267 END DO 270 268 END DO 271 END DO 272 273 DO jp = 0, 1 274 DO kp = 0, 1 275 DO jj = 1, jpjm1 276 DO ji = 1, fs_jpim1 277 ze2vr = 1._wp / e2v(ji,jj) 278 zdyt = zdjt(ji,jj,jk) * ze2vr 279 ze3wr = 1._wp / fse3w(ji,jj+jp,jk+kp) 280 zdzt = zdkt(ji,jj+jp,kp) * ze3wr 281 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 282 zslope_iso = triadj(ji,jj+jp,jk,1-jp,kp) 283 zbv = 0.25_wp * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) 284 zah = fsahtv(ji,jj,jk) !*vmask(ji,jj,jk+kp) !fsaht(ji,jj+jp,jk) 285 zah_slp = zah * zslope_iso 286 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew !fsaeit(ji,jj+jp,jk)*zslope_skew 287 zftv(ji,jj,jk) = zftv(ji,jj,jk) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr 288 ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - (zah_slp + zaei_slp) * zdyt * zbv * ze3wr 269 270 DO jp = 0, 1 271 DO kp = 0, 1 272 DO jj = 1, jpjm1 273 DO ji = 1, fs_jpim1 274 ze2vr = 1._wp / e2v(ji,jj) 275 zdyt = zdjt(ji,jj,jk) * ze2vr 276 ze3wr = 1._wp / fse3w(ji,jj+jp,jk+kp) 277 zdzt = zdkt3d(ji,jj+jp,kp) * ze3wr 278 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 279 zslope_iso = triadj(ji,jj+jp,jk,1-jp,kp) 280 zbv = 0.25_wp * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) 281 ! ln_botmix_grif is .T. don't mask zah for bottom half cells 282 zah = fsahtv(ji,jj,jk) !*vmask(ji,jj,jk+kp) ! fsaht(ji,jj+jp,jk) 283 zah_slp = zah * zslope_iso 284 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew ! fsaeit(ji,jj+jp,jk)*zslope_skew 285 zftv(ji,jj,jk) = zftv(ji,jj,jk) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr 286 ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - (zah_slp + zaei_slp) * zdyt * zbv * ze3wr 287 END DO 289 288 END DO 290 289 END DO 291 290 END DO 292 END DO 293 294 ! !== divergence and add to the general trend ==! 291 ELSE 292 DO ip = 0, 1 !== Horizontal & vertical fluxes 293 DO kp = 0, 1 294 DO jj = 1, jpjm1 295 DO ji = 1, fs_jpim1 296 ze1ur = 1._wp / e1u(ji,jj) 297 zdxt = zdit(ji,jj,jk) * ze1ur 298 ze3wr = 1._wp / fse3w(ji+ip,jj,jk+kp) 299 zdzt = zdkt3d(ji+ip,jj,kp) * ze3wr 300 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 301 zslope_iso = triadi(ji+ip,jj,jk,1-ip,kp) 302 303 zbu = 0.25_wp * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) 304 ! ln_botmix_grif is .F. mask zah for bottom half cells 305 zah = fsahtu(ji,jj,jk) * umask(ji,jj,jk+kp) ! fsaht(ji+ip,jj,jk) ===>> ???? 306 zah_slp = zah * zslope_iso 307 zaei_slp = fsaeiw(ji+ip,jj,jk) * zslope_skew ! fsaeit(ji+ip,jj,jk)*zslope_skew 308 zftu(ji,jj,jk) = zftu(ji,jj,jk) - ( zah * zdxt + (zah_slp - zaei_slp) * zdzt ) * zbu * ze1ur 309 ztfw(ji+ip,jj,jk+kp) = ztfw(ji+ip,jj,jk+kp) - (zah_slp + zaei_slp) * zdxt * zbu * ze3wr 310 END DO 311 END DO 312 END DO 313 END DO 314 315 DO jp = 0, 1 316 DO kp = 0, 1 317 DO jj = 1, jpjm1 318 DO ji = 1, fs_jpim1 319 ze2vr = 1._wp / e2v(ji,jj) 320 zdyt = zdjt(ji,jj,jk) * ze2vr 321 ze3wr = 1._wp / fse3w(ji,jj+jp,jk+kp) 322 zdzt = zdkt3d(ji,jj+jp,kp) * ze3wr 323 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 324 zslope_iso = triadj(ji,jj+jp,jk,1-jp,kp) 325 zbv = 0.25_wp * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) 326 ! ln_botmix_grif is .F. mask zah for bottom half cells 327 zah = fsahtv(ji,jj,jk) * vmask(ji,jj,jk+kp) ! fsaht(ji,jj+jp,jk) 328 zah_slp = zah * zslope_iso 329 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew ! fsaeit(ji,jj+jp,jk)*zslope_skew 330 zftv(ji,jj,jk) = zftv(ji,jj,jk) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr 331 ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - (zah_slp + zaei_slp) * zdyt * zbv * ze3wr 332 END DO 333 END DO 334 END DO 335 END DO 336 END IF 337 ! !== divergence and add to the general trend ==! 295 338 DO jj = 2 , jpjm1 296 DO ji = fs_2, fs_jpim1 339 DO ji = fs_2, fs_jpim1 ! vector opt. 297 340 zbtr = 1._wp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 298 341 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + zbtr * ( zftu(ji-1,jj,jk) - zftu(ji,jj,jk) & … … 303 346 END DO 304 347 ! 305 DO jk = 1, jpkm1 !== Divergence of vertical fluxes added to the general tracer trend348 DO jk = 1, jpkm1 !== Divergence of vertical fluxes added to the general tracer trend 306 349 DO jj = 2, jpjm1 307 DO ji = fs_2, fs_jpim1 350 DO ji = fs_2, fs_jpim1 ! vector opt. 308 351 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ( ztfw(ji,jj,jk+1) - ztfw(ji,jj,jk) ) & 309 352 & / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) … … 312 355 END DO 313 356 ! 314 ! ! "Poleward" diffusive heat or salt transports (T-S case only)357 ! ! "Poleward" diffusive heat or salt transports (T-S case only) 315 358 IF( cdtype == 'TRA' .AND. ln_diaptr .AND. ( MOD( kt, nn_fptr ) == 0 ) ) THEN 316 359 IF( jn == jp_tem) htr_ldf(:) = ptr_vj( zftv(:,:,:) ) ! 3.3 names … … 320 363 #if defined key_diaar5 321 364 IF( cdtype == 'TRA' .AND. jn == jp_tem ) THEN 322 z2d(:,:) = 0._wp 323 zztmp = rau0 * rcp 365 z2d(:,:) = 0._wp 366 zztmp = rau0 * rcp 324 367 DO jk = 1, jpkm1 325 368 DO jj = 2, jpjm1 326 369 DO ji = fs_2, fs_jpim1 ! vector opt. 327 z2d(ji,jj) = z2d(ji,jj) + zftu(ji,jj,jk) 370 z2d(ji,jj) = z2d(ji,jj) + zftu(ji,jj,jk) 328 371 END DO 329 372 END DO … … 332 375 CALL lbc_lnk( z2d, 'U', -1. ) 333 376 CALL iom_put( "udiff_heattr", z2d ) ! heat transport in i-direction 334 z2d(:,:) = 0._wp 377 z2d(:,:) = 0._wp 335 378 DO jk = 1, jpkm1 336 379 DO jj = 2, jpjm1 337 380 DO ji = fs_2, fs_jpim1 ! vector opt. 338 z2d(ji,jj) = z2d(ji,jj) + zftv(ji,jj,jk) 381 z2d(ji,jj) = z2d(ji,jj) + zftv(ji,jj,jk) 339 382 END DO 340 383 END DO … … 342 385 z2d(:,:) = zztmp * z2d(:,:) 343 386 CALL lbc_lnk( z2d, 'V', -1. ) 344 CALL iom_put( "vdiff_heattr", z2d ) ! heat transport in i-direction387 CALL iom_put( "vdiff_heattr", z2d ) ! heat transport in j-direction 345 388 END IF 346 389 #endif -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/NEMO/TOP_SRC/TRP/trcadv.F90
r2715 r2980 103 103 104 104 ! ! add the eiv transport (if necessary) 105 IF( lk_traldf_eiv ) CALL tra_adv_eiv( kt, zun, zvn, zwn, 'TRC' ) 105 IF( lk_traldf_eiv .AND. .NOT. ln_traldf_grif ) & 106 & CALL tra_adv_eiv( kt, zun, zvn, zwn, 'TRC' ) ! add the eiv transport (if necessary) 106 107 ! 107 108 SELECT CASE ( nadv ) !== compute advection trend and add it to general trend ==! -
branches/2011/dev_NOC_2011_MERGE/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90
r2715 r2980 2 2 !!====================================================================== 3 3 !! *** MODULE trcldf *** 4 !! Ocean Passive tracers : lateral diffusive trends 4 !! Ocean Passive tracers : lateral diffusive trends 5 5 !!===================================================================== 6 6 !! History : 9.0 ! 2005-11 (G. Madec) Original code 7 !! NEMO 3.0 ! 2008-01 (C. Ethe, G. Madec) merge TRC-TRA 7 !! NEMO 3.0 ! 2008-01 (C. Ethe, G. Madec) merge TRC-TRA 8 8 !!---------------------------------------------------------------------- 9 9 #if defined key_top … … 23 23 USE traldf_bilap ! lateral mixing (tra_ldf_bilap routine) 24 24 USE traldf_iso ! lateral mixing (tra_ldf_iso routine) 25 USE traldf_iso_grif ! lateral mixing (tra_ldf_iso_grif routine) 25 26 USE traldf_lap ! lateral mixing (tra_ldf_lap routine) 26 27 USE trdmod_oce … … 31 32 PRIVATE 32 33 33 PUBLIC trc_ldf ! called by step.F90 34 PUBLIC trc_ldf ! called by step.F90 34 35 ! !!: ** lateral mixing namelist (nam_trcldf) ** 35 36 INTEGER :: nldf = 0 ! type of lateral diffusion used defined from ln_trcldf_... namlist logicals) … … 39 40 !!---------------------------------------------------------------------- 40 41 !! NEMO/TOP 3.3 , NEMO Consortium (2010) 41 !! $Id$ 42 !! $Id$ 42 43 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 43 44 !!---------------------------------------------------------------------- … … 48 49 !!---------------------------------------------------------------------- 49 50 !! *** ROUTINE tra_ldf *** 50 !! 51 !! 51 52 !! ** Purpose : compute the lateral ocean tracer physics. 52 53 !! … … 61 62 IF( kt == nit000 ) CALL ldf_ctl ! initialisation & control of options 62 63 63 IF( l_trdtrc ) THEN 64 IF( l_trdtrc ) THEN 64 65 ALLOCATE( ztrtrd(jpi,jpj,jpk,jptra) ) ! temporary save of trends 65 66 ztrtrd(:,:,:,:) = tra(:,:,:,:) … … 68 69 SELECT CASE ( nldf ) ! compute lateral mixing trend and add it to the general trend 69 70 CASE ( 0 ) ; CALL tra_ldf_lap ( kt, 'TRC', gtru, gtrv, trb, tra, jptra ) ! iso-level laplacian 70 CASE ( 1 ) ; CALL tra_ldf_iso ( kt, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 ) ! rotated laplacian 71 CASE ( 1 ) 72 IF( ln_traldf_grif ) THEN 73 CALL tra_ldf_iso_grif( kt, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 ) 74 ELSE 75 CALL tra_ldf_iso ( kt, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 ) 76 ENDIF 71 77 CASE ( 2 ) ; CALL tra_ldf_bilap ( kt, 'TRC', gtru, gtrv, trb, tra, jptra ) ! iso-level bilaplacian 72 78 CASE ( 3 ) ; CALL tra_ldf_bilapg( kt, 'TRC', trb, tra, jptra ) ! s-coord. horizontal bilaplacian … … 76 82 WRITE(charout, FMT="('ldf0 ')") ; CALL prt_ctl_trc_info(charout) 77 83 CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 78 CALL tra_ldf_iso ( kt, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 ) 84 IF( ln_traldf_grif ) THEN 85 CALL tra_ldf_iso_grif( kt, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 ) 86 ELSE 87 CALL tra_ldf_iso ( kt, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 ) 88 ENDIF 79 89 WRITE(charout, FMT="('ldf1 ')") ; CALL prt_ctl_trc_info(charout) 80 90 CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) … … 92 102 CALL trd_tra( kt, 'TRC', jn, jptra_trd_ldf, ztrtrd(:,:,:,jn) ) 93 103 END DO 94 DEALLOCATE( ztrtrd ) 104 DEALLOCATE( ztrtrd ) 95 105 ENDIF 96 106 ! ! print mean trends (used for debugging) … … 106 116 !!---------------------------------------------------------------------- 107 117 !! *** ROUTINE ldf_ctl *** 108 !! 118 !! 109 119 !! ** Purpose : Choice of the operator for the lateral tracer diffusion 110 120 !! 111 121 !! ** Method : set nldf from the namtra_ldf logicals 112 !! nldf == -2 No lateral diffusion 122 !! nldf == -2 No lateral diffusion 113 123 !! nldf == -1 ESOPA test: ALL operators are used 114 124 !! nldf == 0 laplacian operator … … 117 127 !! nldf == 3 Rotated bilaplacian 118 128 !!---------------------------------------------------------------------- 119 INTEGER :: ioptio, ierr ! temporary integers 129 INTEGER :: ioptio, ierr ! temporary integers 120 130 !!---------------------------------------------------------------------- 121 131 122 132 ! Define the lateral mixing oparator for tracers 123 133 ! =============================================== 124 134 125 135 ! ! control the input 126 136 ioptio = 0 … … 163 173 ENDIF 164 174 IF ( ln_zps ) THEN ! z-coordinate 165 IF ( ln_trcldf_level ) ierr = 1 ! iso-level not allowed 175 IF ( ln_trcldf_level ) ierr = 1 ! iso-level not allowed 166 176 IF ( ln_trcldf_hor ) nldf = 2 ! horizontal (no rotation) 167 177 IF ( ln_trcldf_iso ) ierr = 2 ! isoneutral ( rotation)
Note: See TracChangeset
for help on using the changeset viewer.
