Changeset 5974 for branches/2015/dev_r5072_UKMO2_OBS_simplification/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
- Timestamp:
- 2015-12-02T11:52:05+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5072_UKMO2_OBS_simplification/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
r5682 r5974 10 10 !! 3.6 ! 2014-11 (P. Mathiot) add ice shelf capability 11 11 !!---------------------------------------------------------------------- 12 !! 'key_vvl' variable volume 13 !!---------------------------------------------------------------------- 12 14 13 !!---------------------------------------------------------------------- 15 14 !! dom_vvl_init : define initial vertical scale factors, depths and column thickness … … 19 18 !! dom_vvl_rst : read/write restart file 20 19 !! dom_vvl_ctl : Check the vvl options 21 !! dom_vvl_orca_fix : Recompute some area-weighted interpolations of vertical scale factors22 !! : to account for manual changes to e[1,2][u,v] in some Straits23 20 !!---------------------------------------------------------------------- 24 !! * Modules used25 21 USE oce ! ocean dynamics and tracers 26 22 USE dom_oce ! ocean space and time domain … … 37 33 PRIVATE 38 34 39 !! * Routine accessibility40 35 PUBLIC dom_vvl_init ! called by domain.F90 41 36 PUBLIC dom_vvl_sf_nxt ! called by step.F90 42 37 PUBLIC dom_vvl_sf_swp ! called by step.F90 43 38 PUBLIC dom_vvl_interpol ! called by dynnxt.F90 44 PRIVATE dom_vvl_orca_fix ! called by dom_vvl_interpol 45 46 !!* Namelist nam_vvl 47 LOGICAL , PUBLIC :: ln_vvl_zstar = .FALSE. ! zstar vertical coordinate 48 LOGICAL , PUBLIC :: ln_vvl_ztilde = .FALSE. ! ztilde vertical coordinate 49 LOGICAL , PUBLIC :: ln_vvl_layer = .FALSE. ! level vertical coordinate 50 LOGICAL , PUBLIC :: ln_vvl_ztilde_as_zstar = .FALSE. ! ztilde vertical coordinate 51 LOGICAL , PUBLIC :: ln_vvl_zstar_at_eqtor = .FALSE. ! ztilde vertical coordinate 52 LOGICAL , PUBLIC :: ln_vvl_kepe = .FALSE. ! kinetic/potential energy transfer 53 ! ! conservation: not used yet 54 REAL(wp) :: rn_ahe3 ! thickness diffusion coefficient 55 REAL(wp) :: rn_rst_e3t ! ztilde to zstar restoration timescale [days] 56 REAL(wp) :: rn_lf_cutoff ! cutoff frequency for low-pass filter [days] 57 REAL(wp) :: rn_zdef_max ! maximum fractional e3t deformation 58 LOGICAL , PUBLIC :: ln_vvl_dbg = .FALSE. ! debug control prints 59 60 !! * Module variables 61 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: un_td, vn_td ! thickness diffusion transport 62 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hdiv_lf ! low frequency part of hz divergence 63 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_b, tilde_e3t_n ! baroclinic scale factors 64 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_a, dtilde_e3t_a ! baroclinic scale factors 65 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: frq_rst_e3t ! retoring period for scale factors 66 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: frq_rst_hdv ! retoring period for low freq. divergence 39 40 ! !!* Namelist nam_vvl 41 LOGICAL , PUBLIC :: ln_vvl_zstar = .FALSE. ! zstar vertical coordinate 42 LOGICAL , PUBLIC :: ln_vvl_ztilde = .FALSE. ! ztilde vertical coordinate 43 LOGICAL , PUBLIC :: ln_vvl_layer = .FALSE. ! level vertical coordinate 44 LOGICAL , PUBLIC :: ln_vvl_ztilde_as_zstar = .FALSE. ! ztilde vertical coordinate 45 LOGICAL , PUBLIC :: ln_vvl_zstar_at_eqtor = .FALSE. ! ztilde vertical coordinate 46 LOGICAL , PUBLIC :: ln_vvl_kepe = .FALSE. ! kinetic/potential energy transfer 47 ! ! conservation: not used yet 48 REAL(wp) :: rn_ahe3 ! thickness diffusion coefficient 49 REAL(wp) :: rn_rst_e3t ! ztilde to zstar restoration timescale [days] 50 REAL(wp) :: rn_lf_cutoff ! cutoff frequency for low-pass filter [days] 51 REAL(wp) :: rn_zdef_max ! maximum fractional e3t deformation 52 LOGICAL , PUBLIC :: ln_vvl_dbg = .FALSE. ! debug control prints 53 54 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: un_td, vn_td ! thickness diffusion transport 55 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hdiv_lf ! low frequency part of hz divergence 56 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_b, tilde_e3t_n ! baroclinic scale factors 57 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_a, dtilde_e3t_a ! baroclinic scale factors 58 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: frq_rst_e3t ! retoring period for scale factors 59 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: frq_rst_hdv ! retoring period for low freq. divergence 67 60 68 61 !! * Substitutions … … 372 365 DO jj = 1, jpjm1 373 366 DO ji = 1, fs_jpim1 ! vector opt. 374 un_td(ji,jj,jk) = rn_ahe3 * umask(ji,jj,jk) * re2u_e1u(ji,jj)&375 &* ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji+1,jj ,jk) )376 vn_td(ji,jj,jk) = rn_ahe3 * vmask(ji,jj,jk) * re1v_e2v(ji,jj)&377 &* ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji ,jj+1,jk) )367 un_td(ji,jj,jk) = rn_ahe3 * umask(ji,jj,jk) * e2_e1u(ji,jj) & 368 & * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji+1,jj ,jk) ) 369 vn_td(ji,jj,jk) = rn_ahe3 * vmask(ji,jj,jk) * e1_e2v(ji,jj) & 370 & * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji ,jj+1,jk) ) 378 371 zwu(ji,jj) = zwu(ji,jj) + un_td(ji,jj,jk) 379 372 zwv(ji,jj) = zwv(ji,jj) + vn_td(ji,jj,jk) … … 394 387 tilde_e3t_a(ji,jj,jk) = tilde_e3t_a(ji,jj,jk) + ( un_td(ji-1,jj ,jk) - un_td(ji,jj,jk) & 395 388 & + vn_td(ji ,jj-1,jk) - vn_td(ji,jj,jk) & 396 & ) * r1_e1 2t(ji,jj)389 & ) * r1_e1e2t(ji,jj) 397 390 END DO 398 391 END DO … … 695 688 !! - vertical interpolation: simple averaging 696 689 !!---------------------------------------------------------------------- 697 !! * Arguments 698 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: pe3_in ! input e3 to be interpolated 699 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( inout ) :: pe3_out ! output interpolated e3 700 CHARACTER(LEN=*), INTENT( in ) :: pout ! grid point of out scale factors 701 ! ! = 'U', 'V', 'W, 'F', 'UW' or 'VW' 702 !! * Local declarations 690 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pe3_in ! input e3 to be interpolated 691 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pe3_out ! output interpolated e3 692 CHARACTER(LEN=*) , INTENT(in ) :: pout ! grid point of out scale factors 693 ! ! = 'U', 'V', 'W, 'F', 'UW' or 'VW' 694 ! 703 695 INTEGER :: ji, jj, jk ! dummy loop indices 704 LOGICAL :: l_is_orca ! local logical705 ! !----------------------------------------------------------------------696 !!---------------------------------------------------------------------- 697 ! 706 698 IF( nn_timing == 1 ) CALL timing_start('dom_vvl_interpol') 707 ! 708 l_is_orca = .FALSE. 709 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) l_is_orca = .TRUE. ! ORCA R2 configuration - will need to correct some locations 710 711 SELECT CASE ( pout ) 712 ! ! ------------------------------------- ! 713 CASE( 'U' ) ! interpolation from T-point to U-point ! 714 ! ! ------------------------------------- ! 715 ! horizontal surface weighted interpolation 699 ! 700 SELECT CASE ( pout ) !== type of interpolation ==! 701 ! 702 CASE( 'U' ) !* from T- to U-point : hor. surface weighted mean 716 703 DO jk = 1, jpk 717 704 DO jj = 1, jpjm1 718 705 DO ji = 1, fs_jpim1 ! vector opt. 719 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * r1_e1 2u(ji,jj) &720 & * ( e1 2t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &721 & + e1 2t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) )706 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * r1_e1e2u(ji,jj) & 707 & * ( e1e2t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) & 708 & + e1e2t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) ) 722 709 END DO 723 710 END DO 724 711 END DO 725 !726 IF( l_is_orca ) CALL dom_vvl_orca_fix( pe3_in, pe3_out, pout )727 ! boundary conditions728 712 CALL lbc_lnk( pe3_out(:,:,:), 'U', 1._wp ) 729 713 pe3_out(:,:,:) = pe3_out(:,:,:) + e3u_0(:,:,:) 730 ! ! ------------------------------------- ! 731 CASE( 'V' ) ! interpolation from T-point to V-point ! 732 ! ! ------------------------------------- ! 733 ! horizontal surface weighted interpolation 714 ! 715 CASE( 'V' ) !* from T- to V-point : hor. surface weighted mean 734 716 DO jk = 1, jpk 735 717 DO jj = 1, jpjm1 736 718 DO ji = 1, fs_jpim1 ! vector opt. 737 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) * r1_e1 2v(ji,jj) &738 & * ( e1 2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) &739 & + e1 2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) )719 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) * r1_e1e2v(ji,jj) & 720 & * ( e1e2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) & 721 & + e1e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) ) 740 722 END DO 741 723 END DO 742 724 END DO 743 !744 IF( l_is_orca ) CALL dom_vvl_orca_fix( pe3_in, pe3_out, pout )745 ! boundary conditions746 725 CALL lbc_lnk( pe3_out(:,:,:), 'V', 1._wp ) 747 726 pe3_out(:,:,:) = pe3_out(:,:,:) + e3v_0(:,:,:) 748 ! ! ------------------------------------- ! 749 CASE( 'F' ) ! interpolation from U-point to F-point ! 750 ! ! ------------------------------------- ! 751 ! horizontal surface weighted interpolation 727 ! 728 CASE( 'F' ) !* from U-point to F-point : hor. surface weighted mean 752 729 DO jk = 1, jpk 753 730 DO jj = 1, jpjm1 754 731 DO ji = 1, fs_jpim1 ! vector opt. 755 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) * r1_e1 2f(ji,jj) &756 & * ( e1 2u(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3u_0(ji,jj ,jk) ) &757 & + e1 2u(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3u_0(ji,jj+1,jk) ) )732 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) * r1_e1e2f(ji,jj) & 733 & * ( e1e2u(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3u_0(ji,jj ,jk) ) & 734 & + e1e2u(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3u_0(ji,jj+1,jk) ) ) 758 735 END DO 759 736 END DO 760 737 END DO 761 !762 IF( l_is_orca ) CALL dom_vvl_orca_fix( pe3_in, pe3_out, pout )763 ! boundary conditions764 738 CALL lbc_lnk( pe3_out(:,:,:), 'F', 1._wp ) 765 739 pe3_out(:,:,:) = pe3_out(:,:,:) + e3f_0(:,:,:) 766 ! ! ------------------------------------- ! 767 CASE( 'W' ) ! interpolation from T-point to W-point ! 768 ! ! ------------------------------------- ! 769 ! vertical simple interpolation 740 ! 741 CASE( 'W' ) !* from T- to W-point : vertical simple mean 742 ! 770 743 pe3_out(:,:,1) = e3w_0(:,:,1) + pe3_in(:,:,1) - e3t_0(:,:,1) 771 ! - ML - The use of mask in this formaula enables the special treatment of the last w- point without indirect adressing 744 ! - ML - The use of mask in this formulea enables the special treatment of the last w-point without indirect adressing 745 !!gm BUG? use here wmask in case of ISF ? to be checked 772 746 DO jk = 2, jpk 773 747 pe3_out(:,:,jk) = e3w_0(:,:,jk) + ( 1.0_wp - 0.5_wp * tmask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3t_0(:,:,jk-1) ) & 774 748 & + 0.5_wp * tmask(:,:,jk) * ( pe3_in(:,:,jk ) - e3t_0(:,:,jk ) ) 775 749 END DO 776 ! ! -------------------------------------- ! 777 CASE( 'UW' ) ! interpolation from U-point to UW-point ! 778 ! ! -------------------------------------- ! 779 ! vertical simple interpolation 750 ! 751 CASE( 'UW' ) !* from U- to UW-point : vertical simple mean 752 ! 780 753 pe3_out(:,:,1) = e3uw_0(:,:,1) + pe3_in(:,:,1) - e3u_0(:,:,1) 781 754 ! - ML - The use of mask in this formaula enables the special treatment of the last w- point without indirect adressing 755 !!gm BUG? use here wumask in case of ISF ? to be checked 782 756 DO jk = 2, jpk 783 757 pe3_out(:,:,jk) = e3uw_0(:,:,jk) + ( 1.0_wp - 0.5_wp * umask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3u_0(:,:,jk-1) ) & 784 758 & + 0.5_wp * umask(:,:,jk) * ( pe3_in(:,:,jk ) - e3u_0(:,:,jk ) ) 785 759 END DO 786 ! ! -------------------------------------- ! 787 CASE( 'VW' ) ! interpolation from V-point to VW-point ! 788 ! ! -------------------------------------- ! 789 ! vertical simple interpolation 760 ! 761 CASE( 'VW' ) !* from V- to VW-point : vertical simple mean 762 ! 790 763 pe3_out(:,:,1) = e3vw_0(:,:,1) + pe3_in(:,:,1) - e3v_0(:,:,1) 791 764 ! - ML - The use of mask in this formaula enables the special treatment of the last w- point without indirect adressing 765 !!gm BUG? use here wvmask in case of ISF ? to be checked 792 766 DO jk = 2, jpk 793 767 pe3_out(:,:,jk) = e3vw_0(:,:,jk) + ( 1.0_wp - 0.5_wp * vmask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3v_0(:,:,jk-1) ) & … … 796 770 END SELECT 797 771 ! 798 799 772 IF( nn_timing == 1 ) CALL timing_stop('dom_vvl_interpol') 800 773 ! 801 774 END SUBROUTINE dom_vvl_interpol 775 802 776 803 777 SUBROUTINE dom_vvl_rst( kt, cdrw ) … … 813 787 !! they are set to 0. 814 788 !!---------------------------------------------------------------------- 815 !! * Arguments816 789 INTEGER , INTENT(in) :: kt ! ocean time-step 817 790 CHARACTER(len=*), INTENT(in) :: cdrw ! "READ"/"WRITE" flag 818 ! ! * Local declarations791 ! 819 792 INTEGER :: jk 820 793 INTEGER :: id1, id2, id3, id4, id5 ! local integers … … 911 884 END IF 912 885 ENDIF 913 886 ! 914 887 ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN ! Create restart file 915 888 ! ! =================== … … 931 904 CALL iom_rstput( kt, nitrst, numrow, 'hdiv_lf', hdiv_lf(:,:,:) ) 932 905 ENDIF 933 934 ENDIF 906 ! 907 ENDIF 908 ! 935 909 IF( nn_timing == 1 ) CALL timing_stop('dom_vvl_rst') 936 910 ! 937 911 END SUBROUTINE dom_vvl_rst 938 912 … … 945 919 !! for vertical coordinate 946 920 !!---------------------------------------------------------------------- 947 INTEGER :: ioptio 948 INTEGER :: ios 949 921 INTEGER :: ioptio, ios 922 !! 950 923 NAMELIST/nam_vvl/ ln_vvl_zstar, ln_vvl_ztilde, ln_vvl_layer, ln_vvl_ztilde_as_zstar, & 951 &ln_vvl_zstar_at_eqtor , rn_ahe3 , rn_rst_e3t , &952 &rn_lf_cutoff , rn_zdef_max , ln_vvl_dbg ! not yet implemented: ln_vvl_kepe924 & ln_vvl_zstar_at_eqtor , rn_ahe3 , rn_rst_e3t , & 925 & rn_lf_cutoff , rn_zdef_max , ln_vvl_dbg ! not yet implemented: ln_vvl_kepe 953 926 !!---------------------------------------------------------------------- 954 927 ! 955 928 REWIND( numnam_ref ) ! Namelist nam_vvl in reference namelist : 956 929 READ ( numnam_ref, nam_vvl, IOSTAT = ios, ERR = 901) 957 930 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_vvl in reference namelist', lwp ) 958 931 ! 959 932 REWIND( numnam_cfg ) ! Namelist nam_vvl in configuration namelist : Parameters of the run 960 933 READ ( numnam_cfg, nam_vvl, IOSTAT = ios, ERR = 902 ) 961 934 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_vvl in configuration namelist', lwp ) 962 935 IF(lwm) WRITE ( numond, nam_vvl ) 963 936 ! 964 937 IF(lwp) THEN ! Namelist print 965 938 WRITE(numout,*) … … 994 967 WRITE(numout,*) ' ln_vvl_dbg = ', ln_vvl_dbg 995 968 ENDIF 996 969 ! 997 970 ioptio = 0 ! Parameter control 998 IF( ln_vvl_ztilde_as_zstar ) ln_vvl_ztilde = .true.999 IF( ln_vvl_zstar ) 1000 IF( ln_vvl_ztilde ) 1001 IF( ln_vvl_layer ) 1002 971 IF( ln_vvl_ztilde_as_zstar ) ln_vvl_ztilde = .true. 972 IF( ln_vvl_zstar ) ioptio = ioptio + 1 973 IF( ln_vvl_ztilde ) ioptio = ioptio + 1 974 IF( ln_vvl_layer ) ioptio = ioptio + 1 975 ! 1003 976 IF( ioptio /= 1 ) CALL ctl_stop( 'Choose ONE vertical coordinate in namelist nam_vvl' ) 1004 977 IF( .NOT. ln_vvl_zstar .AND. nn_isf .NE. 0) CALL ctl_stop( 'Only vvl_zstar has been tested with ice shelf cavity' ) 1005 978 ! 1006 979 IF(lwp) THEN ! Print the choice 1007 980 WRITE(numout,*) … … 1014 987 ! IF( .NOT. ln_vvl_kepe ) WRITE(numout,*) ' kinetic to potential energy transfer : option not used' 1015 988 ENDIF 1016 989 ! 1017 990 #if defined key_agrif 1018 991 IF (.NOT.Agrif_Root()) CALL ctl_stop( 'AGRIF not implemented with non-linear free surface (key_vvl)' ) 1019 992 #endif 1020 993 ! 1021 994 END SUBROUTINE dom_vvl_ctl 1022 1023 SUBROUTINE dom_vvl_orca_fix( pe3_in, pe3_out, pout )1024 !!---------------------------------------------------------------------1025 !! *** ROUTINE dom_vvl_orca_fix ***1026 !!1027 !! ** Purpose : Correct surface weighted, horizontally interpolated,1028 !! scale factors at locations that have been individually1029 !! modified in domhgr. Such modifications break the1030 !! relationship between e12t and e1u*e2u etc.1031 !! Recompute some scale factors ignoring the modified metric.1032 !!----------------------------------------------------------------------1033 !! * Arguments1034 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: pe3_in ! input e3 to be interpolated1035 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( inout ) :: pe3_out ! output interpolated e31036 CHARACTER(LEN=*), INTENT( in ) :: pout ! grid point of out scale factors1037 ! ! = 'U', 'V', 'W, 'F', 'UW' or 'VW'1038 !! * Local declarations1039 INTEGER :: ji, jj, jk ! dummy loop indices1040 INTEGER :: ij0, ij1, ii0, ii1 ! dummy loop indices1041 INTEGER :: isrow ! index for ORCA1 starting row1042 !! acc1043 !! Hmm with the time splitting these "fixes" seem to do more harm than good. Temporarily disabled for1044 !! the ORCA2 tests (by changing jp_cfg test from 2 to 3) pending further investigations1045 !!1046 ! ! =====================1047 IF( cp_cfg == "orca" .AND. jp_cfg == 3 ) THEN ! ORCA R2 configuration1048 ! ! =====================1049 !! acc1050 IF( nn_cla == 0 ) THEN1051 !1052 ii0 = 139 ; ii1 = 140 ! Gibraltar Strait (e2u was modified)1053 ij0 = 102 ; ij1 = 1021054 DO jk = 1, jpkm11055 DO jj = mj0(ij0), mj1(ij1)1056 DO ji = mi0(ii0), mi1(ii1)1057 SELECT CASE ( pout )1058 CASE( 'U' )1059 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) &1060 & * ( e1t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &1061 & + e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) &1062 & ) / e1u(ji,jj) + e3u_0(ji,jj,jk)1063 CASE( 'F' )1064 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) &1065 & * ( e1u(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3u_0(ji ,jj,jk) ) &1066 & + e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) &1067 & ) / e1f(ji,jj) + e3f_0(ji,jj,jk)1068 END SELECT1069 END DO1070 END DO1071 END DO1072 !1073 ii0 = 160 ; ii1 = 160 ! Bab el Mandeb (e2u and e1v were modified)1074 ij0 = 88 ; ij1 = 881075 DO jk = 1, jpkm11076 DO jj = mj0(ij0), mj1(ij1)1077 DO ji = mi0(ii0), mi1(ii1)1078 SELECT CASE ( pout )1079 CASE( 'U' )1080 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) &1081 & * ( e1t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &1082 & + e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) &1083 & ) / e1u(ji,jj) + e3u_0(ji,jj,jk)1084 CASE( 'V' )1085 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) &1086 & * ( e2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) &1087 & + e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) &1088 & ) / e2v(ji,jj) + e3v_0(ji,jj,jk)1089 CASE( 'F' )1090 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) &1091 & * ( e1u(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3u_0(ji ,jj,jk) ) &1092 & + e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) &1093 & ) / e1f(ji,jj) + e3f_0(ji,jj,jk)1094 END SELECT1095 END DO1096 END DO1097 END DO1098 ENDIF1099 1100 ii0 = 145 ; ii1 = 146 ! Danish Straits (e2u was modified)1101 ij0 = 116 ; ij1 = 1161102 DO jk = 1, jpkm11103 DO jj = mj0(ij0), mj1(ij1)1104 DO ji = mi0(ii0), mi1(ii1)1105 SELECT CASE ( pout )1106 CASE( 'U' )1107 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) &1108 & * ( e1t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &1109 & + e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) &1110 & ) / e1u(ji,jj) + e3u_0(ji,jj,jk)1111 CASE( 'F' )1112 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) &1113 & * ( e1u(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3u_0(ji ,jj,jk) ) &1114 & + e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) &1115 & ) / e1f(ji,jj) + e3f_0(ji,jj,jk)1116 END SELECT1117 END DO1118 END DO1119 END DO1120 ENDIF1121 !1122 ! ! =====================1123 IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN ! ORCA R1 configuration1124 ! ! =====================1125 ! This dirty section will be suppressed by simplification process:1126 ! all this will come back in input files1127 ! Currently these hard-wired indices relate to configuration with1128 ! extend grid (jpjglo=332)1129 ! which had a grid-size of 362x292.1130 isrow = 332 - jpjglo1131 !1132 ii0 = 282 ; ii1 = 283 ! Gibraltar Strait (e2u was modified)1133 ij0 = 241 - isrow ; ij1 = 241 - isrow1134 DO jk = 1, jpkm11135 DO jj = mj0(ij0), mj1(ij1)1136 DO ji = mi0(ii0), mi1(ii1)1137 SELECT CASE ( pout )1138 CASE( 'U' )1139 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) &1140 & * ( e1t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &1141 & + e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) &1142 & ) / e1u(ji,jj) + e3u_0(ji,jj,jk)1143 CASE( 'F' )1144 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) &1145 & * ( e1u(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3u_0(ji ,jj,jk) ) &1146 & + e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) &1147 & ) / e1f(ji,jj) + e3f_0(ji,jj,jk)1148 END SELECT1149 END DO1150 END DO1151 END DO1152 !1153 ii0 = 314 ; ii1 = 315 ! Bhosporus Strait (e2u was modified)1154 ij0 = 248 - isrow ; ij1 = 248 - isrow1155 DO jk = 1, jpkm11156 DO jj = mj0(ij0), mj1(ij1)1157 DO ji = mi0(ii0), mi1(ii1)1158 SELECT CASE ( pout )1159 CASE( 'U' )1160 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) &1161 & * ( e1t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &1162 & + e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) &1163 & ) / e1u(ji,jj) + e3u_0(ji,jj,jk)1164 CASE( 'F' )1165 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) &1166 & * ( e1u(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3u_0(ji ,jj,jk) ) &1167 & + e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) &1168 & ) / e1f(ji,jj) + e3f_0(ji,jj,jk)1169 END SELECT1170 END DO1171 END DO1172 END DO1173 !1174 ii0 = 44 ; ii1 = 44 ! Lombok Strait (e1v was modified)1175 ij0 = 164 - isrow ; ij1 = 165 - isrow1176 DO jk = 1, jpkm11177 DO jj = mj0(ij0), mj1(ij1)1178 DO ji = mi0(ii0), mi1(ii1)1179 SELECT CASE ( pout )1180 CASE( 'V' )1181 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) &1182 & * ( e2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) &1183 & + e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) &1184 & ) / e2v(ji,jj) + e3v_0(ji,jj,jk)1185 END SELECT1186 END DO1187 END DO1188 END DO1189 !1190 ii0 = 48 ; ii1 = 48 ! Sumba Strait (e1v was modified) [closed from bathy_11 on]1191 ij0 = 164 - isrow ; ij1 = 165 - isrow1192 DO jk = 1, jpkm11193 DO jj = mj0(ij0), mj1(ij1)1194 DO ji = mi0(ii0), mi1(ii1)1195 SELECT CASE ( pout )1196 CASE( 'V' )1197 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) &1198 & * ( e2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) &1199 & + e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) &1200 & ) / e2v(ji,jj) + e3v_0(ji,jj,jk)1201 END SELECT1202 END DO1203 END DO1204 END DO1205 !1206 ii0 = 53 ; ii1 = 53 ! Ombai Strait (e1v was modified)1207 ij0 = 164 - isrow ; ij1 = 165 - isrow1208 DO jk = 1, jpkm11209 DO jj = mj0(ij0), mj1(ij1)1210 DO ji = mi0(ii0), mi1(ii1)1211 SELECT CASE ( pout )1212 CASE( 'V' )1213 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) &1214 & * ( e2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) &1215 & + e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) &1216 & ) / e2v(ji,jj) + e3v_0(ji,jj,jk)1217 END SELECT1218 END DO1219 END DO1220 END DO1221 !1222 ii0 = 56 ; ii1 = 56 ! Timor Passage (e1v was modified)1223 ij0 = 164 - isrow ; ij1 = 165 - isrow1224 DO jk = 1, jpkm11225 DO jj = mj0(ij0), mj1(ij1)1226 DO ji = mi0(ii0), mi1(ii1)1227 SELECT CASE ( pout )1228 CASE( 'V' )1229 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) &1230 & * ( e2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) &1231 & + e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) &1232 & ) / e2v(ji,jj) + e3v_0(ji,jj,jk)1233 END SELECT1234 END DO1235 END DO1236 END DO1237 !1238 ii0 = 55 ; ii1 = 55 ! West Halmahera Strait (e1v was modified)1239 ij0 = 181 - isrow ; ij1 = 182 - isrow1240 DO jk = 1, jpkm11241 DO jj = mj0(ij0), mj1(ij1)1242 DO ji = mi0(ii0), mi1(ii1)1243 SELECT CASE ( pout )1244 CASE( 'V' )1245 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) &1246 & * ( e2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) &1247 & + e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) &1248 & ) / e2v(ji,jj) + e3v_0(ji,jj,jk)1249 END SELECT1250 END DO1251 END DO1252 END DO1253 !1254 ii0 = 58 ; ii1 = 58 ! East Halmahera Strait (e1v was modified)1255 ij0 = 181 - isrow ; ij1 = 182 - isrow1256 DO jk = 1, jpkm11257 DO jj = mj0(ij0), mj1(ij1)1258 DO ji = mi0(ii0), mi1(ii1)1259 SELECT CASE ( pout )1260 CASE( 'V' )1261 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) &1262 & * ( e2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) &1263 & + e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) &1264 & ) / e2v(ji,jj) + e3v_0(ji,jj,jk)1265 END SELECT1266 END DO1267 END DO1268 END DO1269 ENDIF1270 ! ! =====================1271 IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN ! ORCA R05 configuration1272 ! ! =====================1273 !1274 ii0 = 563 ; ii1 = 564 ! Gibraltar Strait (e2u was modified)1275 ij0 = 327 ; ij1 = 3271276 DO jk = 1, jpkm11277 DO jj = mj0(ij0), mj1(ij1)1278 DO ji = mi0(ii0), mi1(ii1)1279 SELECT CASE ( pout )1280 CASE( 'U' )1281 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) &1282 & * ( e1t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &1283 & + e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) &1284 & ) / e1u(ji,jj) + e3u_0(ji,jj,jk)1285 CASE( 'F' )1286 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) &1287 & * ( e1u(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3u_0(ji ,jj,jk) ) &1288 & + e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) &1289 & ) / e1f(ji,jj) + e3f_0(ji,jj,jk)1290 END SELECT1291 END DO1292 END DO1293 END DO1294 !1295 ii0 = 627 ; ii1 = 628 ! Bosphorus Strait (e2u was modified)1296 ij0 = 343 ; ij1 = 3431297 DO jk = 1, jpkm11298 DO jj = mj0(ij0), mj1(ij1)1299 DO ji = mi0(ii0), mi1(ii1)1300 SELECT CASE ( pout )1301 CASE( 'U' )1302 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) &1303 & * ( e1t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &1304 & + e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) &1305 & ) / e1u(ji,jj) + e3u_0(ji,jj,jk)1306 CASE( 'F' )1307 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) &1308 & * ( e1u(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3u_0(ji ,jj,jk) ) &1309 & + e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) &1310 & ) / e1f(ji,jj) + e3f_0(ji,jj,jk)1311 END SELECT1312 END DO1313 END DO1314 END DO1315 !1316 ii0 = 93 ; ii1 = 94 ! Sumba Strait (e2u was modified)1317 ij0 = 232 ; ij1 = 2321318 DO jk = 1, jpkm11319 DO jj = mj0(ij0), mj1(ij1)1320 DO ji = mi0(ii0), mi1(ii1)1321 SELECT CASE ( pout )1322 CASE( 'U' )1323 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) &1324 & * ( e1t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &1325 & + e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) &1326 & ) / e1u(ji,jj) + e3u_0(ji,jj,jk)1327 CASE( 'F' )1328 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) &1329 & * ( e1u(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3u_0(ji ,jj,jk) ) &1330 & + e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) &1331 & ) / e1f(ji,jj) + e3f_0(ji,jj,jk)1332 END SELECT1333 END DO1334 END DO1335 END DO1336 !1337 ii0 = 103 ; ii1 = 103 ! Ombai Strait (e2u was modified)1338 ij0 = 232 ; ij1 = 2321339 DO jk = 1, jpkm11340 DO jj = mj0(ij0), mj1(ij1)1341 DO ji = mi0(ii0), mi1(ii1)1342 SELECT CASE ( pout )1343 CASE( 'U' )1344 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) &1345 & * ( e1t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &1346 & + e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) &1347 & ) / e1u(ji,jj) + e3u_0(ji,jj,jk)1348 CASE( 'F' )1349 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) &1350 & * ( e1u(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3u_0(ji ,jj,jk) ) &1351 & + e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) &1352 & ) / e1f(ji,jj) + e3f_0(ji,jj,jk)1353 END SELECT1354 END DO1355 END DO1356 END DO1357 !1358 ii0 = 15 ; ii1 = 15 ! Palk Strait (e2u was modified)1359 ij0 = 270 ; ij1 = 2701360 DO jk = 1, jpkm11361 DO jj = mj0(ij0), mj1(ij1)1362 DO ji = mi0(ii0), mi1(ii1)1363 SELECT CASE ( pout )1364 CASE( 'U' )1365 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) &1366 & * ( e1t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) &1367 & + e1t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) &1368 & ) / e1u(ji,jj) + e3u_0(ji,jj,jk)1369 CASE( 'F' )1370 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) &1371 & * ( e1u(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3u_0(ji ,jj,jk) ) &1372 & + e1u(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3u_0(ji+1,jj,jk) ) &1373 & ) / e1f(ji,jj) + e3f_0(ji,jj,jk)1374 END SELECT1375 END DO1376 END DO1377 END DO1378 !1379 ii0 = 87 ; ii1 = 87 ! Lombok Strait (e1v was modified)1380 ij0 = 232 ; ij1 = 2331381 DO jk = 1, jpkm11382 DO jj = mj0(ij0), mj1(ij1)1383 DO ji = mi0(ii0), mi1(ii1)1384 SELECT CASE ( pout )1385 CASE( 'V' )1386 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) &1387 & * ( e2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) &1388 & + e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) &1389 & ) / e2v(ji,jj) + e3v_0(ji,jj,jk)1390 END SELECT1391 END DO1392 END DO1393 END DO1394 !1395 ii0 = 662 ; ii1 = 662 ! Bab el Mandeb (e1v was modified)1396 ij0 = 276 ; ij1 = 2761397 DO jk = 1, jpkm11398 DO jj = mj0(ij0), mj1(ij1)1399 DO ji = mi0(ii0), mi1(ii1)1400 SELECT CASE ( pout )1401 CASE( 'V' )1402 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) &1403 & * ( e2t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) &1404 & + e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) &1405 & ) / e2v(ji,jj) + e3v_0(ji,jj,jk)1406 END SELECT1407 END DO1408 END DO1409 END DO1410 ENDIF1411 END SUBROUTINE dom_vvl_orca_fix1412 995 1413 996 !!======================================================================
Note: See TracChangeset
for help on using the changeset viewer.