Changeset 5842 for branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/OPA_SRC/DOM
- Timestamp:
- 2015-10-30T16:34:24+01:00 (9 years ago)
- Location:
- branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/OPA_SRC/DOM
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90
r5123 r5842 10 10 !! 3.5 ! 2012 (S. Mocavero, I. Epicoco) Add arrays associated 11 11 !! to the optimization of BDY communications 12 !! 3.6.?! 2014 (H. Liu) Add arrays associated with Wetting and Drying 12 13 !!---------------------------------------------------------------------- 13 14 … … 270 271 INTEGER, PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: nicoa, njcoa !: ??? 271 272 #endif 273 274 !!---------------------------------------------------------------------- 275 !! critical depths,filters, limiters,and masks for Wetting and Drying 276 !! --------------------------------------------------------------------- 277 278 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wduflt, wdvflt !: u- and v- filter 279 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wdmask !: u- and v- limiter 280 281 LOGICAL, PUBLIC, SAVE :: ln_wd !: key to turn on/off wetting/drying (T: on, F: off) 282 REAL(wp), PUBLIC, SAVE :: rn_wdmin1 !: minimum water depth on dried cells 283 REAL(wp), PUBLIC, SAVE :: rn_wdmin2 !: tolerrance of minimum water depth on dried cells 284 REAL(wp), PUBLIC, SAVE :: rn_wdld !: land elevation below which wetting/drying will be considered 285 INTEGER , PUBLIC, SAVE :: nn_wdit !: maximum number of iteration for W/D limiter 286 287 !LOGICAL, PUBLIC, SAVE :: ln_wd = .FALSE. !: turn on wetting/drying (T: on, F: off) 288 !REAL(wp), PUBLIC, SAVE :: rn_wdmin1 = 0.10_wp !: minimum water depth on dried cells 289 !REAL(wp), PUBLIC, SAVE :: rn_wdmin2 = 0.01_wp !: tolerrance of minimum water depth on dried cells 290 !REAL(wp), PUBLIC, SAVE :: rn_wdld = 20.0_wp !: land elevation below which wetting/drying will be considered 291 !INTEGER , PUBLIC, SAVE :: nn_wdit = 10 !: maximum number of iteration for W/D limiter 272 292 273 293 !!---------------------------------------------------------------------- … … 407 427 ALLOCATE( npcoa(4,jpk), nicoa(2*(jpi+jpj),4,jpk), njcoa(2*(jpi+jpj),4,jpk), STAT=ierr(12) ) 408 428 #endif 429 430 IF(ln_wd) & 431 ALLOCATE( wduflt(jpi,jpj), wdvflt(jpi,jpj), wdmask(jpi,jpj), STAT=ierr(12) ) 409 432 ! 410 433 dom_oce_alloc = MAXVAL(ierr) -
branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
r5506 r5842 147 147 fse3t_a(:,:,jpk) = e3t_0(:,:,jpk) 148 148 149 !IF(ln_wd) THEN 150 ! DO jj = 1, jpj 151 ! DO ji = 1, jpi 152 ! IF(mbathy(ji,jj) == 2 .AND. e3t_0(ji,jj,1) <= 0.5_wp * rn_wdmin1) THEN 153 ! fse3t_a(ji,jj,1:2) = 0.5_wp * rn_wdmin1 154 ! fse3t_n(ji,jj,1:2) = 0.5_wp * rn_wdmin1 155 ! fse3t_b(ji,jj,1:2) = 0.5_wp * rn_wdmin1 156 ! END IF 157 ! ENDDO 158 ! ENDDO 159 !END IF 160 149 161 ! Reconstruction of all vertical scale factors at now and before time steps 150 162 ! ============================================================================= … … 701 713 ! ! = 'U', 'V', 'W, 'F', 'UW' or 'VW' 702 714 !! * Local declarations 715 REAL(wp) :: zwad ! = 1.0 when ln_wd = .true. 716 ! = 0.0 when ln_wd = .false. 717 ! 703 718 INTEGER :: ji, jj, jk ! dummy loop indices 704 719 LOGICAL :: l_is_orca ! local logical … … 706 721 IF( nn_timing == 1 ) CALL timing_start('dom_vvl_interpol') 707 722 ! 723 IF(ln_wd) THEN 724 zwad = 1.0_wp 725 ELSE 726 zwad = 0.0_wp 727 END IF 728 708 729 l_is_orca = .FALSE. 709 730 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) l_is_orca = .TRUE. ! ORCA R2 configuration - will need to correct some locations … … 717 738 DO jj = 1, jpjm1 718 739 DO ji = 1, fs_jpim1 ! vector opt. 719 pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * r1_e12u(ji,jj) & 740 !pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * r1_e12u(ji,jj) & 741 pe3_out(ji,jj,jk) = 0.5_wp * (umask(ji,jj,jk) * (1.0_wp - zwad) + zwad) * r1_e12u(ji,jj) & 720 742 & * ( e12t(ji ,jj) * ( pe3_in(ji ,jj,jk) - e3t_0(ji ,jj,jk) ) & 721 743 & + e12t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) ) … … 735 757 DO jj = 1, jpjm1 736 758 DO ji = 1, fs_jpim1 ! vector opt. 737 pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) * r1_e12v(ji,jj) & 759 !pe3_out(ji,jj,jk) = 0.5_wp * vmask(ji,jj,jk) * r1_e12v(ji,jj) & 760 pe3_out(ji,jj,jk) = 0.5_wp * (vmask(ji,jj,jk) * (1.0_wp - zwad) + zwad) * r1_e12v(ji,jj) & 738 761 & * ( e12t(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3t_0(ji,jj ,jk) ) & 739 762 & + e12t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) ) … … 753 776 DO jj = 1, jpjm1 754 777 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_e12f(ji,jj) & 778 !pe3_out(ji,jj,jk) = 0.5_wp * umask(ji,jj,jk) * umask(ji,jj+1,jk) * r1_e12f(ji,jj) & 779 pe3_out(ji,jj,jk) = 0.5_wp * (umask(ji,jj,jk) * umask(ji,jj+1,jk) * (1.0_wp - zwad) + zwad) & 780 & * r1_e12f(ji,jj) & 756 781 & * ( e12u(ji,jj ) * ( pe3_in(ji,jj ,jk) - e3u_0(ji,jj ,jk) ) & 757 782 & + e12u(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3u_0(ji,jj+1,jk) ) ) … … 771 796 ! - ML - The use of mask in this formaula enables the special treatment of the last w- point without indirect adressing 772 797 DO jk = 2, jpk 773 pe3_out(:,:,jk) = e3w_0(:,:,jk) + ( 1.0_wp - 0.5_wp * tmask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3t_0(:,:,jk-1) ) & 774 & + 0.5_wp * tmask(:,:,jk) * ( pe3_in(:,:,jk ) - e3t_0(:,:,jk ) ) 798 !pe3_out(:,:,jk) = e3w_0(:,:,jk) + ( 1.0_wp - 0.5_wp * tmask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3t_0(:,:,jk-1) ) & 799 ! & + 0.5_wp * tmask(:,:,jk) * ( pe3_in(:,:,jk ) - e3t_0(:,:,jk ) ) 800 pe3_out(:,:,jk) = e3w_0(:,:,jk) + ( 1.0_wp - 0.5_wp * (tmask(:,:,jk) * (1.0_wp - zwad) + zwad) ) & 801 & * ( pe3_in(:,:,jk-1) - e3t_0(:,:,jk-1) ) & 802 & + 0.5_wp * (tmask(:,:,jk) * (1.0_wp - zwad) + zwad) & 803 & * ( pe3_in(:,:,jk ) - e3t_0(:,:,jk ) ) 775 804 END DO 776 805 ! ! -------------------------------------- ! … … 781 810 ! - ML - The use of mask in this formaula enables the special treatment of the last w- point without indirect adressing 782 811 DO jk = 2, jpk 783 pe3_out(:,:,jk) = e3uw_0(:,:,jk) + ( 1.0_wp - 0.5_wp * umask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3u_0(:,:,jk-1) ) & 784 & + 0.5_wp * umask(:,:,jk) * ( pe3_in(:,:,jk ) - e3u_0(:,:,jk ) ) 812 !pe3_out(:,:,jk) = e3uw_0(:,:,jk) + ( 1.0_wp - 0.5_wp * umask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3u_0(:,:,jk-1) ) & 813 ! & + 0.5_wp * umask(:,:,jk) * ( pe3_in(:,:,jk ) - e3u_0(:,:,jk ) ) 814 pe3_out(:,:,jk) = e3uw_0(:,:,jk) + ( 1.0_wp - 0.5_wp * (umask(:,:,jk) * (1.0_wp - zwad) + zwad) ) & 815 & * ( pe3_in(:,:,jk-1) - e3u_0(:,:,jk-1) ) & 816 & + 0.5_wp * (umask(:,:,jk) * (1.0_wp - zwad) + zwad) & 817 & * ( pe3_in(:,:,jk ) - e3u_0(:,:,jk ) ) 785 818 END DO 786 819 ! ! -------------------------------------- ! … … 791 824 ! - ML - The use of mask in this formaula enables the special treatment of the last w- point without indirect adressing 792 825 DO jk = 2, jpk 793 pe3_out(:,:,jk) = e3vw_0(:,:,jk) + ( 1.0_wp - 0.5_wp * vmask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3v_0(:,:,jk-1) ) & 794 & + 0.5_wp * vmask(:,:,jk) * ( pe3_in(:,:,jk ) - e3v_0(:,:,jk ) ) 826 !pe3_out(:,:,jk) = e3vw_0(:,:,jk) + ( 1.0_wp - 0.5_wp * vmask(:,:,jk) ) * ( pe3_in(:,:,jk-1) - e3v_0(:,:,jk-1) ) & 827 ! & + 0.5_wp * vmask(:,:,jk) * ( pe3_in(:,:,jk ) - e3v_0(:,:,jk ) ) 828 pe3_out(:,:,jk) = e3vw_0(:,:,jk) + ( 1.0_wp - 0.5_wp * (vmask(:,:,jk) * (1.0_wp - zwad) + zwad) ) & 829 & * ( pe3_in(:,:,jk-1) - e3v_0(:,:,jk-1) ) & 830 & + 0.5_wp * (vmask(:,:,jk) * (1.0_wp - zwad) + zwad) & 831 & * ( pe3_in(:,:,jk ) - e3v_0(:,:,jk ) ) 795 832 END DO 796 833 END SELECT … … 817 854 CHARACTER(len=*), INTENT(in) :: cdrw ! "READ"/"WRITE" flag 818 855 !! * Local declarations 819 INTEGER :: j k856 INTEGER :: ji, jj, jk 820 857 INTEGER :: id1, id2, id3, id4, id5 ! local integers 821 858 !!---------------------------------------------------------------------- … … 905 942 fse3t_n(:,:,:) = e3t_0(:,:,:) 906 943 sshn(:,:) = 0.0_wp 944 945 IF(ln_wd) THEN 946 DO jj = 1, jpj 947 DO ji = 1, jpi 948 !IF(e3t_0(ji,jj,1) < 0._wp) THEN 949 !IF(mbathy(ji,jj) == 2 .AND. e3t_0(ji,jj,1) <= 0.5_wp * rn_wdmin1) THEN 950 IF( e3t_0(ji,jj,1) <= 0.5_wp * rn_wdmin1) THEN 951 fse3t_b(ji,jj,:) = 0.5_wp * rn_wdmin1 952 fse3t_n(ji,jj,:) = 0.5_wp * rn_wdmin1 953 fse3t_a(ji,jj,:) = 0.5_wp * rn_wdmin1 954 sshb(ji,jj) = rn_wdmin1 - bathy(ji,jj) 955 sshn(ji,jj) = rn_wdmin1 - bathy(ji,jj) 956 ssha(ji,jj) = rn_wdmin1 - bathy(ji,jj) 957 ENDIF 958 ENDDO 959 ENDDO 960 END IF 961 907 962 IF( ln_vvl_ztilde .OR. ln_vvl_layer) THEN 908 963 tilde_e3t_b(:,:,:) = 0.0_wp -
branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90
r5656 r5842 18 18 !! 3.4 ! 2012-12 (R. Bourdalle-Badie and G. Reffray) modify C1D case 19 19 !! 3.6 ! 2014-11 (P. Mathiot and C. Harris) add ice shelf capabilitye 20 !! 3.? ! 2015-11 (H. Liu) Modifications for Wetting/Drying 20 21 !!---------------------------------------------------------------------- 21 22 … … 1810 1811 REAL(wp) :: zrmax, ztaper ! temporary scalars 1811 1812 REAL(wp) :: zrfact 1813 REAL(wp) :: zsmth 1812 1814 ! 1813 1815 REAL(wp), POINTER, DIMENSION(:,: ) :: ztmpi1, ztmpi2, ztmpj1, ztmpj2 … … 1865 1867 bathy(:,:) = MIN( rn_sbot_max, bathy(:,:) ) 1866 1868 1869 IF(.NOT.ln_wd) THEN 1867 1870 DO jj = 1, jpj 1868 1871 DO ji = 1, jpi … … 1870 1873 END DO 1871 1874 END DO 1875 END IF 1872 1876 ! ! ============================= 1873 1877 ! ! Define the envelop bathymetry (hbatt) … … 1876 1880 zenv(:,:) = bathy(:,:) 1877 1881 ! 1882 IF(.NOT.ln_wd) THEN 1878 1883 ! set first land point adjacent to a wet cell to sbot_min as this needs to be included in smoothing 1879 1884 DO jj = 1, jpj … … 1891 1896 END DO 1892 1897 END DO 1898 END IF 1899 1893 1900 ! apply lateral boundary condition CAUTION: keep the value when the lbc field is zero 1894 1901 CALL lbc_lnk( zenv, 'T', 1._wp, 'no0' ) … … 1930 1937 zri(:,:) = 0._wp 1931 1938 zrj(:,:) = 0._wp 1939 1940 !IF(ln_wd) THEN !extend the smoothed region to cover the W/D zones 1941 ! zsmth = -rn_wdld 1942 !ELSE 1943 zsmth = 0._wp ! The original form (only smooth ocean points) 1944 !ENDIF 1945 1932 1946 DO jj = 1, nlcj 1933 1947 DO ji = 1, nlci 1934 1948 iip1 = MIN( ji+1, nlci ) ! force zri = 0 on last line (ji=ncli+1 to jpi) 1935 1949 ijp1 = MIN( jj+1, nlcj ) ! force zrj = 0 on last raw (jj=nclj+1 to jpj) 1936 IF( (zenv(ji,jj) > 0._wp) .AND. (zenv(iip1,jj) > 0._wp)) THEN1950 IF( (zenv(ji,jj) > zsmth) .AND. (zenv(iip1,jj) > zsmth)) THEN 1937 1951 zri(ji,jj) = ( zenv(iip1,jj ) - zenv(ji,jj) ) / ( zenv(iip1,jj ) + zenv(ji,jj) ) 1938 1952 END IF 1939 IF( (zenv(ji,jj) > 0._wp) .AND. (zenv(ji,ijp1) > 0._wp)) THEN1953 IF( (zenv(ji,jj) > zsmth) .AND. (zenv(ji,ijp1) > zsmth)) THEN 1940 1954 zrj(ji,jj) = ( zenv(ji ,ijp1) - zenv(ji,jj) ) / ( zenv(ji ,ijp1) + zenv(ji,jj) ) 1941 1955 END IF … … 1960 1974 END DO ! End loop ! 1961 1975 ! ! ================ ! 1962 DO jj = 1, jpj 1963 DO ji = 1, jpi 1964 zenv(ji,jj) = MAX( zenv(ji,jj), rn_sbot_min ) ! set all points to avoid undefined scale value warnings 1965 END DO 1966 END DO 1976 IF(ln_wd) THEN 1977 DO jj = 1, jpj 1978 DO ji = 1, jpi 1979 !zenv(ji,jj) = MAX( zenv(ji,jj), -rn_wdld ) ! filt out land bathy data 1980 END DO 1981 END DO 1982 ELSE 1983 DO jj = 1, jpj 1984 DO ji = 1, jpi 1985 zenv(ji,jj) = MAX( zenv(ji,jj), rn_sbot_min ) ! set all points to avoid undefined scale value warnings 1986 END DO 1987 END DO 1988 ENDIF 1967 1989 ! 1968 1990 ! Envelope bathymetry saved in hbatt … … 1994 2016 IF(lwp) THEN 1995 2017 WRITE(numout,*) 1996 WRITE(numout,*) ' zgr_sco: minimum depth of the envelop topography set to : ', rn_sbot_min 2018 IF(.NOT.ln_wd) THEN 2019 WRITE(numout,*) ' zgr_sco: minimum depth of the envelop topography set to : ', rn_sbot_min 2020 ELSE 2021 WRITE(numout,*) ' zgr_sco: minimum positive depth of the envelop topography set to : ', rn_sbot_min 2022 WRITE(numout,*) ' zgr_sco: minimum negative depth of the envelop topography set to : ', -rn_wdld 2023 ENDIF 1997 2024 ENDIF 1998 2025 hbatu(:,:) = rn_sbot_min … … 2007 2034 END DO 2008 2035 END DO 2036 2037 IF(ln_wd) THEN !avoid the zero depth on T- (U-,V-,F-) points 2038 DO jj = 1, jpj 2039 DO ji = 1, jpi 2040 IF(ABS(hbatt(ji,jj)) < rn_wdmin1) & 2041 & hbatt(ji,jj) = SIGN(1._wp, hbatt(ji,jj)) * rn_wdmin1 2042 IF(ABS(hbatu(ji,jj)) < rn_wdmin1) & 2043 & hbatu(ji,jj) = SIGN(1._wp, hbatu(ji,jj)) * rn_wdmin1 2044 IF(ABS(hbatv(ji,jj)) < rn_wdmin1) & 2045 & hbatv(ji,jj) = SIGN(1._wp, hbatv(ji,jj)) * rn_wdmin1 2046 IF(ABS(hbatf(ji,jj)) < rn_wdmin1) & 2047 & hbatf(ji,jj) = SIGN(1._wp, hbatf(ji,jj)) * rn_wdmin1 2048 END DO 2049 END DO 2050 END IF 2009 2051 ! 2010 2052 ! Apply lateral boundary condition … … 2014 2056 DO ji = 1, jpi 2015 2057 IF( hbatu(ji,jj) == 0._wp ) THEN 2058 !No worries about the following line when ln_wd == .true. 2016 2059 IF( zhbat(ji,jj) == 0._wp ) hbatu(ji,jj) = rn_sbot_min 2017 2060 IF( zhbat(ji,jj) /= 0._wp ) hbatu(ji,jj) = zhbat(ji,jj) … … 2039 2082 2040 2083 !!bug: key_helsinki a verifer 2041 hift(:,:) = MIN( hift(:,:), hbatt(:,:) ) 2042 hifu(:,:) = MIN( hifu(:,:), hbatu(:,:) ) 2043 hifv(:,:) = MIN( hifv(:,:), hbatv(:,:) ) 2044 hiff(:,:) = MIN( hiff(:,:), hbatf(:,:) ) 2084 IF(.NOT.ln_wd) THEN 2085 hift(:,:) = MIN( hift(:,:), hbatt(:,:) ) 2086 hifu(:,:) = MIN( hifu(:,:), hbatu(:,:) ) 2087 hifv(:,:) = MIN( hifv(:,:), hbatv(:,:) ) 2088 hiff(:,:) = MIN( hiff(:,:), hbatf(:,:) ) 2089 END IF 2045 2090 2046 2091 IF( nprint == 1 .AND. lwp ) THEN … … 2087 2132 fsde3w(:,:,:) = gdep3w_0(:,:,:) 2088 2133 ! 2089 where (e3t_0 (:,:,:).eq.0.0) e3t_0(:,:,:) = 1.0 2090 where (e3u_0 (:,:,:).eq.0.0) e3u_0(:,:,:) = 1.0 2091 where (e3v_0 (:,:,:).eq.0.0) e3v_0(:,:,:) = 1.0 2092 where (e3f_0 (:,:,:).eq.0.0) e3f_0(:,:,:) = 1.0 2093 where (e3w_0 (:,:,:).eq.0.0) e3w_0(:,:,:) = 1.0 2094 where (e3uw_0 (:,:,:).eq.0.0) e3uw_0(:,:,:) = 1.0 2095 where (e3vw_0 (:,:,:).eq.0.0) e3vw_0(:,:,:) = 1.0 2134 IF(.NOT.ln_wd) THEN 2135 where (e3t_0 (:,:,:).eq.0.0) e3t_0(:,:,:) = 1.0 2136 where (e3u_0 (:,:,:).eq.0.0) e3u_0(:,:,:) = 1.0 2137 where (e3v_0 (:,:,:).eq.0.0) e3v_0(:,:,:) = 1.0 2138 where (e3f_0 (:,:,:).eq.0.0) e3f_0(:,:,:) = 1.0 2139 where (e3w_0 (:,:,:).eq.0.0) e3w_0(:,:,:) = 1.0 2140 where (e3uw_0 (:,:,:).eq.0.0) e3uw_0(:,:,:) = 1.0 2141 where (e3vw_0 (:,:,:).eq.0.0) e3vw_0(:,:,:) = 1.0 2142 END IF 2096 2143 2097 2144 #if defined key_agrif … … 2135 2182 IF( scobot(ji,jj) >= fsdept(ji,jj,jk) ) mbathy(ji,jj) = MAX( 2, jk ) 2136 2183 END DO 2137 IF( scobot(ji,jj) == 0._wp ) mbathy(ji,jj) = 0 2138 END DO 2139 END DO 2184 2185 IF(ln_wd) THEN 2186 IF( scobot(ji,jj) <= -(rn_wdld - rn_wdmin2)) THEN 2187 mbathy(ji,jj) = 0 2188 ELSEIF(scobot(ji,jj) <= rn_wdmin1) THEN 2189 mbathy(ji,jj) = 2 2190 ENDIF 2191 ELSE 2192 IF( scobot(ji,jj) == 0._wp ) mbathy(ji,jj) = 0 2193 ENDIF 2194 END DO 2195 END DO 2196 2140 2197 IF( nprint == 1 .AND. lwp ) WRITE(numout,*) ' MIN val mbathy h90 ', MINVAL( mbathy(:,:) ), & 2141 2198 & ' MAX ', MAXVAL( mbathy(:,:) ) … … 2257 2314 INTEGER :: ji, jj, jk ! dummy loop argument 2258 2315 REAL(wp) :: zcoeft, zcoefw ! temporary scalars 2316 REAL(wp) :: ztmpu, ztmpv, ztmpf 2317 REAL(wp) :: ztmpu1, ztmpv1, ztmpf1 2259 2318 ! 2260 2319 REAL(wp), POINTER, DIMENSION(:,:,:) :: z_gsigw3, z_gsigt3, z_gsi3w3 … … 2310 2369 DO ji = 1, jpim1 2311 2370 DO jj = 1, jpjm1 2371 ! extended for Wetting/Drying case 2372 ztmpu = hbatt(ji,jj)+hbatt(ji+1,jj) 2373 ztmpv = hbatt(ji,jj)+hbatt(ji,jj+1) 2374 ztmpf = hbatt(ji,jj)+hbatt(ji+1,jj)+hbatt(ji,jj+1)+hbatt(ji+1,jj+1) 2375 ztmpu1 = hbatt(ji,jj)*hbatt(ji+1,jj) 2376 ztmpv1 = hbatt(ji,jj)*hbatt(ji,jj+1) 2377 ztmpf1 = MIN(hbatt(ji,jj), hbatt(ji+1,jj), hbatt(ji,jj+1), hbatt(ji+1,jj+1)) * & 2378 & MAX(hbatt(ji,jj), hbatt(ji+1,jj), hbatt(ji,jj+1), hbatt(ji+1,jj+1)) 2312 2379 DO jk = 1, jpk 2313 z_esigtu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) ) & 2314 & / ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 2315 z_esigtv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk) ) & 2316 & / ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 2317 z_esigtf3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) & 2318 & + hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk)+hbatt(ji+1,jj+1)*z_esigt3(ji+1,jj+1,jk) ) & 2319 & / ( hbatt(ji,jj)+hbatt(ji+1,jj)+hbatt(ji,jj+1)+hbatt(ji+1,jj+1) ) 2320 z_esigwu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigw3(ji+1,jj,jk) ) & 2321 & / ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 2322 z_esigwv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigw3(ji,jj+1,jk) ) & 2323 & / ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 2380 IF((ztmpu1 < 0._wp.OR.ABS(ztmpu) < rn_wdmin1).AND.ln_wd) THEN 2381 z_esigtu3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji+1,jj,jk) ) 2382 z_esigwu3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji+1,jj,jk) ) 2383 ELSE 2384 z_esigtu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) ) & 2385 & / ztmpu 2386 z_esigwu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigw3(ji+1,jj,jk) ) & 2387 & / ztmpu 2388 END IF 2389 2390 IF((ztmpv1 < 0._wp.OR.ABS(ztmpv) < rn_wdmin1).AND.ln_wd) THEN 2391 z_esigtv3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji,jj+1,jk) ) 2392 z_esigwv3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji,jj+1,jk) ) 2393 ELSE 2394 z_esigtv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk) ) & 2395 & / ztmpv 2396 z_esigwv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigw3(ji,jj+1,jk) ) & 2397 & / ztmpv 2398 END IF 2399 2400 IF((ztmpf1 < 0._wp.OR.ABS(ztmpf) < rn_wdmin1).AND.ln_wd) THEN 2401 z_esigtf3(ji,jj,jk) = 0.25_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji+1,jj,jk) & 2402 & + z_esigt3(ji,jj+1,jk) + z_esigt3(ji+1,jj+1,jk) ) 2403 ELSE 2404 z_esigtf3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) & 2405 & + hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk)+hbatt(ji+1,jj+1)*z_esigt3(ji+1,jj+1,jk) ) & 2406 & / ztmpf 2407 END IF 2408 2409 ! 2410 !z_esigtu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) ) & 2411 ! & / ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 2412 !z_esigtv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk) ) & 2413 ! & / ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 2414 !z_esigtf3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) & 2415 ! & + hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk)+hbatt(ji+1,jj+1)*z_esigt3(ji+1,jj+1,jk) ) & 2416 ! & / ( hbatt(ji,jj)+hbatt(ji+1,jj)+hbatt(ji,jj+1)+hbatt(ji+1,jj+1) ) 2417 !z_esigwu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigw3(ji+1,jj,jk) ) & 2418 ! & / ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 2419 !z_esigwv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigw3(ji,jj+1,jk) ) & 2420 ! & / ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 2324 2421 ! 2325 2422 e3t_0(ji,jj,jk) = ( (hbatt(ji,jj)-rn_hc)*z_esigt3 (ji,jj,jk) + rn_hc/REAL(jpkm1,wp) ) … … 2361 2458 REAL(wp) :: zsmth ! smoothing around critical depth 2362 2459 REAL(wp) :: zzs, zzb ! Surface and bottom cell thickness in sigma space 2460 REAL(wp) :: ztmpu, ztmpv, ztmpf 2461 REAL(wp) :: ztmpu1, ztmpv1, ztmpf1 2363 2462 ! 2364 2463 REAL(wp), POINTER, DIMENSION(:,:,:) :: z_gsigw3, z_gsigt3, z_gsi3w3 … … 2439 2538 DO jj=1,jpj-1 2440 2539 2441 DO jk = 1, jpk 2442 z_esigtu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) ) / & 2443 ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 2444 z_esigtv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk) ) / & 2445 ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 2446 z_esigtf3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) + & 2447 hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk)+hbatt(ji+1,jj+1)*z_esigt3(ji+1,jj+1,jk) ) / & 2448 ( hbatt(ji,jj)+hbatt(ji+1,jj)+hbatt(ji,jj+1)+hbatt(ji+1,jj+1) ) 2449 z_esigwu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigw3(ji+1,jj,jk) ) / & 2450 ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 2451 z_esigwv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigw3(ji,jj+1,jk) ) / & 2452 ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 2540 ! extend to suit for Wetting/Drying case 2541 ztmpu = hbatt(ji,jj)+hbatt(ji+1,jj) 2542 ztmpv = hbatt(ji,jj)+hbatt(ji,jj+1) 2543 ztmpf = hbatt(ji,jj)+hbatt(ji+1,jj)+hbatt(ji,jj+1)+hbatt(ji+1,jj+1) 2544 ztmpu1 = hbatt(ji,jj)*hbatt(ji+1,jj) 2545 ztmpv1 = hbatt(ji,jj)*hbatt(ji,jj+1) 2546 ztmpf1 = MIN(hbatt(ji,jj), hbatt(ji+1,jj), hbatt(ji,jj+1), hbatt(ji+1,jj+1)) * & 2547 & MAX(hbatt(ji,jj), hbatt(ji+1,jj), hbatt(ji,jj+1), hbatt(ji+1,jj+1)) 2548 DO jk = 1, jpk 2549 IF((ztmpu1 < 0._wp.OR.ABS(ztmpu) < rn_wdmin1).AND.ln_wd) THEN 2550 z_esigtu3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji+1,jj,jk) ) 2551 z_esigwu3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji+1,jj,jk) ) 2552 ELSE 2553 z_esigtu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) ) & 2554 & / ztmpu 2555 z_esigwu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigw3(ji+1,jj,jk) ) & 2556 & / ztmpu 2557 END IF 2558 2559 IF((ztmpv1 < 0._wp.OR.ABS(ztmpv) < rn_wdmin1).AND.ln_wd) THEN 2560 z_esigtv3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji,jj+1,jk) ) 2561 z_esigwv3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji,jj+1,jk) ) 2562 ELSE 2563 z_esigtv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk) ) & 2564 & / ztmpv 2565 z_esigwv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigw3(ji,jj+1,jk) ) & 2566 & / ztmpv 2567 END IF 2568 2569 IF((ztmpf1 < 0._wp.OR.ABS(ztmpf) < rn_wdmin1).AND.ln_wd) THEN 2570 z_esigtf3(ji,jj,jk) = 0.25_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji+1,jj,jk) & 2571 & + z_esigt3(ji,jj+1,jk) + z_esigt3(ji+1,jj+1,jk) ) 2572 ELSE 2573 z_esigtf3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) & 2574 & + hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk)+hbatt(ji+1,jj+1)*z_esigt3(ji+1,jj+1,jk) ) & 2575 & / ztmpf 2576 END IF 2577 2578 !z_esigtu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) ) / & 2579 ! ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 2580 !z_esigtv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk) ) / & 2581 ! ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 2582 !z_esigtf3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) + & 2583 ! hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk)+hbatt(ji+1,jj+1)*z_esigt3(ji+1,jj+1,jk) ) / & 2584 ! ( hbatt(ji,jj)+hbatt(ji+1,jj)+hbatt(ji,jj+1)+hbatt(ji+1,jj+1) ) 2585 !z_esigwu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigw3(ji+1,jj,jk) ) / & 2586 ! ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 2587 !z_esigwv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigw3(ji,jj+1,jk) ) / & 2588 ! ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 2453 2589 2454 2590 e3t_0(ji,jj,jk)=(scosrf(ji,jj)+hbatt(ji,jj))*z_esigt3(ji,jj,jk)
Note: See TracChangeset
for help on using the changeset viewer.