- Timestamp:
- 2020-06-09T18:45:11+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/OCE/SBC/sbccpl.F90
r11715 r13080 48 48 USE lib_mpp ! distribued memory computing library 49 49 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 50 51 #if defined key_oasis3 52 USE mod_oasis, ONLY : OASIS_Sent, OASIS_ToRest, OASIS_SentOut, OASIS_ToRestOut 53 #endif 50 54 51 55 IMPLICIT NONE … … 152 156 INTEGER, PARAMETER :: jps_wlev = 32 ! water level 153 157 INTEGER, PARAMETER :: jps_fice1 = 33 ! first-order ice concentration (for semi-implicit coupling of atmos-ice fluxes) 154 INTEGER, PARAMETER :: jps_a_p = 34 ! meltpond area 158 INTEGER, PARAMETER :: jps_a_p = 34 ! meltpond area fraction 155 159 INTEGER, PARAMETER :: jps_ht_p = 35 ! meltpond thickness 156 160 INTEGER, PARAMETER :: jps_kice = 36 ! sea ice effective conductivity … … 159 163 160 164 INTEGER, PARAMETER :: jpsnd = 38 ! total number of fields sent 165 166 #if ! defined key_oasis3 167 ! Dummy variables to enable compilation when oasis3 is not being used 168 INTEGER :: OASIS_Sent = -1 169 INTEGER :: OASIS_SentOut = -1 170 INTEGER :: OASIS_ToRest = -1 171 INTEGER :: OASIS_ToRestOut = -1 172 #endif 161 173 162 174 ! !!** namelist namsbc_cpl ** … … 184 196 LOGICAL :: ln_usecplmask ! use a coupling mask file to merge data received from several models 185 197 ! -> file cplmask.nc with the float variable called cplmask (jpi,jpj,nn_cplmodel) 198 LOGICAL :: ln_scale_ice_fluxes ! Scale sea ice fluxes by the sea ice fractions at the previous coupling point 186 199 TYPE :: DYNARR 187 200 REAL(wp), POINTER, DIMENSION(:,:,:) :: z3 … … 250 263 NAMELIST/namsbc_cpl/ sn_snd_temp , sn_snd_alb , sn_snd_thick, sn_snd_crt , sn_snd_co2 , & 251 264 & sn_snd_ttilyr, sn_snd_cond , sn_snd_mpnd , sn_snd_sstfrz, sn_snd_thick1, & 265 & ln_scale_ice_fluxes, & 252 266 & sn_snd_ifrac , sn_snd_crtw , sn_snd_wlev , sn_rcv_hsig , sn_rcv_phioc, & 253 267 & sn_rcv_w10m , sn_rcv_taumod, sn_rcv_tau , sn_rcv_dqnsdt, sn_rcv_qsr , & … … 279 293 ENDIF 280 294 IF( lwp .AND. ln_cpl ) THEN ! control print 295 WRITE(numout,*)' ln_scale_ice_fluxes = ', ln_scale_ice_fluxes 281 296 WRITE(numout,*)' received fields (mutiple ice categogies)' 282 297 WRITE(numout,*)' 10m wind module = ', TRIM(sn_rcv_w10m%cldes ), ' (', TRIM(sn_rcv_w10m%clcat ), ')' … … 815 830 END SELECT 816 831 832 ! Initialise ice fractions from last coupling time to zero 833 a_i_last_couple(:,:,:) = 0._wp 834 835 817 836 ! ! ------------------------- ! 818 837 ! ! Ice Meltponds ! … … 1243 1262 IF( srcv(jpr_co2)%laction ) atm_co2(:,:) = frcv(jpr_co2)%z3(:,:,1) 1244 1263 ! 1245 ! ! ================== !1246 ! ! ice skin temp. !1247 ! ! ================== !1248 #if defined key_si31249 ! needed by Met Office1250 IF( srcv(jpr_ts_ice)%laction ) THEN1251 WHERE ( frcv(jpr_ts_ice)%z3(:,:,:) > 0.0 ) ; tsfc_ice(:,:,:) = 0.01252 ELSEWHERE( frcv(jpr_ts_ice)%z3(:,:,:) < -60. ) ; tsfc_ice(:,:,:) = -60.1253 ELSEWHERE ; tsfc_ice(:,:,:) = frcv(jpr_ts_ice)%z3(:,:,:)1254 END WHERE1255 ENDIF1256 #endif1257 1264 ! ! ========================= ! 1258 1265 ! ! Mean Sea Level Pressure ! (taum) … … 1630 1637 !! sprecip solid precipitation over the ocean 1631 1638 !!---------------------------------------------------------------------- 1632 REAL(wp), INTENT(in) , DIMENSION(:,:) :: picefr ! ice fraction [0 to 1]1633 ! !! ! optional arguments, used only in 'mixed oce-ice' case1634 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: palbi ! all skies ice albedo1635 REAL(wp), INTENT(in) , DIMENSION(:,: ), OPTIONAL :: psst ! sea surface temperature [Celsius]1636 REAL(wp), INTENT(in ), DIMENSION(:,:,:), OPTIONAL :: pist ! ice surface temperature [Kelvin]1637 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: phs ! snow depth [m]1638 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: phi ! ice thickness [m]1639 REAL(wp), INTENT(in) , DIMENSION(:,:) :: picefr ! ice fraction [0 to 1] 1640 ! !! ! optional arguments, used only in 'mixed oce-ice' case or for Met-Office coupling 1641 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: palbi ! all skies ice albedo 1642 REAL(wp), INTENT(in) , DIMENSION(:,: ), OPTIONAL :: psst ! sea surface temperature [Celsius] 1643 REAL(wp), INTENT(inout), DIMENSION(:,:,:), OPTIONAL :: pist ! ice surface temperature [Kelvin] => inout for Met-Office 1644 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: phs ! snow depth [m] 1645 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: phi ! ice thickness [m] 1639 1646 ! 1640 1647 INTEGER :: ji, jj, jl ! dummy loop index … … 1643 1650 REAL(wp), DIMENSION(jpi,jpj) :: zemp_tot, zemp_ice, zemp_oce, ztprecip, zsprecip , zevap_oce, zdevap_ice 1644 1651 REAL(wp), DIMENSION(jpi,jpj) :: zqns_tot, zqns_oce, zqsr_tot, zqsr_oce, zqprec_ice, zqemp_oce, zqemp_ice 1645 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zqns_ice, zqsr_ice, zdqns_ice, zqevap_ice, zevap_ice !!gm , zfrqsr_tr_i 1652 REAL(wp), DIMENSION(jpi,jpj) :: zevap_ice_total 1653 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zqns_ice, zqsr_ice, zdqns_ice, zqevap_ice, zevap_ice, zqtr_ice_top, ztsu 1646 1654 !!---------------------------------------------------------------------- 1647 1655 ! … … 1663 1671 ztprecip(:,:) = frcv(jpr_rain)%z3(:,:,1) + zsprecip(:,:) ! May need to ensure positive here 1664 1672 zemp_tot(:,:) = frcv(jpr_tevp)%z3(:,:,1) - ztprecip(:,:) 1665 zemp_ice(:,:) = ( frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1) ) * picefr(:,:)1666 1673 CASE( 'oce and ice' ) ! received fields: jpr_sbpr, jpr_semp, jpr_oemp, jpr_ievp 1667 1674 zemp_tot(:,:) = ziceld(:,:) * frcv(jpr_oemp)%z3(:,:,1) + picefr(:,:) * frcv(jpr_sbpr)%z3(:,:,1) 1668 1675 zemp_ice(:,:) = frcv(jpr_semp)%z3(:,:,1) * picefr(:,:) 1669 1676 zsprecip(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_semp)%z3(:,:,1) 1670 ztprecip(:,:) = frcv(jpr_semp)%z3(:,:,1) - frcv(jpr_sbpr)%z3(:,:,1) + zsprecip(:,:) 1677 ztprecip(:,:) = frcv(jpr_semp)%z3(:,:,1) - frcv(jpr_sbpr)%z3(:,:,1) + zsprecip(:,:) 1671 1678 CASE( 'none' ) ! Not available as for now: needs additional coding below when computing zevap_oce 1672 1679 ! ! since fields received are not defined with none option … … 1675 1682 1676 1683 #if defined key_si3 1684 1685 ! --- evaporation over ice (kg/m2/s) --- ! 1686 zevap_ice_total(:,:) = 0._wp 1687 IF (sn_rcv_emp%clcat == 'yes') THEN 1688 DO jl=1,jpl 1689 IF (ln_scale_ice_fluxes) THEN 1690 zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl) * a_i_last_couple(:,:,jl) 1691 ELSE 1692 zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl) 1693 ENDIF 1694 zevap_ice_total(:,:) = zevap_ice_total(:,:) + zevap_ice(:,:,jl) 1695 ENDDO 1696 ELSE 1697 zevap_ice(:,:,1) = frcv(jpr_ievp)%z3(:,:,1 ) 1698 zevap_ice_total(:,:) = zevap_ice(:,:,1) 1699 ENDIF 1700 1701 IF ( TRIM( sn_rcv_emp%cldes ) == 'conservative' ) THEN 1702 ! For conservative case zemp_ice has not been defined yet. Do it now. 1703 zemp_ice(:,:) = zevap_ice_total(:,:) - frcv(jpr_snow)%z3(:,:,1) * picefr(:,:) 1704 END IF 1705 1677 1706 ! zsnw = snow fraction over ice after wind blowing (=picefr if no blowing) 1678 1707 zsnw(:,:) = 0._wp ; CALL ice_thd_snwblow( ziceld, zsnw ) … … 1683 1712 1684 1713 ! --- evaporation over ocean (used later for qemp) --- ! 1685 zevap_oce(:,:) = frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) 1686 1687 ! --- evaporation over ice (kg/m2/s) --- ! 1688 DO jl=1,jpl 1689 IF (sn_rcv_emp%clcat == 'yes') THEN ; zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl) 1690 ELSE ; zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,1 ) ; ENDIF 1691 ENDDO 1714 zevap_oce(:,:) = frcv(jpr_tevp)%z3(:,:,1) - zevap_ice_total(:,:) 1715 1716 1717 1692 1718 1693 1719 ! since the sensitivity of evap to temperature (devap/dT) is not prescribed by the atmosphere, we set it to 0 … … 1727 1753 sprecip (:,:) = zsprecip (:,:) 1728 1754 tprecip (:,:) = ztprecip (:,:) 1729 evap_ice(:,:,:) = zevap_ice(:,:,:) 1755 IF (ln_scale_ice_fluxes) THEN 1756 ! Convert from grid box means to sea ice means 1757 WHERE( a_i(:,:,:) > 0.0_wp ) evap_ice(:,:,:) = zevap_ice(:,:,:) / a_i(:,:,:) 1758 WHERE( a_i(:,:,:) <= 0.0_wp ) evap_ice(:,:,:) = 0.0 1759 ELSE 1760 evap_ice(:,:,:) = zevap_ice(:,:,:) 1761 ENDIF 1730 1762 DO jl = 1, jpl 1731 1763 devap_ice(:,:,jl) = zdevap_ice(:,:) … … 1774 1806 IF( iom_use('snow_ao_cea') ) CALL iom_put( 'snow_ao_cea' , sprecip(:,:) * ( 1._wp - zsnw(:,:) ) ) ! Snow over ice-free ocean (cell average) 1775 1807 IF( iom_use('snow_ai_cea') ) CALL iom_put( 'snow_ai_cea' , sprecip(:,:) * zsnw(:,:) ) ! Snow over sea-ice (cell average) 1776 IF( iom_use('subl_ai_cea') ) CALL iom_put( 'subl_ai_cea' , frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * tmask(:,:,1) ) ! Sublimation over sea-ice (cell average)1808 IF( iom_use('subl_ai_cea') ) CALL iom_put( 'subl_ai_cea' , zevap_ice_total(:,:) * tmask(:,:,1) ) ! Sublimation over sea-ice (cell average) 1777 1809 IF( iom_use('evap_ao_cea') ) CALL iom_put( 'evap_ao_cea' , ( frcv(jpr_tevp)%z3(:,:,1) & 1778 & - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) ) * tmask(:,:,1) ) ! ice-free oce evap (cell average)1810 & - zevap_ice_total(:,:) ) * tmask(:,:,1) ) ! ice-free oce evap (cell average) 1779 1811 ! note: runoff output is done in sbcrnf (which includes icebergs too) and iceshelf output is done in sbcisf 1780 1812 ! … … 1784 1816 CASE( 'oce only' ) ! the required field is directly provided 1785 1817 zqns_tot(:,:) = frcv(jpr_qnsoce)%z3(:,:,1) 1818 ! For Met Office sea ice non-solar fluxes are already delt with by JULES so setting to zero 1819 ! here so the only flux is the ocean only one. 1820 zqns_ice(:,:,:) = 0._wp 1786 1821 CASE( 'conservative' ) ! the required fields are directly provided 1787 1822 zqns_tot(:,:) = frcv(jpr_qnsmix)%z3(:,:,1) … … 1847 1882 & + zsprecip(:,:) * ( 1._wp - zsnw ) * ( zcptsnw (:,:) - rLfus ) ! solid precip over ocean + snow melting 1848 1883 zqemp_ice(:,:) = zsprecip(:,:) * zsnw * ( zcptsnw (:,:) - rLfus ) ! solid precip over ice (qevap_ice=0 since atm. does not take it into account) 1849 !! zqemp_ice(:,:) = - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * zcptsnw (:,:) & ! ice evap1850 !! & + zsprecip(:,:) * zsnw * zqprec_ice(:,:) * r1_rhos ! solid precip over ice1851 1884 1852 1885 ! --- total non solar flux (including evap/precip) --- ! … … 1900 1933 IF ( srcv(jpr_icb)%laction ) CALL iom_put('hflx_icb_cea' , - frcv(jpr_icb)%z3(:,:,1) * rLfus ) ! latent heat from icebergs melting 1901 1934 IF ( iom_use('hflx_rain_cea') ) CALL iom_put('hflx_rain_cea' , ( tprecip(:,:) - sprecip(:,:) ) * zcptrain(:,:) ) ! heat flux from rain (cell average) 1902 IF ( iom_use('hflx_evap_cea') ) CALL iom_put('hflx_evap_cea' , ( frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1)&1903 & * picefr(:,:) )* zcptn(:,:) * tmask(:,:,1) ) ! heat flux from evap (cell average)1935 IF ( iom_use('hflx_evap_cea') ) CALL iom_put('hflx_evap_cea' , ( frcv(jpr_tevp)%z3(:,:,1) - zevap_ice_total(:,:) ) & 1936 * zcptn(:,:) * tmask(:,:,1) ) ! heat flux from evap (cell average) 1904 1937 IF ( iom_use('hflx_snow_cea') ) CALL iom_put('hflx_snow_cea' , sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) ) ! heat flux from snow (cell average) 1905 1938 IF ( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea', sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) & … … 1914 1947 CASE( 'oce only' ) 1915 1948 zqsr_tot(:,: ) = MAX( 0._wp , frcv(jpr_qsroce)%z3(:,:,1) ) 1949 ! For Met Office sea ice solar fluxes are already delt with by JULES so setting to zero 1950 ! here so the only flux is the ocean only one. 1951 zqsr_ice(:,:,:) = 0._wp 1916 1952 CASE( 'conservative' ) 1917 1953 zqsr_tot(:,: ) = frcv(jpr_qsrmix)%z3(:,:,1) … … 1992 2028 ENDDO 1993 2029 ENDIF 2030 CASE( 'none' ) 2031 zdqns_ice(:,:,:) = 0._wp 1994 2032 END SELECT 1995 2033 … … 2007 2045 ! ! ========================= ! 2008 2046 CASE ('coupled') 2009 qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:) 2010 qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:) 2047 IF (ln_scale_ice_fluxes) THEN 2048 WHERE( a_i(:,:,:) > 0.0_wp ) qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) 2049 WHERE( a_i(:,:,:) <= 0.0_wp ) qml_ice(:,:,:) = 0.0_wp 2050 WHERE( a_i(:,:,:) > 0.0_wp ) qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) 2051 WHERE( a_i(:,:,:) <= 0.0_wp ) qcn_ice(:,:,:) = 0.0_wp 2052 ELSE 2053 qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:) 2054 qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:) 2055 ENDIF 2011 2056 END SELECT 2012 2057 ! … … 2017 2062 ! 2018 2063 ! ! ===> used prescribed cloud fraction representative for polar oceans in summer (0.81) 2019 ztri = 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ! surface transmission parameter(Grenfell Maykut 77)2064 ztri = 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ! surface transmission when hi>10cm (Grenfell Maykut 77) 2020 2065 ! 2021 qtr_ice_top(:,:,:) = ztri * qsr_ice(:,:,:) 2022 WHERE( phs(:,:,:) >= 0.0_wp ) qtr_ice_top(:,:,:) = 0._wp ! snow fully opaque 2023 WHERE( phi(:,:,:) <= 0.1_wp ) qtr_ice_top(:,:,:) = qsr_ice(:,:,:) ! thin ice transmits all solar radiation 2066 WHERE ( phs(:,:,:) <= 0._wp .AND. phi(:,:,:) < 0.1_wp ) ! linear decrease from hi=0 to 10cm 2067 zqtr_ice_top(:,:,:) = qsr_ice(:,:,:) * ( ztri + ( 1._wp - ztri ) * ( 1._wp - phi(:,:,:) * 10._wp ) ) 2068 ELSEWHERE( phs(:,:,:) <= 0._wp .AND. phi(:,:,:) >= 0.1_wp ) ! constant (ztri) when hi>10cm 2069 zqtr_ice_top(:,:,:) = qsr_ice(:,:,:) * ztri 2070 ELSEWHERE ! zero when hs>0 2071 zqtr_ice_top(:,:,:) = 0._wp 2072 END WHERE 2024 2073 ! 2025 2074 ELSEIF( ln_cndflx .AND. .NOT.ln_cndemulate ) THEN !== conduction flux as surface forcing ==! … … 2027 2076 ! ! ===> here we must receive the qtr_ice_top array from the coupler 2028 2077 ! for now just assume zero (fully opaque ice) 2029 qtr_ice_top(:,:,:) = 0._wp 2078 zqtr_ice_top(:,:,:) = 0._wp 2079 ! 2080 ENDIF 2081 ! 2082 IF( ln_mixcpl ) THEN 2083 DO jl=1,jpl 2084 qtr_ice_top(:,:,jl) = qtr_ice_top(:,:,jl) * xcplmask(:,:,0) + zqtr_ice_top(:,:,jl) * zmsk(:,:) 2085 ENDDO 2086 ELSE 2087 qtr_ice_top(:,:,:) = zqtr_ice_top(:,:,:) 2088 ENDIF 2089 ! ! ================== ! 2090 ! ! ice skin temp. ! 2091 ! ! ================== ! 2092 ! needed by Met Office 2093 IF( srcv(jpr_ts_ice)%laction ) THEN 2094 WHERE ( frcv(jpr_ts_ice)%z3(:,:,:) > 0.0 ) ; ztsu(:,:,:) = 0.0 + rt0 2095 ELSEWHERE( frcv(jpr_ts_ice)%z3(:,:,:) < -60. ) ; ztsu(:,:,:) = -60. + rt0 2096 ELSEWHERE ; ztsu(:,:,:) = frcv(jpr_ts_ice)%z3(:,:,:) + rt0 2097 END WHERE 2098 ! 2099 IF( ln_mixcpl ) THEN 2100 DO jl=1,jpl 2101 pist(:,:,jl) = pist(:,:,jl) * xcplmask(:,:,0) + ztsu(:,:,jl) * zmsk(:,:) 2102 ENDDO 2103 ELSE 2104 pist(:,:,:) = ztsu(:,:,:) 2105 ENDIF 2030 2106 ! 2031 2107 ENDIF … … 2047 2123 INTEGER, INTENT(in) :: kt 2048 2124 ! 2125 ! 2049 2126 INTEGER :: ji, jj, jl ! dummy loop indices 2050 2127 INTEGER :: isec, info ! local integer 2051 2128 REAL(wp) :: zumax, zvmax 2052 2129 REAL(wp), DIMENSION(jpi,jpj) :: zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1 2053 REAL(wp), DIMENSION(jpi,jpj,jpl) :: ztmp3, ztmp4 2130 REAL(wp), DIMENSION(jpi,jpj,jpl) :: ztmp3, ztmp4 2131 2054 2132 !!---------------------------------------------------------------------- 2055 2133 ! … … 2193 2271 ENDIF 2194 2272 2273 ! If this coupling was successful then save ice fraction for use between coupling points. 2274 ! This is needed for some calculations where the ice fraction at the last coupling point 2275 ! is needed. 2276 IF( info == OASIS_Sent .OR. info == OASIS_ToRest .OR. & 2277 & info == OASIS_SentOut .OR. info == OASIS_ToRestOut ) THEN 2278 IF ( sn_snd_thick%clcat == 'yes' ) THEN 2279 a_i_last_couple(:,:,1:jpl) = a_i(:,:,1:jpl) 2280 ENDIF 2281 ENDIF 2282 2195 2283 IF( ssnd(jps_fice1)%laction ) THEN 2196 2284 SELECT CASE( sn_snd_thick1%clcat ) … … 2250 2338 ! ! Ice melt ponds ! 2251 2339 ! ! ------------------------- ! 2252 ! needed by Met Office 2340 ! needed by Met Office - 1) fraction of ponded ice; 2) local/actual pond depth 2253 2341 IF( ssnd(jps_a_p)%laction .OR. ssnd(jps_ht_p)%laction ) THEN 2254 2342 SELECT CASE( sn_snd_mpnd%cldes) … … 2256 2344 SELECT CASE( sn_snd_mpnd%clcat ) 2257 2345 CASE( 'yes' ) 2258 ztmp3(:,:,1:jpl) = a_ip(:,:,1:jpl) 2259 ztmp4(:,:,1:jpl) = v_ip(:,:,1:jpl) 2346 2347 ztmp3(:,:,1:jpl) = a_ip_eff(:,:,1:jpl) 2348 ztmp4(:,:,1:jpl) = h_ip(:,:,1:jpl) 2349 2260 2350 CASE( 'no' ) 2261 2351 ztmp3(:,:,:) = 0.0 2262 2352 ztmp4(:,:,:) = 0.0 2263 2353 DO jl=1,jpl 2264 ztmp3(:,:,1) = ztmp3(:,:,1) + a_ip (:,:,jpl)2265 ztmp4(:,:,1) = ztmp4(:,:,1) + v_ip(:,:,jpl)2354 ztmp3(:,:,1) = ztmp3(:,:,1) + a_ip_frac(:,:,jpl) 2355 ztmp4(:,:,1) = ztmp4(:,:,1) + h_ip(:,:,jpl) 2266 2356 ENDDO 2267 2357 CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_mpnd%clcat' )
Note: See TracChangeset
for help on using the changeset viewer.