Changeset 2969
- Timestamp:
- 2011-10-20T18:00:23+02:00 (13 years ago)
- Location:
- branches/2011/dev_r2782_NOCS_Griffies
- Files:
-
- 9 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) (18 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/traldf_iso_grif.F90 (modified) (19 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/dev_r2782_NOCS_Griffies/DOC/TexFiles/Namelist/namtra_ldf
r2540 r2969 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_r2782_NOCS_Griffies/NEMOGCM/CONFIG/GYRE/EXP00/namelist
r2715 r2969 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_r2782_NOCS_Griffies/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist
r2715 r2969 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_r2782_NOCS_Griffies/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist
r2715 r2969 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_r2782_NOCS_Griffies/NEMOGCM/CONFIG/POMME/EXP00/namelist
r2650 r2969 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_r2782_NOCS_Griffies/NEMOGCM/NEMO/OPA_SRC/LDF/ldfslp.F90
r2840 r2969 402 402 !!---------------------------------------------------------------------- 403 403 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 404 USE oce , ONLY: zalbet => ua ! use ua as workspace404 USE oce , ONLY: zalbet => ua ! use ua as workspace 405 405 USE wrk_nemo, ONLY: z1_mlbw => wrk_2d_1 406 406 !! 407 INTEGER, INTENT( in ) :: kt ! ocean time-step index407 INTEGER, INTENT( in ) :: kt ! ocean time-step index 408 408 !! 409 409 INTEGER :: ji, jj, jk, jl, ip, jp, kp ! dummy loop indices … … 425 425 !--------------------------------! 426 426 ! 427 CALL eos_alpbet( tsb, zalbet, zbeta0 ) 428 ! 429 DO jl = 0, 1 !== unmasked before density i- j-, k-gradients ==!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 430 ! 431 ip = jl ; jp = jl ! guaranteed nonzero gradients ( absolute value larger than repsln)432 DO jk = 1, jpkm1 433 DO jj = 1, jpjm1 434 DO ji = 1, fs_jpim1 ! vector opt.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 435 zdit = ( tsb(ji+1,jj,jk,jp_tem) - tsb(ji,jj,jk,jp_tem) ) ! i-gradient of T & S at u-point 436 436 zdis = ( tsb(ji+1,jj,jk,jp_sal) - tsb(ji,jj,jk,jp_sal) ) … … 445 445 END DO 446 446 ! 447 IF( ln_zps.and.l n_grad_zps ) THEN! partial steps: correction of i- & j-grad on bottom447 IF( ln_zps.and.l_grad_zps ) THEN ! partial steps: correction of i- & j-grad on bottom 448 448 # if defined key_vectopt_loop 449 449 DO jj = 1, 1 450 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling)450 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 451 451 # else 452 452 DO jj = 1, jpjm1 … … 466 466 END DO 467 467 468 DO kp = 0, 1 !== unmasked before density i- j-, k-gradients ==!469 DO jk = 1, jpkm1 470 DO jj = 1, jpj 471 DO ji = 1, jpi ! vector opt.472 IF( jk+kp > 1 ) THEN 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 473 zdkt = ( tsb(ji,jj,jk+kp-1,jp_tem) - tsb(ji,jj,jk+kp,jp_tem) ) 474 474 zdks = ( tsb(ji,jj,jk+kp-1,jp_sal) - tsb(ji,jj,jk+kp,jp_sal) ) … … 484 484 END DO 485 485 ! 486 DO jj = 1, jpj !== Reciprocal depth of the w-point below ML base ==!486 DO jj = 1, jpj !== Reciprocal depth of the w-point below ML base ==! 487 487 DO ji = 1, jpi 488 488 jk = MIN( nmln(ji,jj), mbkt(ji,jj) ) + 1 ! MIN in case ML depth is the ocean depth … … 491 491 END DO 492 492 ! 493 ! !== intialisations to zero ==!494 ! 495 wslp2 (:,:,:) = 0._wp 493 ! !== intialisations to zero ==! 494 ! 495 wslp2 (:,:,:) = 0._wp ! wslp2 will be cumulated 3D field set to zero 496 496 triadi_g(:,:,1,:,:) = 0._wp ; triadi_g(:,:,jpk,:,:) = 0._wp ! set surface and bottom slope to zero 497 497 triadj_g(:,:,1,:,:) = 0._wp ; triadj_g(:,:,jpk,:,:) = 0._wp 498 !!gm _iso set to zero missing498 !!gm _iso set to zero missing 499 499 triadi (:,:,1,:,:) = 0._wp ; triadj (:,:,jpk,:,:) = 0._wp ! set surface and bottom slope to zero 500 500 triadj (:,:,1,:,:) = 0._wp ; triadj (:,:,jpk,:,:) = 0._wp … … 504 504 !-------------------------------------! 505 505 ! 506 DO jl = 0, 1 ! calculate slope of the 4 triads immediately ONE level below mixed-layer base507 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 508 508 DO jj = 1, jpjm1 509 509 DO ji = 1, fs_jpim1 … … 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 … … 575 575 zti_lim = ( zti_g_lim + zti_coord ) * umask(ji,jj,jk+kp) ! remove coordinate slope => relative to coordinate surfaces 576 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 <<==577 zti_lim2 = zti_lim * zti_lim ! square of limited slopes ! masked <<== 578 578 ztj_lim2 = ztj_lim * ztj_lim 579 579 ! … … 596 596 zbtj = e1t(ji,jj+jp) * e2t(ji,jj+jp) * fse3w(ji,jj+jp,jk+kp) 597 597 ! 598 !!gm this may inhibit vectorization on Vect Computers, and even on scalar computers.... ==> to be checked599 ! agn may need to change to using ztj_g_lim**2, as wslp2 just used for eddy growth rate, needs *real* slope598 !!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 600 wslp2 (ji+ip,jj,jk+kp) = wslp2(ji+ip,jj,jk+kp) + 0.25_wp * zbu / zbti * zti_lim2 ! masked 601 601 wslp2 (ji,jj+jp,jk+kp) = wslp2(ji,jj+jp,jk+kp) + 0.25_wp * zbv / zbtj * ztj_lim2 … … 637 637 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: p_dzr ! z-gradient of density (T-point) 638 638 !! 639 INTEGER :: ji , jj , jk ! dummy loop indices640 INTEGER :: iku, ikv, ik, ikm1 ! local integers639 INTEGER :: ji , jj , jk ! dummy loop indices 640 INTEGER :: iku, ikv, ik, ikm1 ! local integers 641 641 REAL(wp) :: zeps, zm1_g, zm1_2g ! local scalars 642 642 REAL(wp) :: zci, zfi, zau, zbu, zai, zbi ! - - … … 654 654 wslpjml(1,:) = 0._wp ; wslpjml(jpi,:) = 0._wp 655 655 ! 656 ! !== surface mixed layer mask !657 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 658 658 # if defined key_vectopt_loop 659 659 DO jj = 1, 1 660 DO ji = 1, jpij ! vector opt. (forced unrolling)660 DO ji = 1, jpij ! vector opt. (forced unrolling) 661 661 # else 662 662 DO jj = 1, jpj … … 689 689 DO ji = 2, jpim1 690 690 # endif 691 ! !== Slope at u- & v-points just below the Mixed Layer ==!691 ! !== Slope at u- & v-points just below the Mixed Layer ==! 692 692 ! 693 ! 693 ! !- vertical density gradient for u- and v-slopes (from dzr at T-point) 694 694 iku = MIN( MAX( 1, nmln(ji,jj) , nmln(ji+1,jj) ) , jpkm1 ) ! ML (MAX of T-pts, bound by jpkm1) 695 695 ikv = MIN( MAX( 1, nmln(ji,jj) , nmln(ji,jj+1) ) , jpkm1 ) ! 696 696 zbu = 0.5_wp * ( p_dzr(ji,jj,iku) + p_dzr(ji+1,jj ,iku) ) 697 697 zbv = 0.5_wp * ( p_dzr(ji,jj,ikv) + p_dzr(ji ,jj+1,ikv) ) 698 ! 698 ! !- horizontal density gradient at u- & v-points 699 699 zau = p_gru(ji,jj,iku) / e1u(ji,jj) 700 700 zav = p_grv(ji,jj,ikv) / e2v(ji,jj) 701 ! 702 ! 701 ! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0 702 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 703 703 zbu = MIN( zbu , -100._wp* ABS( zau ) , -7.e+3_wp/fse3u(ji,jj,iku)* ABS( zau ) ) 704 704 zbv = MIN( zbv , -100._wp* ABS( zav ) , -7.e+3_wp/fse3v(ji,jj,ikv)* ABS( zav ) ) 705 ! 705 ! !- Slope at u- & v-points (uslpml, vslpml) 706 706 uslpml(ji,jj) = zau / ( zbu - zeps ) * umask(ji,jj,iku) 707 707 vslpml(ji,jj) = zav / ( zbv - zeps ) * vmask(ji,jj,ikv) 708 708 ! 709 ! !== i- & j-slopes at w-points just below the Mixed Layer ==!709 ! !== i- & j-slopes at w-points just below the Mixed Layer ==! 710 710 ! 711 711 ik = MIN( nmln(ji,jj) + 1, jpk ) 712 712 ikm1 = MAX( 1, ik-1 ) 713 ! 713 ! !- vertical density gradient for w-slope (from N^2) 714 714 zbw = zm1_2g * pn2 (ji,jj,ik) * ( prd (ji,jj,ik) + prd (ji,jj,ikm1) + 2. ) 715 ! 715 ! !- horizontal density i- & j-gradient at w-points 716 716 zci = MAX( umask(ji-1,jj,ik ) + umask(ji,jj,ik ) & 717 717 & + umask(ji-1,jj,ikm1) + umask(ji,jj,ikm1) , zeps ) * e1t(ji,jj) … … 722 722 zaj = ( p_grv(ji,jj-1,ik ) + p_grv(ji,jj,ik ) & 723 723 & + p_grv(ji,jj-1,ikm1) + p_grv(ji,jj,ikm1) ) / zcj * tmask(ji,jj,ik) 724 ! 725 ! 724 ! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0. 725 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 726 726 zbi = MIN( zbw , -100._wp* ABS( zai ) , -7.e+3_wp/fse3w(ji,jj,ik)* ABS( zai ) ) 727 727 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/fse3w(ji,jj,ik)* ABS( zaj ) ) 728 ! 728 ! !- i- & j-slope at w-points (wslpiml, wslpjml) 729 729 wslpiml(ji,jj) = zai / ( zbi - zeps ) * tmask (ji,jj,ik) 730 730 wslpjml(ji,jj) = zaj / ( zbj - zeps ) * tmask (ji,jj,ik) 731 731 END DO 732 732 END DO 733 !!gm this lbc_lnk should be useless....733 !!gm this lbc_lnk should be useless.... 734 734 CALL lbc_lnk( uslpml , 'U', -1. ) ; CALL lbc_lnk( vslpml , 'V', -1. ) ! lateral boundary cond. (sign change) 735 735 CALL lbc_lnk( wslpiml, 'W', -1. ) ; CALL lbc_lnk( wslpjml, 'W', -1. ) ! lateral boundary conditions … … 765 765 IF( ln_dynldf_iso ) CALL ctl_stop( 'ldf_slp_init: Griffies operator on momentum not supported' ) 766 766 ! 767 ! IF( ( ln_traldf_hor .OR. ln_dynldf_hor ) .AND. ln_sco ) &768 ! CALL ctl_stop( 'ldf_slp_init: horizontal Griffies operator in s-coordinate not supported' )769 !770 767 ELSE ! Madec operator : slopes at u-, v-, and w-points 771 768 ALLOCATE( uslp(jpi,jpj,jpk) , vslp(jpi,jpj,jpk) , wslpi(jpi,jpj,jpk) , wslpj(jpi,jpj,jpk) , & … … 780 777 wslpj(:,:,:) = 0._wp ; wslpjml(:,:) = 0._wp 781 778 782 !!gm I no longer understand this.....779 !!gm I no longer understand this..... 783 780 IF( (ln_traldf_hor .OR. ln_dynldf_hor) .AND. .NOT. (lk_vvl .AND. ln_rstart) ) THEN 784 781 IF(lwp) WRITE(numout,*) ' Horizontal mixing in s-coordinate: slope = slope of s-surfaces' … … 813 810 LOGICAL, PUBLIC, PARAMETER :: lk_ldfslp = .FALSE. !: slopes flag 814 811 CONTAINS 815 SUBROUTINE ldf_slp( kt, prd, pn2 ) 812 SUBROUTINE ldf_slp( kt, prd, pn2 ) ! Dummy routine 816 813 INTEGER, INTENT(in) :: kt 817 814 REAL, DIMENSION(:,:,:), INTENT(in) :: prd, pn2 … … 822 819 WRITE(*,*) 'ldf_slp_grif: You should not have seen this print! error?', kt 823 820 END SUBROUTINE ldf_slp_grif 824 SUBROUTINE ldf_slp_init ! Dummy routine821 SUBROUTINE ldf_slp_init ! Dummy routine 825 822 END SUBROUTINE ldf_slp_init 826 823 #endif -
branches/2011/dev_r2782_NOCS_Griffies/NEMOGCM/NEMO/OPA_SRC/LDF/ldftra.F90
r2840 r2969 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, ln_grad_zps,&69 & ln_triad_iso , ln_botmix_grif, & 70 70 & rn_aht_0 , rn_ahtb_0 , rn_aeiv_0, & 71 71 & rn_slpmax … … 95 95 WRITE(numout,*) ' maximum isoppycnal slope rn_slpmax = ', rn_slpmax 96 96 WRITE(numout,*) ' + griffies operator internal controls not set via the namelist (experimental): ' 97 WRITE(numout,*) ' calculate triads twice ln_triad_iso = ', ln_triad_iso 98 WRITE(numout,*) ' special Tgradts w ptl steps ln_grad_zps = ', ln_grad_zps 99 WRITE(numout,*) ' GM -->lat mixing on bottom ln_botmix = ', ln_botmix 97 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 100 99 WRITE(numout,*) 101 100 ENDIF -
branches/2011/dev_r2782_NOCS_Griffies/NEMOGCM/NEMO/OPA_SRC/LDF/ldftra_oce.F90
r2840 r2969 32 32 33 33 REAL(wp), PUBLIC :: aht0, ahtb0, aeiv0 !!: OLD namelist names 34 LOGICAL , PUBLIC :: ln_triad_iso 35 LOGICAL , PUBLIC :: ln_ grad_zps = .FALSE. !: special treatment for Horz Tgradients w partial steps36 LOGICAL , PUBLIC :: l n_botmix = .FALSE. !: mixing on bottom34 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 37 37 38 38 #if defined key_traldf_c3d -
branches/2011/dev_r2782_NOCS_Griffies/NEMOGCM/NEMO/OPA_SRC/TRA/traldf_iso_grif.F90
r2958 r2969 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 … … 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 in124 ! Fortran 2003 and up. We would be OK if code was written to use125 ! 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 122 129 123 IF( kt == nit000.AND..NOT.ALLOCATED(ah_wslp2) ) THEN … … 131 125 IF(lwp) WRITE(numout,*) 'tra_ldf_iso_grif : rotated laplacian diffusion operator on ', cdtype 132 126 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 133 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 ) 134 128 IF( lk_mpp ) CALL mpp_sum ( ierr ) 135 129 IF( ierr > 0 ) CALL ctl_stop('STOP', 'tra_ldf_iso_grif: unable to allocate arrays') … … 145 139 !!---------------------------------------------------------------------- 146 140 147 !!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 148 142 149 143 ah_wslp2(:,:,:) = 0._wp … … 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 159 ! Subtract s-coordinate slope at t-points to give slope rel to s surfaces 166 ! 160 ! (do this by *adding* gradient of depth) 167 161 zslope2 = zslope_skew + ( fsdept(ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * ze1ur * umask(ji,jj,jk+kp) 168 162 zslope2 = zslope2 *zslope2 … … 170 164 & + zah * ( zbu * ze3wr / ( e1t(ji+ip,jj) * e2t(ji+ip,jj) ) ) * zslope2 171 165 IF( ln_traldf_gdia ) THEN 172 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 173 167 psix_eiv(ji,jj,jk+kp) = psix_eiv(ji,jj,jk+kp) + 0.25_wp * zaei_slp 174 168 ENDIF … … 187 181 ze3wr = 1.0_wp / fse3w(ji,jj+jp,jk+kp) 188 182 zbv = 0.25_wp * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) 189 zah = fsahtv(ji,jj,jk) !fsaht(ji,jj+jp,jk)183 zah = fsahtv(ji,jj,jk) ! fsaht(ji,jj+jp,jk) 190 184 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 191 185 ! Subtract s-coordinate slope at t-points to give slope rel to s surfaces … … 196 190 & + zah * ( zbv * ze3wr / ( e1t(ji,jj+jp) * e2t(ji,jj+jp) ) ) * zslope2 197 191 IF( ln_traldf_gdia ) THEN 198 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 199 193 psiy_eiv(ji,jj,jk+kp) = psiy_eiv(ji,jj,jk+kp) + 0.25_wp * zaei_slp 200 194 ENDIF … … 221 215 END DO 222 216 END DO 223 IF( ln_zps.and.l n_grad_zps ) THEN ! partial steps: correction at the last level217 IF( ln_zps.and.l_grad_zps ) THEN ! partial steps: correction at the last level 224 218 # if defined key_vectopt_loop 225 219 DO jj = 1, 1 … … 242 236 ! 243 237 ! !== Vertical tracer gradient at level jk and jk+1 244 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) 245 239 ! 246 ! ! surface boundary condition: zdkt(jk=1)=zdkt(jk=2)247 IF( jk == 1 ) THEN ; zdkt (:,:,0) = zdkt(:,:,1)248 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) 249 243 ENDIF 250 244 251 245 252 IF (ln_botmix ) THEN246 IF (ln_botmix_grif) THEN 253 247 DO ip = 0, 1 !== Horizontal & vertical fluxes 254 248 DO kp = 0, 1 … … 256 250 DO ji = 1, fs_jpim1 257 251 ze1ur = 1._wp / e1u(ji,jj) 258 zdxt = zdit(ji,jj,jk) * ze1ur252 zdxt = zdit(ji,jj,jk) * ze1ur 259 253 ze3wr = 1._wp / fse3w(ji+ip,jj,jk+kp) 260 zdzt = zdkt (ji+ip,jj,kp) * ze3wr254 zdzt = zdkt3d(ji+ip,jj,kp) * ze3wr 261 255 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 262 256 zslope_iso = triadi(ji+ip,jj,jk,1-ip,kp) 263 257 264 258 zbu = 0.25_wp * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) 265 ! ln_botmix is .T. don't mask zah for bottom half cells259 ! ln_botmix_grif is .T. don't mask zah for bottom half cells 266 260 zah = fsahtu(ji,jj,jk) !*umask(ji,jj,jk+kp) !fsaht(ji+ip,jj,jk) ===>> ???? 267 261 zah_slp = zah * zslope_iso … … 279 273 DO ji = 1, fs_jpim1 280 274 ze2vr = 1._wp / e2v(ji,jj) 281 zdyt = zdjt(ji,jj,jk) * ze2vr275 zdyt = zdjt(ji,jj,jk) * ze2vr 282 276 ze3wr = 1._wp / fse3w(ji,jj+jp,jk+kp) 283 zdzt = zdkt(ji,jj+jp,kp) * ze3wr277 zdzt = zdkt3d(ji,jj+jp,kp) * ze3wr 284 278 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 285 279 zslope_iso = triadj(ji,jj+jp,jk,1-jp,kp) 286 280 zbv = 0.25_wp * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) 287 ! ln_botmix is .T. don't mask zah for bottom half cells288 zah = fsahtv(ji,jj,jk) !*vmask(ji,jj,jk+kp) !fsaht(ji,jj+jp,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) 289 283 zah_slp = zah * zslope_iso 290 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew !fsaeit(ji,jj+jp,jk)*zslope_skew284 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew ! fsaeit(ji,jj+jp,jk)*zslope_skew 291 285 zftv(ji,jj,jk) = zftv(ji,jj,jk) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr 292 286 ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - (zah_slp + zaei_slp) * zdyt * zbv * ze3wr … … 301 295 DO ji = 1, fs_jpim1 302 296 ze1ur = 1._wp / e1u(ji,jj) 303 zdxt = zdit(ji,jj,jk) * ze1ur297 zdxt = zdit(ji,jj,jk) * ze1ur 304 298 ze3wr = 1._wp / fse3w(ji+ip,jj,jk+kp) 305 zdzt = zdkt (ji+ip,jj,kp) * ze3wr299 zdzt = zdkt3d(ji+ip,jj,kp) * ze3wr 306 300 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 307 301 zslope_iso = triadi(ji+ip,jj,jk,1-ip,kp) 308 302 309 303 zbu = 0.25_wp * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) 310 ! ln_botmix is .F. mask zah for bottom half cells311 zah = fsahtu(ji,jj,jk) * umask(ji,jj,jk+kp) ! fsaht(ji+ip,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) ===>> ???? 312 306 zah_slp = zah * zslope_iso 313 zaei_slp = fsaeiw(ji+ip,jj,jk) * zslope_skew !fsaeit(ji+ip,jj,jk)*zslope_skew307 zaei_slp = fsaeiw(ji+ip,jj,jk) * zslope_skew ! fsaeit(ji+ip,jj,jk)*zslope_skew 314 308 zftu(ji,jj,jk) = zftu(ji,jj,jk) - ( zah * zdxt + (zah_slp - zaei_slp) * zdzt ) * zbu * ze1ur 315 309 ztfw(ji+ip,jj,jk+kp) = ztfw(ji+ip,jj,jk+kp) - (zah_slp + zaei_slp) * zdxt * zbu * ze3wr … … 324 318 DO ji = 1, fs_jpim1 325 319 ze2vr = 1._wp / e2v(ji,jj) 326 zdyt = zdjt(ji,jj,jk) * ze2vr320 zdyt = zdjt(ji,jj,jk) * ze2vr 327 321 ze3wr = 1._wp / fse3w(ji,jj+jp,jk+kp) 328 zdzt = zdkt(ji,jj+jp,kp) * ze3wr322 zdzt = zdkt3d(ji,jj+jp,kp) * ze3wr 329 323 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 330 324 zslope_iso = triadj(ji,jj+jp,jk,1-jp,kp) 331 325 zbv = 0.25_wp * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) 332 ! ln_botmix is .F. mask zah for bottom half cells333 zah = fsahtv(ji,jj,jk) * vmask(ji,jj,jk+kp) ! fsaht(ji,jj+jp,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) 334 328 zah_slp = zah * zslope_iso 335 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew !fsaeit(ji,jj+jp,jk)*zslope_skew329 zaei_slp = fsaeiw(ji,jj+jp,jk) * zslope_skew ! fsaeit(ji,jj+jp,jk)*zslope_skew 336 330 zftv(ji,jj,jk) = zftv(ji,jj,jk) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr 337 331 ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - (zah_slp + zaei_slp) * zdyt * zbv * ze3wr … … 341 335 END DO 342 336 END IF 343 ! !== divergence and add to the general trend ==!337 ! !== divergence and add to the general trend ==! 344 338 DO jj = 2 , jpjm1 345 DO ji = fs_2, fs_jpim1 339 DO ji = fs_2, fs_jpim1 ! vector opt. 346 340 zbtr = 1._wp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 347 341 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + zbtr * ( zftu(ji-1,jj,jk) - zftu(ji,jj,jk) & … … 352 346 END DO 353 347 ! 354 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 355 349 DO jj = 2, jpjm1 356 DO ji = fs_2, fs_jpim1 350 DO ji = fs_2, fs_jpim1 ! vector opt. 357 351 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ( ztfw(ji,jj,jk+1) - ztfw(ji,jj,jk) ) & 358 352 & / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) … … 361 355 END DO 362 356 ! 363 ! ! "Poleward" diffusive heat or salt transports (T-S case only)357 ! ! "Poleward" diffusive heat or salt transports (T-S case only) 364 358 IF( cdtype == 'TRA' .AND. ln_diaptr .AND. ( MOD( kt, nn_fptr ) == 0 ) ) THEN 365 359 IF( jn == jp_tem) htr_ldf(:) = ptr_vj( zftv(:,:,:) ) ! 3.3 names … … 391 385 z2d(:,:) = zztmp * z2d(:,:) 392 386 CALL lbc_lnk( z2d, 'V', -1. ) 393 CALL iom_put( "vdiff_heattr", z2d ) ! heat transport in i-direction387 CALL iom_put( "vdiff_heattr", z2d ) ! heat transport in j-direction 394 388 END IF 395 389 #endif
Note: See TracChangeset
for help on using the changeset viewer.
