Changeset 14037 for NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/LDF/ldftra.F90
- Timestamp:
- 2020-12-03T12:20:38+01:00 (3 years ago)
- Location:
- NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG
- Property svn:externals
-
old new 8 8 9 9 # SETTE 10 ^/utils/CI/sette @13292sette10 ^/utils/CI/sette_wave@13990 sette
-
- Property svn:externals
-
NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/LDF/ldftra.F90
r13295 r14037 246 246 ENDIF 247 247 ! 248 IF( ln_ldfeiv .AND. .NOT.( ln_traldf_iso .OR. ln_traldf_triad ) ) & 249 & CALL ctl_stop( 'ln_ldfeiv=T requires iso-neutral laplacian diffusion' ) 250 IF( ln_isfcav .AND. ln_traldf_triad ) & 251 & CALL ctl_stop( ' ice shelf cavity and traldf_triad not tested' ) 248 IF( ln_isfcav .AND. ln_traldf_triad ) CALL ctl_stop( ' ice shelf cavity and traldf_triad not tested' ) 252 249 ! 253 250 IF( nldf_tra == np_lap_i .OR. nldf_tra == np_lap_it .OR. & … … 430 427 zaht_min = 0.2_wp * aht0 ! minimum value for aht 431 428 zDaht = aht0 - zaht_min 432 DO_2D( 1, 1, 1, 1 ) 429 ! NOTE: [tiling-comms-merge] Change needed to preserve results with respect to the trunk 430 DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 433 431 !!gm CAUTION : here we assume lat/lon grid in 20deg N/S band (like all ORCA cfg) 434 432 !! ==>>> The Coriolis value is identical for t- & u_points, and for v- and f-points … … 541 539 IF( ln_traldf_blp ) CALL ctl_stop( 'ldf_eiv_init: eddy induced velocity ONLY with laplacian diffusivity' ) 542 540 ! 541 IF( .NOT.( ln_traldf_iso .OR. ln_traldf_triad ) ) & 542 & CALL ctl_stop( 'ln_ldfeiv=T requires iso-neutral laplacian diffusion' ) 543 543 ! != allocate the aei arrays 544 544 ALLOCATE( aeiu(jpi,jpj,jpk), aeiv(jpi,jpj,jpk), STAT=ierr ) … … 694 694 CALL lbc_lnk( 'ldftra', zaeiw(:,:), 'W', 1.0_wp ) ! lateral boundary condition 695 695 ! 696 DO_2D( 0, 0, 0, 0 ) 696 DO_2D( 0, 0, 0, 0 ) !== aei at u- and v-points ==! 697 697 paeiu(ji,jj,1) = 0.5_wp * ( zaeiw(ji,jj) + zaeiw(ji+1,jj ) ) * umask(ji,jj,1) 698 698 paeiv(ji,jj,1) = 0.5_wp * ( zaeiw(ji,jj) + zaeiw(ji ,jj+1) ) * vmask(ji,jj,1) … … 726 726 !! ** Action : pu, pv increased by the eiv transport 727 727 !!---------------------------------------------------------------------- 728 INTEGER 729 INTEGER 730 INTEGER 731 CHARACTER(len=3) 732 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pu! in : 3 ocean transport components [m3/s]733 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pv! out: 3 ocean transport components [m3/s]734 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pw! increased by the eiv [m3/s]728 INTEGER , INTENT(in ) :: kt ! ocean time-step index 729 INTEGER , INTENT(in ) :: kit000 ! first time step index 730 INTEGER , INTENT(in ) :: Kmm, Krhs ! ocean time level indices 731 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 732 REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(inout) :: pu ! in : 3 ocean transport components [m3/s] 733 REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(inout) :: pv ! out: 3 ocean transport components [m3/s] 734 REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(inout) :: pw ! increased by the eiv [m3/s] 735 735 !! 736 736 INTEGER :: ji, jj, jk ! dummy loop indices 737 737 REAL(wp) :: zuwk, zuwk1, zuwi, zuwi1 ! local scalars 738 738 REAL(wp) :: zvwk, zvwk1, zvwj, zvwj1 ! - - 739 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpsi_uw, zpsi_vw 740 !!---------------------------------------------------------------------- 741 ! 742 IF( kt == kit000 ) THEN 743 IF(lwp) WRITE(numout,*) 744 IF(lwp) WRITE(numout,*) 'ldf_eiv_trp : eddy induced advection on ', cdtype,' :' 745 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ add to velocity fields the eiv component' 739 REAL(wp), DIMENSION(A2D(nn_hls),jpk) :: zpsi_uw, zpsi_vw 740 !!---------------------------------------------------------------------- 741 ! 742 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 743 IF( kt == kit000 ) THEN 744 IF(lwp) WRITE(numout,*) 745 IF(lwp) WRITE(numout,*) 'ldf_eiv_trp : eddy induced advection on ', cdtype,' :' 746 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ add to velocity fields the eiv component' 747 ENDIF 746 748 ENDIF 747 749 … … 782 784 !! 783 785 !!---------------------------------------------------------------------- 784 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: psi_uw, psi_vw ! streamfunction [m3/s]785 INTEGER 786 REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(inout) :: psi_uw, psi_vw ! streamfunction [m3/s] 787 INTEGER , INTENT(in ) :: Kmm ! ocean time level indices 786 788 ! 787 789 INTEGER :: ji, jj, jk ! dummy loop indices 788 790 REAL(wp) :: zztmp ! local scalar 789 REAL(wp), DIMENSION( jpi,jpj) :: zw2d ! 2D workspace790 REAL(wp), DIMENSION( jpi,jpj,jpk) :: zw3d ! 3D workspace791 REAL(wp), DIMENSION(A2D(nn_hls)) :: zw2d ! 2D workspace 792 REAL(wp), DIMENSION(A2D(nn_hls),jpk) :: zw3d ! 3D workspace 791 793 !!---------------------------------------------------------------------- 792 794 ! … … 794 796 !!gm to be redesigned.... 795 797 ! !== eiv stream function: output ==! 796 CALL lbc_lnk_multi( 'ldftra', psi_uw, 'U', -1.0_wp , psi_vw, 'V', -1.0_wp )797 !798 798 !!gm CALL iom_put( "psi_eiv_uw", psi_uw ) ! output 799 799 !!gm CALL iom_put( "psi_eiv_vw", psi_vw ) … … 803 803 zw3d(:,:,jpk) = 0._wp ! bottom value always 0 804 804 ! 805 DO jk = 1, jpkm1! e2u e3u u_eiv = -dk[psi_uw]806 zw3d( :,:,jk) = ( psi_uw(:,:,jk+1) - psi_uw(:,:,jk) ) / ( e2u(:,:) * e3u(:,:,jk,Kmm) )807 END DO805 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! e2u e3u u_eiv = -dk[psi_uw] 806 zw3d(ji,jj,jk) = ( psi_uw(ji,jj,jk+1) - psi_uw(ji,jj,jk) ) / ( e2u(ji,jj) * e3u(ji,jj,jk,Kmm) ) 807 END_3D 808 808 CALL iom_put( "uoce_eiv", zw3d ) 809 809 ! 810 DO jk = 1, jpkm1! e1v e3v v_eiv = -dk[psi_vw]811 zw3d( :,:,jk) = ( psi_vw(:,:,jk+1) - psi_vw(:,:,jk) ) / ( e1v(:,:) * e3v(:,:,jk,Kmm) )812 END DO810 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! e1v e3v v_eiv = -dk[psi_vw] 811 zw3d(ji,jj,jk) = ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj,jk) ) / ( e1v(ji,jj) * e3v(ji,jj,jk,Kmm) ) 812 END_3D 813 813 CALL iom_put( "voce_eiv", zw3d ) 814 814 ! 815 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 815 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! e1 e2 w_eiv = dk[psix] + dk[psix] 816 816 zw3d(ji,jj,jk) = ( psi_vw(ji,jj,jk) - psi_vw(ji ,jj-1,jk) & 817 817 & + psi_uw(ji,jj,jk) - psi_uw(ji-1,jj ,jk) ) / e1e2t(ji,jj) 818 818 END_3D 819 CALL lbc_lnk( 'ldftra', zw3d, 'T', 1.0_wp ) ! lateral boundary condition820 819 CALL iom_put( "woce_eiv", zw3d ) 821 820 ! 822 821 IF( iom_use('weiv_masstr') ) THEN ! vertical mass transport & its square value 823 zw2d(:,:) = rho0 * e1e2t(:,:) 822 DO_2D( 0, 0, 0, 0 ) 823 zw2d(ji,jj) = rho0 * e1e2t(ji,jj) 824 END_2D 824 825 DO jk = 1, jpk 825 826 zw3d(:,:,jk) = zw3d(:,:,jk) * zw2d(:,:) … … 845 846 zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) 846 847 END_3D 847 CALL lbc_lnk( 'ldftra', zw2d, 'U', -1.0_wp )848 CALL lbc_lnk( 'ldftra', zw3d, 'U', -1.0_wp )849 848 CALL iom_put( "ueiv_heattr" , zztmp * zw2d ) ! heat transport in i-direction 850 849 CALL iom_put( "ueiv_heattr3d", zztmp * zw3d ) ! heat transport in i-direction … … 866 865 zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) 867 866 END_3D 868 CALL lbc_lnk( 'ldftra', zw2d, 'V', -1.0_wp ) 869 CALL iom_put( "veiv_heattr", zztmp * zw2d ) ! heat transport in j-direction 870 CALL iom_put( "veiv_heattr", zztmp * zw3d ) ! heat transport in j-direction 867 CALL iom_put( "veiv_heattr" , zztmp * zw2d ) ! heat transport in j-direction 868 CALL iom_put( "veiv_heattr3d", zztmp * zw3d ) ! heat transport in j-direction 871 869 ! 872 870 IF( iom_use( 'sophteiv' ) ) CALL dia_ptr_hst( jp_tem, 'eiv', 0.5 * zw3d ) … … 881 879 zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) 882 880 END_3D 883 CALL lbc_lnk( 'ldftra', zw2d, 'U', -1.0_wp )884 CALL lbc_lnk( 'ldftra', zw3d, 'U', -1.0_wp )885 881 CALL iom_put( "ueiv_salttr", zztmp * zw2d ) ! salt transport in i-direction 886 882 CALL iom_put( "ueiv_salttr3d", zztmp * zw3d ) ! salt transport in i-direction … … 893 889 zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) 894 890 END_3D 895 CALL lbc_lnk( 'ldftra', zw2d, 'V', -1.0_wp ) 896 CALL iom_put( "veiv_salttr", zztmp * zw2d ) ! salt transport in j-direction 897 CALL iom_put( "veiv_salttr", zztmp * zw3d ) ! salt transport in j-direction 891 CALL iom_put( "veiv_salttr" , zztmp * zw2d ) ! salt transport in j-direction 892 CALL iom_put( "veiv_salttr3d", zztmp * zw3d ) ! salt transport in j-direction 898 893 ! 899 894 IF( iom_use( 'sopsteiv' ) ) CALL dia_ptr_hst( jp_sal, 'eiv', 0.5 * zw3d )
Note: See TracChangeset
for help on using the changeset viewer.