Changeset 5144 for LMDZ6/branches/Amaury_dev/libf/phylmd/wx_pbl_mod.F90
- Timestamp:
- Jul 29, 2024, 11:01:04 PM (3 months ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/branches/Amaury_dev/libf/phylmd/wx_pbl_mod.F90
r5143 r5144 1 1 MODULE wx_pbl_mod 2 2 3 ! Split Planetary Boundary Layer4 5 ! This module manages the splitting of the boundary layer between two regions; the (w) 6 ! region (inside cold pools) and the (x) region (outside cold pools)3 ! Split Planetary Boundary Layer 4 5 ! This module manages the splitting of the boundary layer between two regions; the (w) 6 ! region (inside cold pools) and the (x) region (outside cold pools) 7 7 8 8 USE dimphy … … 12 12 CONTAINS 13 13 14 !****************************************************************************************15 16 SUBROUTINE wx_pbl0_merge(knon, ypplay, ypaprs,&17 sigw, dTs_forcing, dqs_forcing,&18 19 20 21 22 23 24 25 26 27 28 29 30 14 !**************************************************************************************** 15 16 SUBROUTINE wx_pbl0_merge(knon, ypplay, ypaprs, & 17 sigw, dTs_forcing, dqs_forcing, & 18 yt_x, yt_w, yq_x, yq_w, & 19 yu_x, yu_w, yv_x, yv_w, & 20 ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w, & 21 ycdragm_x, ycdragm_w, & 22 AcoefT_x, AcoefT_w, AcoefQ_x, AcoefQ_w, & 23 AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, & 24 BcoefT_x, BcoefT_w, BcoefQ_x, BcoefQ_w, & 25 BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, & 26 AcoefT, AcoefQ, AcoefU, AcoefV, & 27 BcoefT, BcoefQ, BcoefU, BcoefV, & 28 ycdragh, ycdragq, ycdragm, & 29 yt1, yq1, yu1, yv1 & 30 ) 31 31 32 32 USE wx_pbl_var_mod 33 33 34 USE lmdz_print_control, ONLY: prt_level, lunout34 USE lmdz_print_control, ONLY: prt_level, lunout 35 35 USE indice_sol_mod, ONLY: is_oce 36 36 USE lmdz_clesphys 37 USE lmdz_ YOETHF37 USE lmdz_yoethf 38 38 USE lmdz_fcttre, ONLY: foeew, foede, qsats, qsatl, dqsats, dqsatl, thermcep 39 40 INCLUDE "YOMCST.h" 41 42 INTEGER, INTENT(IN) :: knon ! number of grid cells 43 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) 44 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) 45 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area 46 REAL, DIMENSION(knon), INTENT(IN) :: dTs_forcing ! forced temperature difference (w)-(x) 47 REAL, DIMENSION(knon), INTENT(IN) :: dqs_forcing ! forced humidity difference (w)-(x) 48 REAL, DIMENSION(knon,klev), INTENT(IN) :: yt_x, yt_w, yq_x, yq_w 49 REAL, DIMENSION(knon,klev), INTENT(IN) :: yu_x, yu_w, yv_x, yv_w 50 REAL, DIMENSION(knon), INTENT(IN) :: ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w 51 REAL, DIMENSION(knon), INTENT(IN) :: ycdragm_x, ycdragm_w 52 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w, AcoefQ_x, AcoefQ_w 53 REAL, DIMENSION(knon), INTENT(IN) :: AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w 54 REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w, BcoefQ_x, BcoefQ_w 55 REAL, DIMENSION(knon), INTENT(IN) :: BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w 56 REAL, DIMENSION(knon), INTENT(OUT) :: AcoefT, AcoefQ, AcoefU, AcoefV 57 REAL, DIMENSION(knon), INTENT(OUT) :: BcoefT, BcoefQ, BcoefU, BcoefV 58 REAL, DIMENSION(knon), INTENT(OUT) :: ycdragh, ycdragq, ycdragm 59 REAL, DIMENSION(knon), INTENT(OUT) :: yt1, yq1, yu1, yv1 ! Apparent T, q, u, v at first level, as 60 !seen by surface modules 61 62 ! Local variables 63 INTEGER :: j 64 REAL :: dd_Kh 65 REAL :: dd_Kq 66 REAL :: dd_Km 67 REAL :: dd_u 68 REAL :: dd_v 69 REAL :: dd_t 70 REAL :: dd_q 71 72 REAL :: LambdaTs, LambdaQs, LambdaUs, LambdaVs 73 74 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region 75 76 77 sigx(1:knon) = 1.-sigw(1:knon) 78 79 80 DO j=1,knon 81 82 83 ! Compute w-x differences 84 dd_t = yt_w(j,1) - yt_x(j,1) 85 dd_q = yq_w(j,1) - yq_x(j,1) 86 dd_u = yu_w(j,1) - yu_x(j,1) 87 dd_v = yv_w(j,1) - yv_x(j,1) 88 89 ! Merged exchange coefficients 90 dd_Kh = Kech_h_w(j) - Kech_h_x(j) 91 dd_Kq = Kech_q_w(j) - Kech_q_x(j) 92 dd_Km = Kech_m_w(j) - Kech_m_x(j) 93 94 LambdaTs = dd_KTp(j)/Kech_Tp(j) 95 LambdaQs = dd_KQs(j)/Kech_Qs(j) 96 LambdaUs = dd_KUp(j)/Kech_Up(j) 97 LambdaVs = dd_KVp(j)/Kech_Vp(j) 98 99 ! Calcul des coef A, B \'equivalents dans la couche 1 100 101 ! The dTs_forcing and dqs_forcing terms are added for diagnostic purpose ; they should be zero in normal operation. 102 AcoefT(j) = AcoefT_x(j) + sigw(j)*(1.+sigx(j)*LambdaTs)*(dd_AT(j) - C_p(j)*dTs_forcing(j)) 103 AcoefQ(j) = AcoefQ_x(j) + sigw(j)*(1.+sigx(j)*LambdaQs)*(dd_AQ(j) - dqs_forcing(j)) 104 AcoefU(j) = AcoefU_x(j) + sigw(j)*(1.+sigx(j)*LambdaUs)*dd_AU(j) 105 AcoefV(j) = AcoefV_x(j) + sigw(j)*(1.+sigx(j)*LambdaVs)*dd_AV(j) 106 107 108 !! BcoefT(j) = (sigw(j)*Kech_h_w(j)*Kech_T_pw(j)*BcoefT_w(j) + & 109 !! sigx(j)*Kech_h_x(j)*Kech_T_px(j)*BcoefT_x(j) )/(Kech_h(j)*Kech_Tp(j)) 110 !! BcoefQ(j) = (sigw(j)*Kech_q_w(j)*Kech_Q_pw(j)*BcoefQ_w(j) + & 111 !! sigx(j)*Kech_q_x(j)*Kech_Q_px(j)*BcoefQ_x(j) )/(Kech_q(j)*Kech_Qp(j)) 112 !! BcoefU(j) = (sigw(j)*Kech_m_w(j)*Kech_U_pw(j)*BcoefU_w(j) + & 113 !! sigx(j)*Kech_m_x(j)*Kech_U_px(j)*BcoefU_x(j) )/(Kech_m(j)*Kech_Up(j)) 114 !! BcoefV(j) = (sigw(j)*Kech_m_w(j)*Kech_V_pw(j)*BcoefV_w(j) + & 115 !! sigx(j)*Kech_m_x(j)*Kech_V_px(j)*BcoefV_x(j) )/(Kech_m(j)*Kech_Vp(j)) 116 117 !! Print *,'YYYYpbl0: BcoefT_x, sigw, sigx, dd_Kh, dd_KTp, Kech_h_w ', & 118 !! BcoefT_x, sigw, sigx, dd_Kh, dd_KTp, Kech_h_w 119 !! Print *,'YYYYpbl0: Kech_T_pw, dd_BT, Kech_h, Kech_Tp ', & 120 !! Kech_T_pw, dd_BT, Kech_h, Kech_Tp 121 BcoefT(j) = BcoefT_x(j) + sigw(j)*(sigx(j)*dd_Kh*dd_KTp(j)*BcoefT_x(j) + & 122 Kech_h_w(j)*Kech_T_pw(j)*dd_BT(j))/(Kech_h(j)*Kech_Tp(j)) 123 BcoefQ(j) = BcoefQ_x(j) + sigw(j)*(sigx(j)*dd_Kq*dd_KQs(j)*BcoefQ_x(j) + & 124 Kech_q_w(j)*Kech_Q_sw(j)*dd_BQ(j))/(Kech_q(j)*Kech_Qs(j)) 125 BcoefU(j) = BcoefU_x(j) + sigw(j)*(sigx(j)*dd_Km*dd_KUp(j)*BcoefU_x(j) + & 126 Kech_m_w(j)*Kech_U_pw(j)*dd_BU(j))/(Kech_m(j)*Kech_Up(j)) 127 BcoefV(j) = BcoefV_x(j) + sigw(j)*(sigx(j)*dd_Km*dd_KVp(j)*BcoefV_x(j) + & 128 Kech_m_w(j)*Kech_V_pw(j)*dd_BV(j))/(Kech_m(j)*Kech_Vp(j)) 129 !>jyg 130 131 132 ! Calcul des cdrag \'equivalents dans la couche 133 134 ycdragm(j) = ycdragm_x(j) + sigw(j)*dd_Cdragm(j) 135 ycdragh(j) = ycdragh_x(j) + sigw(j)*dd_Cdragh(j) 136 ycdragq(j) = ycdragq_x(j) + sigw(j)*dd_Cdragq(j) 137 138 ! Calcul de T, q, u et v \'equivalents dans la couche 1 139 !! yt1(j) = yt_x(j,1) + sigw(j)*dd_t*(1.+sigx(j)*dd_Kh/KCT) 140 !! yq1(j) = yq_x(j,1) + sigw(j)*dd_q*(1.+sigx(j)*dd_Kh/KCQ) 141 !! yu1(j) = yu_x(j,1) + sigw(j)*dd_u*(1.+sigx(j)*dd_Km/KCU) 142 !! yv1(j) = yv_x(j,1) + sigw(j)*dd_v*(1.+sigx(j)*dd_Km/KCV) 143 yt1(j) = yt_x(j,1) + sigw(j)*dd_t 144 yq1(j) = yq_x(j,1) + sigw(j)*dd_q 145 yu1(j) = yu_x(j,1) + sigw(j)*dd_u 146 yv1(j) = yv_x(j,1) + sigw(j)*dd_v 147 148 149 ENDDO 150 151 152 153 END SUBROUTINE wx_pbl0_merge 154 155 SUBROUTINE wx_pbl_dts_merge(knon, dtime, ypplay, ypaprs, & 156 sigw, beta, wcstar, wdens, & 157 AT_x, AT_w, & 158 BT_x, BT_w, & 159 AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, & 160 AcoefT, AcoefQ, BcoefT, BcoefQ, & 161 HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, & 162 phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0, & 163 g_T, g_Q, & 164 Gamma_phiT, Gamma_phiQ, & 165 dTs_ins, dqsatsrf_ins & 166 ) 39 USE lmdz_yomcst 40 41 IMPLICIT NONE 42 43 INTEGER, INTENT(IN) :: knon ! number of grid cells 44 REAL, DIMENSION(knon, klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) 45 REAL, DIMENSION(knon, klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) 46 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area 47 REAL, DIMENSION(knon), INTENT(IN) :: dTs_forcing ! forced temperature difference (w)-(x) 48 REAL, DIMENSION(knon), INTENT(IN) :: dqs_forcing ! forced humidity difference (w)-(x) 49 REAL, DIMENSION(knon, klev), INTENT(IN) :: yt_x, yt_w, yq_x, yq_w 50 REAL, DIMENSION(knon, klev), INTENT(IN) :: yu_x, yu_w, yv_x, yv_w 51 REAL, DIMENSION(knon), INTENT(IN) :: ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w 52 REAL, DIMENSION(knon), INTENT(IN) :: ycdragm_x, ycdragm_w 53 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w, AcoefQ_x, AcoefQ_w 54 REAL, DIMENSION(knon), INTENT(IN) :: AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w 55 REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w, BcoefQ_x, BcoefQ_w 56 REAL, DIMENSION(knon), INTENT(IN) :: BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w 57 REAL, DIMENSION(knon), INTENT(OUT) :: AcoefT, AcoefQ, AcoefU, AcoefV 58 REAL, DIMENSION(knon), INTENT(OUT) :: BcoefT, BcoefQ, BcoefU, BcoefV 59 REAL, DIMENSION(knon), INTENT(OUT) :: ycdragh, ycdragq, ycdragm 60 REAL, DIMENSION(knon), INTENT(OUT) :: yt1, yq1, yu1, yv1 ! Apparent T, q, u, v at first level, as 61 !seen by surface modules 62 63 ! Local variables 64 INTEGER :: j 65 REAL :: dd_Kh 66 REAL :: dd_Kq 67 REAL :: dd_Km 68 REAL :: dd_u 69 REAL :: dd_v 70 REAL :: dd_t 71 REAL :: dd_q 72 73 REAL :: LambdaTs, LambdaQs, LambdaUs, LambdaVs 74 75 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region 76 77 sigx(1:knon) = 1. - sigw(1:knon) 78 79 DO j = 1, knon 80 81 82 ! Compute w-x differences 83 dd_t = yt_w(j, 1) - yt_x(j, 1) 84 dd_q = yq_w(j, 1) - yq_x(j, 1) 85 dd_u = yu_w(j, 1) - yu_x(j, 1) 86 dd_v = yv_w(j, 1) - yv_x(j, 1) 87 88 ! Merged exchange coefficients 89 dd_Kh = Kech_h_w(j) - Kech_h_x(j) 90 dd_Kq = Kech_q_w(j) - Kech_q_x(j) 91 dd_Km = Kech_m_w(j) - Kech_m_x(j) 92 93 LambdaTs = dd_KTp(j) / Kech_Tp(j) 94 LambdaQs = dd_KQs(j) / Kech_Qs(j) 95 LambdaUs = dd_KUp(j) / Kech_Up(j) 96 LambdaVs = dd_KVp(j) / Kech_Vp(j) 97 98 ! Calcul des coef A, B \'equivalents dans la couche 1 99 100 ! The dTs_forcing and dqs_forcing terms are added for diagnostic purpose ; they should be zero in normal operation. 101 AcoefT(j) = AcoefT_x(j) + sigw(j) * (1. + sigx(j) * LambdaTs) * (dd_AT(j) - C_p(j) * dTs_forcing(j)) 102 AcoefQ(j) = AcoefQ_x(j) + sigw(j) * (1. + sigx(j) * LambdaQs) * (dd_AQ(j) - dqs_forcing(j)) 103 AcoefU(j) = AcoefU_x(j) + sigw(j) * (1. + sigx(j) * LambdaUs) * dd_AU(j) 104 AcoefV(j) = AcoefV_x(j) + sigw(j) * (1. + sigx(j) * LambdaVs) * dd_AV(j) 105 106 107 !! BcoefT(j) = (sigw(j)*Kech_h_w(j)*Kech_T_pw(j)*BcoefT_w(j) + & 108 !! sigx(j)*Kech_h_x(j)*Kech_T_px(j)*BcoefT_x(j) )/(Kech_h(j)*Kech_Tp(j)) 109 !! BcoefQ(j) = (sigw(j)*Kech_q_w(j)*Kech_Q_pw(j)*BcoefQ_w(j) + & 110 !! sigx(j)*Kech_q_x(j)*Kech_Q_px(j)*BcoefQ_x(j) )/(Kech_q(j)*Kech_Qp(j)) 111 !! BcoefU(j) = (sigw(j)*Kech_m_w(j)*Kech_U_pw(j)*BcoefU_w(j) + & 112 !! sigx(j)*Kech_m_x(j)*Kech_U_px(j)*BcoefU_x(j) )/(Kech_m(j)*Kech_Up(j)) 113 !! BcoefV(j) = (sigw(j)*Kech_m_w(j)*Kech_V_pw(j)*BcoefV_w(j) + & 114 !! sigx(j)*Kech_m_x(j)*Kech_V_px(j)*BcoefV_x(j) )/(Kech_m(j)*Kech_Vp(j)) 115 116 !! Print *,'YYYYpbl0: BcoefT_x, sigw, sigx, dd_Kh, dd_KTp, Kech_h_w ', & 117 !! BcoefT_x, sigw, sigx, dd_Kh, dd_KTp, Kech_h_w 118 !! Print *,'YYYYpbl0: Kech_T_pw, dd_BT, Kech_h, Kech_Tp ', & 119 !! Kech_T_pw, dd_BT, Kech_h, Kech_Tp 120 BcoefT(j) = BcoefT_x(j) + sigw(j) * (sigx(j) * dd_Kh * dd_KTp(j) * BcoefT_x(j) + & 121 Kech_h_w(j) * Kech_T_pw(j) * dd_BT(j)) / (Kech_h(j) * Kech_Tp(j)) 122 BcoefQ(j) = BcoefQ_x(j) + sigw(j) * (sigx(j) * dd_Kq * dd_KQs(j) * BcoefQ_x(j) + & 123 Kech_q_w(j) * Kech_Q_sw(j) * dd_BQ(j)) / (Kech_q(j) * Kech_Qs(j)) 124 BcoefU(j) = BcoefU_x(j) + sigw(j) * (sigx(j) * dd_Km * dd_KUp(j) * BcoefU_x(j) + & 125 Kech_m_w(j) * Kech_U_pw(j) * dd_BU(j)) / (Kech_m(j) * Kech_Up(j)) 126 BcoefV(j) = BcoefV_x(j) + sigw(j) * (sigx(j) * dd_Km * dd_KVp(j) * BcoefV_x(j) + & 127 Kech_m_w(j) * Kech_V_pw(j) * dd_BV(j)) / (Kech_m(j) * Kech_Vp(j)) 128 !>jyg 129 130 131 ! Calcul des cdrag \'equivalents dans la couche 132 133 ycdragm(j) = ycdragm_x(j) + sigw(j) * dd_Cdragm(j) 134 ycdragh(j) = ycdragh_x(j) + sigw(j) * dd_Cdragh(j) 135 ycdragq(j) = ycdragq_x(j) + sigw(j) * dd_Cdragq(j) 136 137 ! Calcul de T, q, u et v \'equivalents dans la couche 1 138 !! yt1(j) = yt_x(j,1) + sigw(j)*dd_t*(1.+sigx(j)*dd_Kh/KCT) 139 !! yq1(j) = yq_x(j,1) + sigw(j)*dd_q*(1.+sigx(j)*dd_Kh/KCQ) 140 !! yu1(j) = yu_x(j,1) + sigw(j)*dd_u*(1.+sigx(j)*dd_Km/KCU) 141 !! yv1(j) = yv_x(j,1) + sigw(j)*dd_v*(1.+sigx(j)*dd_Km/KCV) 142 yt1(j) = yt_x(j, 1) + sigw(j) * dd_t 143 yq1(j) = yq_x(j, 1) + sigw(j) * dd_q 144 yu1(j) = yu_x(j, 1) + sigw(j) * dd_u 145 yv1(j) = yv_x(j, 1) + sigw(j) * dd_v 146 147 ENDDO 148 149 END SUBROUTINE wx_pbl0_merge 150 151 SUBROUTINE wx_pbl_dts_merge(knon, dtime, ypplay, ypaprs, & 152 sigw, beta, wcstar, wdens, & 153 AT_x, AT_w, & 154 BT_x, BT_w, & 155 AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, & 156 AcoefT, AcoefQ, BcoefT, BcoefQ, & 157 HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, & 158 phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0, & 159 g_T, g_Q, & 160 Gamma_phiT, Gamma_phiQ, & 161 dTs_ins, dqsatsrf_ins & 162 ) 167 163 168 164 USE wx_pbl_var_mod 169 165 170 USE lmdz_print_control, ONLY: prt_level, lunout171 USE lmdz_ YOETHF166 USE lmdz_print_control, ONLY: prt_level, lunout 167 USE lmdz_yoethf 172 168 USE lmdz_fcttre, ONLY: foeew, foede, qsats, qsatl, dqsats, dqsatl, thermcep 173 174 INCLUDE "YOMCST.h" 175 176 INTEGER, INTENT(IN) :: knon ! number of grid cells 177 REAL, INTENT(IN) :: dtime ! time step size (s) 178 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) 179 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) 180 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pool fractional area 181 REAL, DIMENSION(knon), INTENT(IN) :: beta ! evaporation by potential evaporation 182 REAL, DIMENSION(knon), INTENT(IN) :: wcstar ! cold pool gust front speed 183 REAL, DIMENSION(knon), INTENT(IN) :: wdens ! cold pool number density 184 REAL, DIMENSION(knon), INTENT(IN) :: AT_x, AT_w 185 REAL, DIMENSION(knon), INTENT(IN) :: BT_x, BT_w 186 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0 187 188 REAL, DIMENSION(knon), INTENT(OUT) :: AcoefT, AcoefQ, BcoefT, BcoefQ 189 REAL, DIMENSION(knon), INTENT(OUT) :: HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn 190 REAL, DIMENSION(knon), INTENT(OUT) :: phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0 191 REAL, DIMENSION(knon), INTENT(OUT) :: g_T, g_Q 192 REAL, DIMENSION(knon), INTENT(OUT) :: Gamma_phiT, Gamma_phiQ 193 REAL, DIMENSION(knon), INTENT(OUT) :: dTs_ins, dqsatsrf_ins 194 195 ! Local variables 196 REAL, DIMENSION(knon) :: qsat_x 197 REAL, DIMENSION(knon) :: qsat_w 198 REAL, DIMENSION(knon) :: dqsatdT_x 199 REAL, DIMENSION(knon) :: dqsatdT_w 200 201 REAL, DIMENSION(knon) :: T10_x 202 REAL, DIMENSION(knon) :: T10_w 203 REAL, DIMENSION(knon) :: phiT0_x 204 REAL, DIMENSION(knon) :: phiT0_w 205 REAL, DIMENSION(knon) :: phiQ0_x 206 REAL, DIMENSION(knon) :: phiQ0_w 207 REAL, DIMENSION(knon) :: Rn0_x 208 REAL, DIMENSION(knon) :: Rn0_w 209 REAL, DIMENSION(knon) :: Rp1_x 210 REAL, DIMENSION(knon) :: Rp1_w 211 REAL, DIMENSION(knon) :: Rps_x 212 REAL, DIMENSION(knon) :: Rps_w 213 214 REAL, DIMENSION(knon) :: HTphiT_x 215 REAL, DIMENSION(knon) :: HTphiT_w 216 REAL, DIMENSION(knon) :: HTphiQ_x 217 REAL, DIMENSION(knon) :: HTphiQ_w 218 REAL, DIMENSION(knon) :: HTRn_x 219 REAL, DIMENSION(knon) :: HTRn_w 220 221 REAL, DIMENSION(knon) :: HQphiT_x 222 REAL, DIMENSION(knon) :: HQphiT_w 223 REAL, DIMENSION(knon) :: HQphiQ_x 224 REAL, DIMENSION(knon) :: HQphiQ_w 225 REAL, DIMENSION(knon) :: HQRn_x 226 REAL, DIMENSION(knon) :: HQRn_w 227 228 REAL, DIMENSION(knon) :: HQphiT_b 229 REAL, DIMENSION(knon) :: dd_HQphiT 230 REAL, DIMENSION(knon) :: HQphiQ_b 231 REAL, DIMENSION(knon) :: dd_HQphiQ 232 REAL, DIMENSION(knon) :: HQRn_b 233 REAL, DIMENSION(knon) :: dd_HQRn 234 235 REAL, DIMENSION(knon) :: sigx 236 237 REAL, DIMENSION(knon) :: Ts, T1 238 !!! REAL, DIMENSION(knon) :: qsat, dqsat_dT 239 !!! REAL, DIMENSION(knon) :: phiT0 240 241 !!! REAL, DIMENSION(knon) :: Cp, Lv 242 REAL, DIMENSION(knon) :: tau, Inert 243 244 REAL :: dd_Kh 245 REAL :: zdelta, zcvm5, zcor 246 REAL :: qsat 247 248 INTEGER :: j 249 250 251 !---------------------------------------------------------------------------- 252 ! Reference state 253 ! --------------- 254 ! dqsat_dT_w = dqsat_dT(Ts0_w) dqsat_dT_x = dqsat_dT(Ts0_x) 255 ! T10_w = (AT_w/Cp - Kech_T_w BT_w dtime Ts0_w)/(1 - Kech_T_w BT_w dtime) 256 ! T10_x = (AT_x/Cp - Kech_T_x BT_x dtime Ts0_x)/(1 - Kech_T_x BT_x dtime) 257 ! phiT0_w = Kech_T_pw (AT_w - Cp Ts0_w) phiT0_x = Kech_T_px (AT_x - Cp Ts0_x) 258 ! phiQ0_w = Kech_Q_sw (beta AQ_w - qsatsrf0_w) phiQ0_x = Kech_Q_sx (beta AQ_x - qsatsrf0_x) 259 ! Rn0_w = eps_1 Rsigma T10_w^4 - Rsigma Ts0_w^4 Rn0_x = eps_1 Rsigma T10_x^4 - Rsigma Ts0_x^4 260 ! Rp1_w = 4 eps_1 Rsigma T10_w^3 Rp1_x = 4 eps_1 Rsigma T10_x^3 261 ! Rps_w = 4 Rsigma Ts0_w^3 Rps_x = 4 Rsigma Ts0_x^3 262 263 ! phiT0_b = sigw phiT0_w + sigx phiT0_x 264 ! dphiT0 = phiT0_w - phiT0_x 265 ! phiQ0_b = sigw phiQ0_w + sigx phiQ0_x 266 ! dphiQ0 = phiQ0_w - phiQ0_x 267 ! Rn0_b = sigw Rn0_w + sigx Rn0_x 169 USE lmdz_yomcst 170 171 IMPLICIT NONE 172 173 INTEGER, INTENT(IN) :: knon ! number of grid cells 174 REAL, INTENT(IN) :: dtime ! time step size (s) 175 REAL, DIMENSION(knon, klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) 176 REAL, DIMENSION(knon, klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) 177 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pool fractional area 178 REAL, DIMENSION(knon), INTENT(IN) :: beta ! evaporation by potential evaporation 179 REAL, DIMENSION(knon), INTENT(IN) :: wcstar ! cold pool gust front speed 180 REAL, DIMENSION(knon), INTENT(IN) :: wdens ! cold pool number density 181 REAL, DIMENSION(knon), INTENT(IN) :: AT_x, AT_w 182 REAL, DIMENSION(knon), INTENT(IN) :: BT_x, BT_w 183 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0 184 185 REAL, DIMENSION(knon), INTENT(OUT) :: AcoefT, AcoefQ, BcoefT, BcoefQ 186 REAL, DIMENSION(knon), INTENT(OUT) :: HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn 187 REAL, DIMENSION(knon), INTENT(OUT) :: phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0 188 REAL, DIMENSION(knon), INTENT(OUT) :: g_T, g_Q 189 REAL, DIMENSION(knon), INTENT(OUT) :: Gamma_phiT, Gamma_phiQ 190 REAL, DIMENSION(knon), INTENT(OUT) :: dTs_ins, dqsatsrf_ins 191 192 ! Local variables 193 REAL, DIMENSION(knon) :: qsat_x 194 REAL, DIMENSION(knon) :: qsat_w 195 REAL, DIMENSION(knon) :: dqsatdT_x 196 REAL, DIMENSION(knon) :: dqsatdT_w 197 198 REAL, DIMENSION(knon) :: T10_x 199 REAL, DIMENSION(knon) :: T10_w 200 REAL, DIMENSION(knon) :: phiT0_x 201 REAL, DIMENSION(knon) :: phiT0_w 202 REAL, DIMENSION(knon) :: phiQ0_x 203 REAL, DIMENSION(knon) :: phiQ0_w 204 REAL, DIMENSION(knon) :: Rn0_x 205 REAL, DIMENSION(knon) :: Rn0_w 206 REAL, DIMENSION(knon) :: Rp1_x 207 REAL, DIMENSION(knon) :: Rp1_w 208 REAL, DIMENSION(knon) :: Rps_x 209 REAL, DIMENSION(knon) :: Rps_w 210 211 REAL, DIMENSION(knon) :: HTphiT_x 212 REAL, DIMENSION(knon) :: HTphiT_w 213 REAL, DIMENSION(knon) :: HTphiQ_x 214 REAL, DIMENSION(knon) :: HTphiQ_w 215 REAL, DIMENSION(knon) :: HTRn_x 216 REAL, DIMENSION(knon) :: HTRn_w 217 218 REAL, DIMENSION(knon) :: HQphiT_x 219 REAL, DIMENSION(knon) :: HQphiT_w 220 REAL, DIMENSION(knon) :: HQphiQ_x 221 REAL, DIMENSION(knon) :: HQphiQ_w 222 REAL, DIMENSION(knon) :: HQRn_x 223 REAL, DIMENSION(knon) :: HQRn_w 224 225 REAL, DIMENSION(knon) :: HQphiT_b 226 REAL, DIMENSION(knon) :: dd_HQphiT 227 REAL, DIMENSION(knon) :: HQphiQ_b 228 REAL, DIMENSION(knon) :: dd_HQphiQ 229 REAL, DIMENSION(knon) :: HQRn_b 230 REAL, DIMENSION(knon) :: dd_HQRn 231 232 REAL, DIMENSION(knon) :: sigx 233 234 REAL, DIMENSION(knon) :: Ts, T1 235 !!! REAL, DIMENSION(knon) :: qsat, dqsat_dT 236 !!! REAL, DIMENSION(knon) :: phiT0 237 238 !!! REAL, DIMENSION(knon) :: Cp, Lv 239 REAL, DIMENSION(knon) :: tau, Inert 240 241 REAL :: dd_Kh 242 REAL :: zdelta, zcvm5, zcor 243 REAL :: qsat 244 245 INTEGER :: j 246 247 248 !---------------------------------------------------------------------------- 249 ! Reference state 250 ! --------------- 251 ! dqsat_dT_w = dqsat_dT(Ts0_w) dqsat_dT_x = dqsat_dT(Ts0_x) 252 ! T10_w = (AT_w/Cp - Kech_T_w BT_w dtime Ts0_w)/(1 - Kech_T_w BT_w dtime) 253 ! T10_x = (AT_x/Cp - Kech_T_x BT_x dtime Ts0_x)/(1 - Kech_T_x BT_x dtime) 254 ! phiT0_w = Kech_T_pw (AT_w - Cp Ts0_w) phiT0_x = Kech_T_px (AT_x - Cp Ts0_x) 255 ! phiQ0_w = Kech_Q_sw (beta AQ_w - qsatsrf0_w) phiQ0_x = Kech_Q_sx (beta AQ_x - qsatsrf0_x) 256 ! Rn0_w = eps_1 Rsigma T10_w^4 - Rsigma Ts0_w^4 Rn0_x = eps_1 Rsigma T10_x^4 - Rsigma Ts0_x^4 257 ! Rp1_w = 4 eps_1 Rsigma T10_w^3 Rp1_x = 4 eps_1 Rsigma T10_x^3 258 ! Rps_w = 4 Rsigma Ts0_w^3 Rps_x = 4 Rsigma Ts0_x^3 259 260 ! phiT0_b = sigw phiT0_w + sigx phiT0_x 261 ! dphiT0 = phiT0_w - phiT0_x 262 ! phiQ0_b = sigw phiQ0_w + sigx phiQ0_x 263 ! dphiQ0 = phiQ0_w - phiQ0_x 264 ! Rn0_b = sigw Rn0_w + sigx Rn0_x 268 265 dRn0 = Rn0_w - Rn0_x 269 266 270 267 271 !---------------------------------------------------------------------------- 272 ! Elementary enthalpy equations 273 ! ----------------------------- 274 ! phiT_w = phiT0_w - HTphiT_w (Ts_w-Ts0_w) phiT_x = phiT0_x - HTphiT_x (Ts_x-Ts0_x) 275 ! phiQ_w = phiQ0_w - HTphiQ_w (Ts_w-Ts0_w) phiQ_x = phiQ0_x - HTphiQ_x (Ts_x-Ts0_x) 276 ! Rn_w = Rn0_w - HTRn_w (Ts_w-Ts0_w) Rn_x = Rn0_x - HTRn_x (Ts_x-Ts0_x) 277 ! DFlux_DT coefficients 278 ! --------------------- 279 ! Heat flux equation 280 ! HTphiT_w = Cp Kech_T_pw HTphiT_x = Cp Kech_T_px 281 ! Moisture flux equation 282 ! HTphiQ_w = beta Kech_Q_sw dqsat_dT_w HTphiQ_x = beta Kech_Q_sx dqsat_dT_x 283 ! Radiation equation 284 ! HTRn_w = Rp1_w Kech_T_pw BcoefT_w dtime + Rps_w HTRn_x = Rp1_x Kech_T_px BcoefT_x dtime + Rps_x 285 286 !---------------------------------------------------------------------------- 287 ! Elementary moisture equations 288 ! ----------------------------- 289 ! beta Ts_w = beta Ts0_w + QQ_w (qsatsrf_w-qsatsrf0_w) beta Ts_x = beta Ts0_x + QQ_x (qsatsrf_x-qsatsrf0_x) 290 ! beta phiT_w = beta phiT0_w - HQphiT_w (qsatsrf_w-qsatsrf0_w) beta phiQ_x = beta phiQ0_x - HTphiQ_x (qsatsrf_x-qsatsrf0_x) 291 ! beta phiQ_w = beta phiQ0_w - HQphiQ_w (qsatsrf_w-qsatsrf0_w) beta phiQ_x = beta phiQ0_x - HTphiQ_x (qsatsrf_x-qsatsrf0_x) 292 ! beta Rn_w = beta Rn0_w - HQRn_w (qsatsrf_w-qsatsrf0_w) beta Rn_x = beta Rn0_x - HTRn_x (qsatsrf_x-qsatsrf0_x) 293 ! DFluxDQ coefficients 294 ! --------------------- 295 ! dqsat_dT equation 296 ! QQ_w = 1. / dqsat_dT_w QQ_x = 1. / dqsat_dT_x 297 ! Heat flux equation 298 ! HQphiT_w = Cp Kech_T_pw QQ_w HQphiT_x = Cp Kech_T_px QQ_x 299 ! Moisture flux equation 300 ! HQphiQ_w = beta Kech_Q_sw HQphiQ_x = beta Kech_Q_sx 301 ! Radiation equation 302 ! HQRn_w = (Rp1_w Kech_T_pw BcoefT_w dtime + Rps_w) QQ_w 303 ! HQRn_x = (Rp1_x Kech_T_px BcoefT_x dtime + Rps_x) QQ_x 304 305 !---------------------------------------------------------------------------- 306 ! Mean values and w-x differences 307 ! ------------------------------- 308 ! HTphiT_b = sigw HTphiT_w + sigx HTphiT_x dd_HTphiT = HTphiT_w - HTphiT_x 309 ! HTphiQ_b = sigw HTphiQ_w + sigx HTphiQ_x dd_HTphiQ = HTphiQ_w - HTphiQ_x 310 ! HTRn_b = sigw HTRn_w + sigx HTRn_x dd_HTRn = HTRn_w - HTRn_x 311 312 ! QQ_b = sigw QQ_w + sigx QQ_x dd_QQ = QQ_w - QQ_x 313 ! HQphiT_b = sigw HQphiT_w + sigx HQphiT_x dd_HQphiT = HQphiT_w - HQphiT_x 314 ! HQphiQ_b = sigw HQphiQ_w + sigx HQphiQ_x dd_HQphiQ = HQphiQ_w - HQphiQ_x 315 ! HQRn_b = sigw HQRn_w + sigx HQRn_x dd_HQRn = HQRn_w - HQRn_x 316 317 !---------------------------------------------------------------------------- 318 ! Equations 319 ! --------- 320 ! (1 - g_T) dTs = dTs_ins + Gamma_phiT phiT 321 ! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ 322 323 ! Feedback Gains 324 ! -------------- 325 ! g_T = - (sqrt(tau)/I) [ HTphiT_b + Lv HTphiQ_b + HTRn_b + & 326 ! (dd_HTphiT + Lv dd_HTphiQ + dd_HTRn) (sigx - sigw - sigw sigx dd_HTphiT/HTphiT_b) ] 327 ! g_Q = - (sqrt(tau)/(I QQ_b)) ( HQphiT_b + Lv HQphiQ_b + HQRn_b ) - & 328 ! (sigx - sigw - sigw sigx dd_HQphiQ/HQphiQ_b) & 329 ! [ dd_QQ/QQ_b + (sqrt(tau)/(I QQ_b))(dd_HQphiT + Lv dd_HQphiQ + dd_HQRn) ] 330 331 ! Ts, qs Coupling coefficients / 332 ! ---------------------------- 333 ! Gamma_phiT = (sqrt(tau)/(I HTphiT_b)) (dd_HTphiT + Lv dd_HTphiQ + dd_HTRn) 334 ! Gamma_phiQ = (1/(HQphiQ_b QQ_b)) [ dd_QQ + (sqrt(tau)/(I )) (dd_HQphiT + Lv dd_HQphiQ + dd_HQRn) ] 335 336 ! Insensitive changes 337 ! ------------------- 338 ! dTs_ins = (1 - g_T) dTs0 - Gamma_phiT phiT0_b 339 ! dqsatsrf_ins = (1 - g_Q) dqsatsrf0 - Gamma_phiQ phiQ0_b 340 341 !---------------------------------------------------------------------------- 342 ! Effective coefficients Acoef and Bcoef 343 ! -------------------------------------- 344 ! Equations 345 ! --------- 346 ! Cp Ta = AcoefT + BcoefT phiT dtime 347 ! qa = AcoefQ + BcoefQ phiQ dtime 348 ! Coefficients 349 ! ------------ 350 ! AcoefT = AcoefT0 - sigw sigx (dd_KTp/Kech_Tp) Cp dTs_ins/(1 - g_T) 351 ! BcoefT = BcoefT0 - sigw sigx (dd_KTp/Kech_Tp) Cp Gamma_phiT/(1 - g_T)/dtime 352 353 ! AcoefQ = AcoefQ0 - sigw sigx (dd_KQp/Kech_Qp) dqs_ins/(1 - g_Q) 354 ! BcoefQ = BcoefQ0 - sigw sigx (dd_KQp/Kech_Qp) Gamma_phiq/(1 - g_Q)/dtime 355 356 !============================================================================== 357 358 359 ! Parameters 360 ! ---------- 361 Inert(1:knon) = 2000. 362 tau(1:knon) = sqrt(sigw(1:knon)/max(rpi*wdens(1:knon)*wcstar(1:knon)**2 , & 363 sigw(1:knon)*1.e-12,smallestreal)) 364 sigx(1:knon) = 1.-sigw(1:knon) 365 !! Compute Cp, Lv, qsat, dqsat_dT. 366 ! C_p(1:knon) = RCpd 367 ! L_v(1:knon) = RLvtt 368 369 ! print *,' AAAA wx_pbl_dTs, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:) 370 371 372 T10_x(1:knon) = (AT_x(1:knon)/C_p(1:knon) - Kech_h_x(1:knon)*BT_x(1:knon)*dtime*Ts0_x(1:knon))/ & 373 (1 - Kech_h_x(1:knon)*BT_x(1:knon)*dtime) 374 T10_w(1:knon) = (AT_w(1:knon)/C_p(1:knon) - Kech_h_w(1:knon)*BT_w(1:knon)*dtime*Ts0_w(1:knon))/ & 375 (1 - Kech_h_w(1:knon)*BT_w(1:knon)*dtime) 376 377 phiT0_x(1:knon) = Kech_T_px(1:knon)*(AT_x(1:knon) - C_p(1:knon)*Ts0_x(1:knon)) 378 phiT0_w(1:knon) = Kech_T_pw(1:knon)*(AT_w(1:knon) - C_p(1:knon)*Ts0_w(1:knon)) 379 380 phiQ0_x(1:knon) = Kech_Q_sx(1:knon)*(beta(1:knon)*AQ_x(1:knon) - qsatsrf0_x(1:knon)) 381 phiQ0_w(1:knon) = Kech_Q_sw(1:knon)*(beta(1:knon)*AQ_w(1:knon) - qsatsrf0_w(1:knon)) 382 383 Rn0_x(1:knon) = eps_1*Rsigma*T10_x(1:knon)**4 - Rsigma*Ts0_x(1:knon)**4 384 Rn0_w(1:knon) = eps_1*Rsigma*T10_w(1:knon)**4 - Rsigma*Ts0_w(1:knon)**4 385 386 Rp1_x(1:knon) = 4*eps_1*Rsigma*T10_x(1:knon)**3 387 Rp1_w(1:knon) = 4*eps_1*Rsigma*T10_w(1:knon)**3 388 389 Rps_x(1:knon) = 4*Rsigma*Ts0_x(1:knon)**3 390 Rps_w(1:knon) = 4*Rsigma*Ts0_w(1:knon)**3 391 392 ! DFlux_DT coefficients 393 ! --------------------- 394 ! Heat flux equation 395 HTphiT_x(1:knon) = C_p(1:knon)*Kech_T_px(1:knon) 396 HTphiT_w(1:knon) = C_p(1:knon)*Kech_T_pw(1:knon) 397 ! Moisture flux equation 398 HTphiQ_x(1:knon) = beta(1:knon)*Kech_Q_sx(1:knon)*dqsatdT0_x(1:knon) 399 HTphiQ_w(1:knon) = beta(1:knon)*Kech_Q_sw(1:knon)*dqsatdT0_w(1:knon) 400 ! Radiation equation 401 HTRn_x(1:knon) = Rp1_x(1:knon)*Kech_T_px(1:knon)*BT_x(1:knon)*dtime + Rps_x(1:knon) 402 HTRn_w(1:knon) = Rp1_w(1:knon)*Kech_T_pw(1:knon)*BT_w(1:knon)*dtime + Rps_w(1:knon) 403 404 ! DFluxDQ coefficients 405 ! --------------------- 406 ! Heat flux equation 407 HQphiT_x(1:knon) = C_p(1:knon)*Kech_T_px(1:knon)*QQ_x(1:knon) 408 HQphiT_w(1:knon) = C_p(1:knon)*Kech_T_pw(1:knon)*QQ_w(1:knon) 409 ! Moisture flux equation 410 HQphiQ_x(1:knon) = beta(1:knon)*Kech_Q_sx(1:knon) 411 HQphiQ_w(1:knon) = beta(1:knon)*Kech_Q_sw(1:knon) 412 ! Radiation equation 413 HQRn_x(1:knon) = (Rp1_x(1:knon)*Kech_T_px(1:knon)*BT_x(1:knon)*dtime + Rps_x(1:knon))*QQ_x(1:knon) 414 HQRn_w(1:knon) = (Rp1_w(1:knon)*Kech_T_pw(1:knon)*BT_w(1:knon)*dtime + Rps_w(1:knon))*QQ_w(1:knon) 415 416 ! Mean values and w-x differences 417 ! ------------------------------- 418 phiT0_b(1:knon) = sigw(1:knon)*phiT0_w(1:knon) + sigx(1:knon)*phiT0_x(1:knon) 419 phiQ0_b(1:knon) = sigw(1:knon)*phiQ0_w(1:knon) + sigx(1:knon)*phiQ0_x(1:knon) 420 Rn0_b(1:knon) = sigw(1:knon)*Rn0_w(1:knon) + sigx(1:knon)*Rn0_x(1:knon) 421 422 dphiT0(1:knon) = phiT0_w(1:knon) - phiT0_x(1:knon) 423 dphiQ0(1:knon) = phiQ0_w(1:knon) - phiQ0_x(1:knon) 424 dRn0(1:knon) = Rn0_w(1:knon) - Rn0_x(1:knon) 425 426 HTphiT_b(1:knon) = sigw(1:knon)*HTphiT_w(1:knon) + sigx(1:knon)*HTphiT_x(1:knon) 427 dd_HTphiT(1:knon) = HTphiT_w(1:knon) - HTphiT_x(1:knon) 428 429 HTphiQ_b(1:knon) = sigw(1:knon)*HTphiQ_w(1:knon) + sigx(1:knon)*HTphiQ_x(1:knon) 430 dd_HTphiQ(1:knon) = HTphiQ_w(1:knon) - HTphiQ_x(1:knon) 431 432 HTRn_b(1:knon) = sigw(1:knon)*HTRn_w(1:knon) + sigx(1:knon)*HTRn_x(1:knon) 433 dd_HTRn(1:knon) = HTRn_w(1:knon) - HTRn_x(1:knon) 434 435 HQphiT_b(1:knon) = sigw(1:knon)*HQphiT_w(1:knon) + sigx(1:knon)*HQphiT_x(1:knon) 436 dd_HQphiT(1:knon) = HQphiT_w(1:knon) - HQphiT_x(1:knon) 437 438 HQphiQ_b(1:knon) = sigw(1:knon)*HQphiQ_w(1:knon) + sigx(1:knon)*HQphiQ_x(1:knon) 439 dd_HQphiQ(1:knon) = HQphiQ_w - HQphiQ_x(1:knon) 440 441 HQRn_b(1:knon) = sigw(1:knon)*HQRn_w(1:knon) + sigx(1:knon)*HQRn_x(1:knon) 442 dd_HQRn(1:knon) = HQRn_w(1:knon) - HQRn_x(1:knon) 443 444 ! Feedback Gains 445 ! -------------- 446 g_T(1:knon) = - (sqrt(tau(1:knon))/Inert(1:knon)) & 447 * (HTphiT_b(1:knon) + L_v(1:knon)*HTphiQ_b(1:knon) + HTRn_b(1:knon) & 448 + (dd_HTphiT(1:knon) + L_v(1:knon)*dd_HTphiQ(1:knon) + dd_HTRn(1:knon)) & 449 * (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HTphiT(1:knon)/HTphiT_b(1:knon)) ) 450 451 !!!! DO j = 1,knon 452 !!!! IF (mod(j,20) .EQ.0) THEN 453 !!!! print *, ' j dd_QQ QQ_b dd_HQphiQ dd_HQphiT dd_HQRn HQphiQ_b HQphiT_b HQRn_b ' 454 !!!! ENDIF 455 !!!! print 1789, j, dd_QQ(j), QQ_b(j), dd_HQphiQ(j), dd_HQphiT(j), dd_HQRn(j), HQphiQ_b(j), HQphiT_b(j), HQRn_b(j) 456 !!!! 1789 FORMAT( I4, 10(1X,E10.2)) 457 !!!! ENDDO 458 g_Q(1:knon) = - (dd_QQ(1:knon)/QQ_b(1:knon)) * & 459 (sigx(1:knon)-sigw(1:knon)-sigw(1:knon)*sigx(1:knon)*dd_KQs(1:knon)/Kech_Qs(1:knon)) & 460 - sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon)) * & 461 ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) + & 462 (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_KQs(1:knon)/Kech_Qs(1:knon)) * & 463 (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) 464 465 !! g_Q(1:knon) = - (dd_QQ(1:knon)/QQ_b(1:knon)) * & 466 !! (sigx(1:knon)-sigw(1:knon)-sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) & 467 !! - sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon)) * & 468 !! ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) + & 469 !! (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) * & 470 !! (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) 471 472 !! g_Q(1:knon) = - (sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon))) * & 473 !! ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) ) & 474 !! - (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) * & 475 !! ( dd_QQ(1:knon)/QQ_b(1:knon) & 476 !! + (sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon))) & 477 !! * (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) 478 479 ! Ts, qs Coupling coefficients / 480 ! ---------------------------- 481 Gamma_phiT(1:knon) = (sqrt(tau(1:knon))/(Inert(1:knon)*HTphiT_b(1:knon))) & 482 * (dd_HTphiT(1:knon) + L_v(1:knon)*dd_HTphiQ(1:knon) + dd_HTRn(1:knon)) 483 484 Gamma_phiQ(1:knon) = (1./(Kech_Qs(1:knon)*QQ_b(1:knon))) * & 485 ( dd_QQ(1:knon) & 486 + (sqrt(tau(1:knon))/(Inert(1:knon))) * & 487 (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) 488 489 !! Gamma_phiQ(1:knon) = (beta(1:knon)/(HQphiQ_b(1:knon)*QQ_b(1:knon))) * & 490 !! ( dd_QQ(1:knon) & 491 !! + (sqrt(tau(1:knon))/(Inert(1:knon))) * & 492 !! (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) 493 494 !! Gamma_phiQ(1:knon) = (1/(HQphiQ_b(1:knon)*QQ_b(1:knon))) & 495 !! * ( dd_QQ(1:knon) & 496 !! + (sqrt(tau(1:knon))/(Inert(1:knon))) & 497 !! * (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) 498 499 ! Insensitive changes 500 ! ------------------- 501 dTs_ins(1:knon) = (sqrt(tau(1:knon))/Inert(1:knon))* & 502 (dphiT0(1:knon) + L_v(1:knon)*dphiQ0(1:knon) + dRn0(1:knon)) 503 504 dqsatsrf_ins(1:knon) = (beta(1:knon)/QQ_b(1:knon))*dTs_ins(1:knon) 505 506 IF (prt_level >= 10) THEN 507 print *,'wx_pbl_merge, tau ', tau 508 print *,'wx_pbl_merge, AcoefT0 ', AcoefT0 509 print *,'wx_pbl_merge, AcoefQ0 ', AcoefQ0 510 print *,'wx_pbl_merge, BcoefT0 ', BcoefT0 511 print *,'wx_pbl_merge, BcoefQ0 ', BcoefQ0 512 print *,'wx_pbl_merge, qsat0_w, qsat0_x ', (qsat0_w(j), qsat0_x(j),j=1,knon) 513 print *,'wx_pbl_merge, dqsatdT0_w, dqsatdT0_x ', (dqsatdT0_w(j), dqsatdT0_x(j),j=1,knon) 514 ENDIF 515 516 !---------------------------------------------------------------------------- 517 518 !------------------------------------------------------------------------------ 519 520 ! Effective coefficients Acoef and Bcoef 521 ! -------------------------------------- 522 DO j = 1,knon 523 AcoefT(j) = AcoefT0(j) - sigw(j)*sigx(j)*(dd_KTp(j)/Kech_Tp(j))*C_p(j)* & 524 (dTs0(j) + (dTs_ins(j)-dTs0(j)-Gamma_phiT(j)*phiT0_b(j))/(1. - g_T(j))) 525 BcoefT(j) = BcoefT0(j) - sigw(j)*sigx(j)*(dd_KTp(j)/Kech_Tp(j))*C_p(j)*Gamma_phiT(j)/(1. - g_T(j))/dtime 526 527 AcoefQ(j) = AcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))* & 528 (dqsatsrf0(j) + (dqsatsrf_ins(j)-(beta(j)/QQ_b(j))*dTs0(j)-Gamma_phiQ(j)*phiQ0_b(j))/(1 - g_Q(j)))/ & 529 max(beta(j),1.e-4) 530 BcoefQ(j) = BcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))*Gamma_phiQ(j)/(1 - g_Q(j))/ & 531 (max(beta(j),1.e-4)*dtime) 532 !! AcoefQ(j) = AcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))* & 533 !! (dqsatsrf0(j) + (dqsatsrf_ins(j)-(beta(j)/QQ_b(j))*dTs0(j)-Gamma_phiQ(j)*phiQ0_b(j))/(1 - g_Q(j)))/ & 534 !! beta(j) 535 !! BcoefQ(j) = BcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))*Gamma_phiQ(j)/(1 - g_Q(j))/(beta(j)*dtime) 536 ENDDO ! j = 1,knon 537 538 IF (prt_level >= 10) THEN 539 print *,'wx_pbl_dts AAAA BcoefQ, BcoefQ0, sigw ', & 540 BcoefQ, BcoefQ0, sigw 541 print *,'wx_pbl_dts_merge, dTs_ins ', dTs_ins 542 print *,'wx_pbl_dts_merge, dqs_ins ', dqsatsrf_ins 543 ENDIF 544 545 546 547 END SUBROUTINE wx_pbl_dts_merge 548 549 SUBROUTINE wx_pbl_split(knon, nsrf, dtime, sigw, beta, iflag_split, & 550 g_T, g_Q, & 551 Gamma_phiT, Gamma_phiQ, & 552 dTs_ins, dqsatsrf_ins, & 553 phiT, phiQ, phiU, phiV, & 554 !!!! HTRn_b, dd_HTRn, HTphiT_b, dd_HTphiT, & 555 phiQ0_b, phiT0_b, & 556 phiT_x, phiT_w, & 557 phiQ_x, phiQ_w, & 558 phiU_x, phiU_w, & 559 phiV_x, phiV_w, & 560 philat_x, philat_w, & 561 !!!! Rn_b, dRn, & 562 dqsatsrf, & 563 dTs, delta_qsurf & 564 ) 268 !---------------------------------------------------------------------------- 269 ! Elementary enthalpy equations 270 ! ----------------------------- 271 ! phiT_w = phiT0_w - HTphiT_w (Ts_w-Ts0_w) phiT_x = phiT0_x - HTphiT_x (Ts_x-Ts0_x) 272 ! phiQ_w = phiQ0_w - HTphiQ_w (Ts_w-Ts0_w) phiQ_x = phiQ0_x - HTphiQ_x (Ts_x-Ts0_x) 273 ! Rn_w = Rn0_w - HTRn_w (Ts_w-Ts0_w) Rn_x = Rn0_x - HTRn_x (Ts_x-Ts0_x) 274 ! DFlux_DT coefficients 275 ! --------------------- 276 ! Heat flux equation 277 ! HTphiT_w = Cp Kech_T_pw HTphiT_x = Cp Kech_T_px 278 ! Moisture flux equation 279 ! HTphiQ_w = beta Kech_Q_sw dqsat_dT_w HTphiQ_x = beta Kech_Q_sx dqsat_dT_x 280 ! Radiation equation 281 ! HTRn_w = Rp1_w Kech_T_pw BcoefT_w dtime + Rps_w HTRn_x = Rp1_x Kech_T_px BcoefT_x dtime + Rps_x 282 283 !---------------------------------------------------------------------------- 284 ! Elementary moisture equations 285 ! ----------------------------- 286 ! beta Ts_w = beta Ts0_w + QQ_w (qsatsrf_w-qsatsrf0_w) beta Ts_x = beta Ts0_x + QQ_x (qsatsrf_x-qsatsrf0_x) 287 ! beta phiT_w = beta phiT0_w - HQphiT_w (qsatsrf_w-qsatsrf0_w) beta phiQ_x = beta phiQ0_x - HTphiQ_x (qsatsrf_x-qsatsrf0_x) 288 ! beta phiQ_w = beta phiQ0_w - HQphiQ_w (qsatsrf_w-qsatsrf0_w) beta phiQ_x = beta phiQ0_x - HTphiQ_x (qsatsrf_x-qsatsrf0_x) 289 ! beta Rn_w = beta Rn0_w - HQRn_w (qsatsrf_w-qsatsrf0_w) beta Rn_x = beta Rn0_x - HTRn_x (qsatsrf_x-qsatsrf0_x) 290 ! DFluxDQ coefficients 291 ! --------------------- 292 ! dqsat_dT equation 293 ! QQ_w = 1. / dqsat_dT_w QQ_x = 1. / dqsat_dT_x 294 ! Heat flux equation 295 ! HQphiT_w = Cp Kech_T_pw QQ_w HQphiT_x = Cp Kech_T_px QQ_x 296 ! Moisture flux equation 297 ! HQphiQ_w = beta Kech_Q_sw HQphiQ_x = beta Kech_Q_sx 298 ! Radiation equation 299 ! HQRn_w = (Rp1_w Kech_T_pw BcoefT_w dtime + Rps_w) QQ_w 300 ! HQRn_x = (Rp1_x Kech_T_px BcoefT_x dtime + Rps_x) QQ_x 301 302 !---------------------------------------------------------------------------- 303 ! Mean values and w-x differences 304 ! ------------------------------- 305 ! HTphiT_b = sigw HTphiT_w + sigx HTphiT_x dd_HTphiT = HTphiT_w - HTphiT_x 306 ! HTphiQ_b = sigw HTphiQ_w + sigx HTphiQ_x dd_HTphiQ = HTphiQ_w - HTphiQ_x 307 ! HTRn_b = sigw HTRn_w + sigx HTRn_x dd_HTRn = HTRn_w - HTRn_x 308 309 ! QQ_b = sigw QQ_w + sigx QQ_x dd_QQ = QQ_w - QQ_x 310 ! HQphiT_b = sigw HQphiT_w + sigx HQphiT_x dd_HQphiT = HQphiT_w - HQphiT_x 311 ! HQphiQ_b = sigw HQphiQ_w + sigx HQphiQ_x dd_HQphiQ = HQphiQ_w - HQphiQ_x 312 ! HQRn_b = sigw HQRn_w + sigx HQRn_x dd_HQRn = HQRn_w - HQRn_x 313 314 !---------------------------------------------------------------------------- 315 ! Equations 316 ! --------- 317 ! (1 - g_T) dTs = dTs_ins + Gamma_phiT phiT 318 ! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ 319 320 ! Feedback Gains 321 ! -------------- 322 ! g_T = - (sqrt(tau)/I) [ HTphiT_b + Lv HTphiQ_b + HTRn_b + & 323 ! (dd_HTphiT + Lv dd_HTphiQ + dd_HTRn) (sigx - sigw - sigw sigx dd_HTphiT/HTphiT_b) ] 324 ! g_Q = - (sqrt(tau)/(I QQ_b)) ( HQphiT_b + Lv HQphiQ_b + HQRn_b ) - & 325 ! (sigx - sigw - sigw sigx dd_HQphiQ/HQphiQ_b) & 326 ! [ dd_QQ/QQ_b + (sqrt(tau)/(I QQ_b))(dd_HQphiT + Lv dd_HQphiQ + dd_HQRn) ] 327 328 ! Ts, qs Coupling coefficients / 329 ! ---------------------------- 330 ! Gamma_phiT = (sqrt(tau)/(I HTphiT_b)) (dd_HTphiT + Lv dd_HTphiQ + dd_HTRn) 331 ! Gamma_phiQ = (1/(HQphiQ_b QQ_b)) [ dd_QQ + (sqrt(tau)/(I )) (dd_HQphiT + Lv dd_HQphiQ + dd_HQRn) ] 332 333 ! Insensitive changes 334 ! ------------------- 335 ! dTs_ins = (1 - g_T) dTs0 - Gamma_phiT phiT0_b 336 ! dqsatsrf_ins = (1 - g_Q) dqsatsrf0 - Gamma_phiQ phiQ0_b 337 338 !---------------------------------------------------------------------------- 339 ! Effective coefficients Acoef and Bcoef 340 ! -------------------------------------- 341 ! Equations 342 ! --------- 343 ! Cp Ta = AcoefT + BcoefT phiT dtime 344 ! qa = AcoefQ + BcoefQ phiQ dtime 345 ! Coefficients 346 ! ------------ 347 ! AcoefT = AcoefT0 - sigw sigx (dd_KTp/Kech_Tp) Cp dTs_ins/(1 - g_T) 348 ! BcoefT = BcoefT0 - sigw sigx (dd_KTp/Kech_Tp) Cp Gamma_phiT/(1 - g_T)/dtime 349 350 ! AcoefQ = AcoefQ0 - sigw sigx (dd_KQp/Kech_Qp) dqs_ins/(1 - g_Q) 351 ! BcoefQ = BcoefQ0 - sigw sigx (dd_KQp/Kech_Qp) Gamma_phiq/(1 - g_Q)/dtime 352 353 !============================================================================== 354 355 356 ! Parameters 357 ! ---------- 358 Inert(1:knon) = 2000. 359 tau(1:knon) = sqrt(sigw(1:knon) / max(rpi * wdens(1:knon) * wcstar(1:knon)**2, & 360 sigw(1:knon) * 1.e-12, smallestreal)) 361 sigx(1:knon) = 1. - sigw(1:knon) 362 !! Compute Cp, Lv, qsat, dqsat_dT. 363 ! C_p(1:knon) = RCpd 364 ! L_v(1:knon) = RLvtt 365 366 ! print *,' AAAA wx_pbl_dTs, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:) 367 368 T10_x(1:knon) = (AT_x(1:knon) / C_p(1:knon) - Kech_h_x(1:knon) * BT_x(1:knon) * dtime * Ts0_x(1:knon)) / & 369 (1 - Kech_h_x(1:knon) * BT_x(1:knon) * dtime) 370 T10_w(1:knon) = (AT_w(1:knon) / C_p(1:knon) - Kech_h_w(1:knon) * BT_w(1:knon) * dtime * Ts0_w(1:knon)) / & 371 (1 - Kech_h_w(1:knon) * BT_w(1:knon) * dtime) 372 373 phiT0_x(1:knon) = Kech_T_px(1:knon) * (AT_x(1:knon) - C_p(1:knon) * Ts0_x(1:knon)) 374 phiT0_w(1:knon) = Kech_T_pw(1:knon) * (AT_w(1:knon) - C_p(1:knon) * Ts0_w(1:knon)) 375 376 phiQ0_x(1:knon) = Kech_Q_sx(1:knon) * (beta(1:knon) * AQ_x(1:knon) - qsatsrf0_x(1:knon)) 377 phiQ0_w(1:knon) = Kech_Q_sw(1:knon) * (beta(1:knon) * AQ_w(1:knon) - qsatsrf0_w(1:knon)) 378 379 Rn0_x(1:knon) = eps_1 * Rsigma * T10_x(1:knon)**4 - Rsigma * Ts0_x(1:knon)**4 380 Rn0_w(1:knon) = eps_1 * Rsigma * T10_w(1:knon)**4 - Rsigma * Ts0_w(1:knon)**4 381 382 Rp1_x(1:knon) = 4 * eps_1 * Rsigma * T10_x(1:knon)**3 383 Rp1_w(1:knon) = 4 * eps_1 * Rsigma * T10_w(1:knon)**3 384 385 Rps_x(1:knon) = 4 * Rsigma * Ts0_x(1:knon)**3 386 Rps_w(1:knon) = 4 * Rsigma * Ts0_w(1:knon)**3 387 388 ! DFlux_DT coefficients 389 ! --------------------- 390 ! Heat flux equation 391 HTphiT_x(1:knon) = C_p(1:knon) * Kech_T_px(1:knon) 392 HTphiT_w(1:knon) = C_p(1:knon) * Kech_T_pw(1:knon) 393 ! Moisture flux equation 394 HTphiQ_x(1:knon) = beta(1:knon) * Kech_Q_sx(1:knon) * dqsatdT0_x(1:knon) 395 HTphiQ_w(1:knon) = beta(1:knon) * Kech_Q_sw(1:knon) * dqsatdT0_w(1:knon) 396 ! Radiation equation 397 HTRn_x(1:knon) = Rp1_x(1:knon) * Kech_T_px(1:knon) * BT_x(1:knon) * dtime + Rps_x(1:knon) 398 HTRn_w(1:knon) = Rp1_w(1:knon) * Kech_T_pw(1:knon) * BT_w(1:knon) * dtime + Rps_w(1:knon) 399 400 ! DFluxDQ coefficients 401 ! --------------------- 402 ! Heat flux equation 403 HQphiT_x(1:knon) = C_p(1:knon) * Kech_T_px(1:knon) * QQ_x(1:knon) 404 HQphiT_w(1:knon) = C_p(1:knon) * Kech_T_pw(1:knon) * QQ_w(1:knon) 405 ! Moisture flux equation 406 HQphiQ_x(1:knon) = beta(1:knon) * Kech_Q_sx(1:knon) 407 HQphiQ_w(1:knon) = beta(1:knon) * Kech_Q_sw(1:knon) 408 ! Radiation equation 409 HQRn_x(1:knon) = (Rp1_x(1:knon) * Kech_T_px(1:knon) * BT_x(1:knon) * dtime + Rps_x(1:knon)) * QQ_x(1:knon) 410 HQRn_w(1:knon) = (Rp1_w(1:knon) * Kech_T_pw(1:knon) * BT_w(1:knon) * dtime + Rps_w(1:knon)) * QQ_w(1:knon) 411 412 ! Mean values and w-x differences 413 ! ------------------------------- 414 phiT0_b(1:knon) = sigw(1:knon) * phiT0_w(1:knon) + sigx(1:knon) * phiT0_x(1:knon) 415 phiQ0_b(1:knon) = sigw(1:knon) * phiQ0_w(1:knon) + sigx(1:knon) * phiQ0_x(1:knon) 416 Rn0_b(1:knon) = sigw(1:knon) * Rn0_w(1:knon) + sigx(1:knon) * Rn0_x(1:knon) 417 418 dphiT0(1:knon) = phiT0_w(1:knon) - phiT0_x(1:knon) 419 dphiQ0(1:knon) = phiQ0_w(1:knon) - phiQ0_x(1:knon) 420 dRn0(1:knon) = Rn0_w(1:knon) - Rn0_x(1:knon) 421 422 HTphiT_b(1:knon) = sigw(1:knon) * HTphiT_w(1:knon) + sigx(1:knon) * HTphiT_x(1:knon) 423 dd_HTphiT(1:knon) = HTphiT_w(1:knon) - HTphiT_x(1:knon) 424 425 HTphiQ_b(1:knon) = sigw(1:knon) * HTphiQ_w(1:knon) + sigx(1:knon) * HTphiQ_x(1:knon) 426 dd_HTphiQ(1:knon) = HTphiQ_w(1:knon) - HTphiQ_x(1:knon) 427 428 HTRn_b(1:knon) = sigw(1:knon) * HTRn_w(1:knon) + sigx(1:knon) * HTRn_x(1:knon) 429 dd_HTRn(1:knon) = HTRn_w(1:knon) - HTRn_x(1:knon) 430 431 HQphiT_b(1:knon) = sigw(1:knon) * HQphiT_w(1:knon) + sigx(1:knon) * HQphiT_x(1:knon) 432 dd_HQphiT(1:knon) = HQphiT_w(1:knon) - HQphiT_x(1:knon) 433 434 HQphiQ_b(1:knon) = sigw(1:knon) * HQphiQ_w(1:knon) + sigx(1:knon) * HQphiQ_x(1:knon) 435 dd_HQphiQ(1:knon) = HQphiQ_w - HQphiQ_x(1:knon) 436 437 HQRn_b(1:knon) = sigw(1:knon) * HQRn_w(1:knon) + sigx(1:knon) * HQRn_x(1:knon) 438 dd_HQRn(1:knon) = HQRn_w(1:knon) - HQRn_x(1:knon) 439 440 ! Feedback Gains 441 ! -------------- 442 g_T(1:knon) = - (sqrt(tau(1:knon)) / Inert(1:knon)) & 443 * (HTphiT_b(1:knon) + L_v(1:knon) * HTphiQ_b(1:knon) + HTRn_b(1:knon) & 444 + (dd_HTphiT(1:knon) + L_v(1:knon) * dd_HTphiQ(1:knon) + dd_HTRn(1:knon)) & 445 * (sigx(1:knon) - sigw(1:knon) - sigw(1:knon) * sigx(1:knon) * dd_HTphiT(1:knon) / HTphiT_b(1:knon))) 446 447 !!!! DO j = 1,knon 448 !!!! IF (mod(j,20) .EQ.0) THEN 449 !!!! print *, ' j dd_QQ QQ_b dd_HQphiQ dd_HQphiT dd_HQRn HQphiQ_b HQphiT_b HQRn_b ' 450 !!!! ENDIF 451 !!!! print 1789, j, dd_QQ(j), QQ_b(j), dd_HQphiQ(j), dd_HQphiT(j), dd_HQRn(j), HQphiQ_b(j), HQphiT_b(j), HQRn_b(j) 452 !!!! 1789 FORMAT( I4, 10(1X,E10.2)) 453 !!!! ENDDO 454 g_Q(1:knon) = - (dd_QQ(1:knon) / QQ_b(1:knon)) * & 455 (sigx(1:knon) - sigw(1:knon) - sigw(1:knon) * sigx(1:knon) * dd_KQs(1:knon) / Kech_Qs(1:knon)) & 456 - sqrt(tau(1:knon)) / (Inert(1:knon) * QQ_b(1:knon)) * & 457 (HQphiT_b(1:knon) + L_v(1:knon) * HQphiQ_b(1:knon) + HQRn_b(1:knon) + & 458 (sigx(1:knon) - sigw(1:knon) - sigw(1:knon) * sigx(1:knon) * dd_KQs(1:knon) / Kech_Qs(1:knon)) * & 459 (dd_HQphiT(1:knon) + L_v(1:knon) * dd_HQphiQ(1:knon) + dd_HQRn(1:knon))) 460 461 !! g_Q(1:knon) = - (dd_QQ(1:knon)/QQ_b(1:knon)) * & 462 !! (sigx(1:knon)-sigw(1:knon)-sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) & 463 !! - sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon)) * & 464 !! ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) + & 465 !! (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) * & 466 !! (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) 467 468 !! g_Q(1:knon) = - (sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon))) * & 469 !! ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) ) & 470 !! - (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) * & 471 !! ( dd_QQ(1:knon)/QQ_b(1:knon) & 472 !! + (sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon))) & 473 !! * (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) 474 475 ! Ts, qs Coupling coefficients / 476 ! ---------------------------- 477 Gamma_phiT(1:knon) = (sqrt(tau(1:knon)) / (Inert(1:knon) * HTphiT_b(1:knon))) & 478 * (dd_HTphiT(1:knon) + L_v(1:knon) * dd_HTphiQ(1:knon) + dd_HTRn(1:knon)) 479 480 Gamma_phiQ(1:knon) = (1. / (Kech_Qs(1:knon) * QQ_b(1:knon))) * & 481 (dd_QQ(1:knon) & 482 + (sqrt(tau(1:knon)) / (Inert(1:knon))) * & 483 (dd_HQphiT(1:knon) + L_v(1:knon) * dd_HQphiQ(1:knon) + dd_HQRn(1:knon))) 484 485 !! Gamma_phiQ(1:knon) = (beta(1:knon)/(HQphiQ_b(1:knon)*QQ_b(1:knon))) * & 486 !! ( dd_QQ(1:knon) & 487 !! + (sqrt(tau(1:knon))/(Inert(1:knon))) * & 488 !! (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) 489 490 !! Gamma_phiQ(1:knon) = (1/(HQphiQ_b(1:knon)*QQ_b(1:knon))) & 491 !! * ( dd_QQ(1:knon) & 492 !! + (sqrt(tau(1:knon))/(Inert(1:knon))) & 493 !! * (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) 494 495 ! Insensitive changes 496 ! ------------------- 497 dTs_ins(1:knon) = (sqrt(tau(1:knon)) / Inert(1:knon)) * & 498 (dphiT0(1:knon) + L_v(1:knon) * dphiQ0(1:knon) + dRn0(1:knon)) 499 500 dqsatsrf_ins(1:knon) = (beta(1:knon) / QQ_b(1:knon)) * dTs_ins(1:knon) 501 502 IF (prt_level >= 10) THEN 503 print *, 'wx_pbl_merge, tau ', tau 504 print *, 'wx_pbl_merge, AcoefT0 ', AcoefT0 505 print *, 'wx_pbl_merge, AcoefQ0 ', AcoefQ0 506 print *, 'wx_pbl_merge, BcoefT0 ', BcoefT0 507 print *, 'wx_pbl_merge, BcoefQ0 ', BcoefQ0 508 print *, 'wx_pbl_merge, qsat0_w, qsat0_x ', (qsat0_w(j), qsat0_x(j), j = 1, knon) 509 print *, 'wx_pbl_merge, dqsatdT0_w, dqsatdT0_x ', (dqsatdT0_w(j), dqsatdT0_x(j), j = 1, knon) 510 ENDIF 511 512 !---------------------------------------------------------------------------- 513 514 !------------------------------------------------------------------------------ 515 516 ! Effective coefficients Acoef and Bcoef 517 ! -------------------------------------- 518 DO j = 1, knon 519 AcoefT(j) = AcoefT0(j) - sigw(j) * sigx(j) * (dd_KTp(j) / Kech_Tp(j)) * C_p(j) * & 520 (dTs0(j) + (dTs_ins(j) - dTs0(j) - Gamma_phiT(j) * phiT0_b(j)) / (1. - g_T(j))) 521 BcoefT(j) = BcoefT0(j) - sigw(j) * sigx(j) * (dd_KTp(j) / Kech_Tp(j)) * C_p(j) * Gamma_phiT(j) / (1. - g_T(j)) / dtime 522 523 AcoefQ(j) = AcoefQ0(j) - sigw(j) * sigx(j) * (dd_KQs(j) / Kech_Qs(j)) * & 524 (dqsatsrf0(j) + (dqsatsrf_ins(j) - (beta(j) / QQ_b(j)) * dTs0(j) - Gamma_phiQ(j) * phiQ0_b(j)) / (1 - g_Q(j))) / & 525 max(beta(j), 1.e-4) 526 BcoefQ(j) = BcoefQ0(j) - sigw(j) * sigx(j) * (dd_KQs(j) / Kech_Qs(j)) * Gamma_phiQ(j) / (1 - g_Q(j)) / & 527 (max(beta(j), 1.e-4) * dtime) 528 !! AcoefQ(j) = AcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))* & 529 !! (dqsatsrf0(j) + (dqsatsrf_ins(j)-(beta(j)/QQ_b(j))*dTs0(j)-Gamma_phiQ(j)*phiQ0_b(j))/(1 - g_Q(j)))/ & 530 !! beta(j) 531 !! BcoefQ(j) = BcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))*Gamma_phiQ(j)/(1 - g_Q(j))/(beta(j)*dtime) 532 ENDDO ! j = 1,knon 533 534 IF (prt_level >= 10) THEN 535 print *, 'wx_pbl_dts AAAA BcoefQ, BcoefQ0, sigw ', & 536 BcoefQ, BcoefQ0, sigw 537 print *, 'wx_pbl_dts_merge, dTs_ins ', dTs_ins 538 print *, 'wx_pbl_dts_merge, dqs_ins ', dqsatsrf_ins 539 ENDIF 540 541 END SUBROUTINE wx_pbl_dts_merge 542 543 SUBROUTINE wx_pbl_split(knon, nsrf, dtime, sigw, beta, iflag_split, & 544 g_T, g_Q, & 545 Gamma_phiT, Gamma_phiQ, & 546 dTs_ins, dqsatsrf_ins, & 547 phiT, phiQ, phiU, phiV, & 548 !!!! HTRn_b, dd_HTRn, HTphiT_b, dd_HTphiT, & 549 phiQ0_b, phiT0_b, & 550 phiT_x, phiT_w, & 551 phiQ_x, phiQ_w, & 552 phiU_x, phiU_w, & 553 phiV_x, phiV_w, & 554 philat_x, philat_w, & 555 !!!! Rn_b, dRn, & 556 dqsatsrf, & 557 dTs, delta_qsurf & 558 ) 565 559 566 560 USE wx_pbl_var_mod 567 561 568 USE lmdz_print_control, ONLY: prt_level, lunout562 USE lmdz_print_control, ONLY: prt_level, lunout 569 563 USE indice_sol_mod, ONLY: is_oce 570 571 INCLUDE "YOMCST.h" 572 573 INTEGER, INTENT(IN) :: knon ! number of grid cells 574 INTEGER, INTENT(IN) :: nsrf ! surface type 575 REAL, INTENT(IN) :: dtime ! time step size (s) 576 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area 577 REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor 578 INTEGER, INTENT(IN) :: iflag_split 579 REAL, DIMENSION(knon), INTENT(IN) :: g_T, g_Q 580 REAL, DIMENSION(knon), INTENT(IN) :: Gamma_phiT, Gamma_phiQ 581 REAL, DIMENSION(knon), INTENT(IN) :: dTs_ins, dqsatsrf_ins 582 REAL, DIMENSION(knon), INTENT(IN) :: phiT, phiQ, phiU, phiV 583 REAL, DIMENSION(knon), INTENT(IN) :: phiQ0_b, phiT0_b 584 585 REAL, DIMENSION(knon), INTENT(OUT) :: phiT_x, phiT_w 586 REAL, DIMENSION(knon), INTENT(OUT) :: phiQ_x, phiQ_w 587 REAL, DIMENSION(knon), INTENT(OUT) :: phiU_x, phiU_w 588 REAL, DIMENSION(knon), INTENT(OUT) :: phiV_x, phiV_w 589 REAL, DIMENSION(knon), INTENT(OUT) :: philat_x, philat_w 590 REAL, DIMENSION(knon), INTENT(OUT) :: dqsatsrf ! beta delta(qsat(Ts)) 591 REAL, DIMENSION(knon), INTENT(OUT) :: dTs ! Temperature difference at surface 592 REAL, DIMENSION(knon), INTENT(OUT) :: delta_qsurf 593 594 !! Local variables 595 INTEGER :: j 596 REAL, DIMENSION(knon) :: dphiT, dphiQ, dphiU, dphiV 597 REAL, DIMENSION(knon) :: q1_x, q1_w 598 599 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region 600 601 !---------------------------------------------------------------------------- 602 ! Equations 603 ! --------- 604 !!!!!! (1 - g_T) dTs = dTs_ins + Gamma_phiT phiT 605 !!!!!! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ 606 !!!!!! dphiT = (dd_KTp/KTp) phiT + ( dd_AT - C_p dTs)*KxKwTp/KTp 607 !!!!!! dphiQ = (dd_KQs/KQs) phiQ + (beta dd_AQ - dqsatsrf )*KxKwQs/KQs 608 !!!!!! dphiU = (dd_KUp/KUp) phiU + ( dd_AU )*KxKwUp/KUp 609 !!!!!! dphiV = (dd_KVp/KVp) phiV + ( dd_AV )*KxKwVp/KVp 610 611 ! (1 - g_T) (dTs-dTs0) = dTs_ins-dTs0 + Gamma_phiT (phiT-phiT0) 612 ! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ 613 ! dphiT = (dd_KTp/KTp) phiT + ( dd_AT - C_p dTs)*KxKwTp/KTp 614 ! dphiQ = (dd_KQs/KQs) phiQ + (beta dd_AQ - dqsatsrf )*KxKwQs/KQs 615 ! dphiU = (dd_KUp/KUp) phiU + ( dd_AU )*KxKwUp/KUp 616 ! dphiV = (dd_KVp/KVp) phiV + ( dd_AV )*KxKwVp/KVp 617 618 !! 619 sigx(:) = 1.-sigw(:) 620 621 ! print *,' AAAA wx_pbl_split, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:) 622 623 IF (iflag_split == 2 .AND. nsrf /= is_oce) THEN 624 625 ! Delta_tsurf and Delta_qsurf computation 626 ! ----------------------------------------- 627 IF (prt_level >=10 ) THEN 628 print *,' wx_pbl_split, dTs_ins, dTs0 , Gamma_phiT, g_T ', dTs_ins, dTs0, Gamma_phiT, g_T 629 print *,' wx_pbl_split, dqsatsrf_ins, Gamma_phiQ, g_q ', dqsatsrf_ins, Gamma_phiQ, g_q 564 USE lmdz_yomcst 565 566 IMPLICIT NONE 567 568 INTEGER, INTENT(IN) :: knon ! number of grid cells 569 INTEGER, INTENT(IN) :: nsrf ! surface type 570 REAL, INTENT(IN) :: dtime ! time step size (s) 571 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area 572 REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor 573 INTEGER, INTENT(IN) :: iflag_split 574 REAL, DIMENSION(knon), INTENT(IN) :: g_T, g_Q 575 REAL, DIMENSION(knon), INTENT(IN) :: Gamma_phiT, Gamma_phiQ 576 REAL, DIMENSION(knon), INTENT(IN) :: dTs_ins, dqsatsrf_ins 577 REAL, DIMENSION(knon), INTENT(IN) :: phiT, phiQ, phiU, phiV 578 REAL, DIMENSION(knon), INTENT(IN) :: phiQ0_b, phiT0_b 579 580 REAL, DIMENSION(knon), INTENT(OUT) :: phiT_x, phiT_w 581 REAL, DIMENSION(knon), INTENT(OUT) :: phiQ_x, phiQ_w 582 REAL, DIMENSION(knon), INTENT(OUT) :: phiU_x, phiU_w 583 REAL, DIMENSION(knon), INTENT(OUT) :: phiV_x, phiV_w 584 REAL, DIMENSION(knon), INTENT(OUT) :: philat_x, philat_w 585 REAL, DIMENSION(knon), INTENT(OUT) :: dqsatsrf ! beta delta(qsat(Ts)) 586 REAL, DIMENSION(knon), INTENT(OUT) :: dTs ! Temperature difference at surface 587 REAL, DIMENSION(knon), INTENT(OUT) :: delta_qsurf 588 589 !! Local variables 590 INTEGER :: j 591 REAL, DIMENSION(knon) :: dphiT, dphiQ, dphiU, dphiV 592 REAL, DIMENSION(knon) :: q1_x, q1_w 593 594 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region 595 596 !---------------------------------------------------------------------------- 597 ! Equations 598 ! --------- 599 !!!!!! (1 - g_T) dTs = dTs_ins + Gamma_phiT phiT 600 !!!!!! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ 601 !!!!!! dphiT = (dd_KTp/KTp) phiT + ( dd_AT - C_p dTs)*KxKwTp/KTp 602 !!!!!! dphiQ = (dd_KQs/KQs) phiQ + (beta dd_AQ - dqsatsrf )*KxKwQs/KQs 603 !!!!!! dphiU = (dd_KUp/KUp) phiU + ( dd_AU )*KxKwUp/KUp 604 !!!!!! dphiV = (dd_KVp/KVp) phiV + ( dd_AV )*KxKwVp/KVp 605 606 ! (1 - g_T) (dTs-dTs0) = dTs_ins-dTs0 + Gamma_phiT (phiT-phiT0) 607 ! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ 608 ! dphiT = (dd_KTp/KTp) phiT + ( dd_AT - C_p dTs)*KxKwTp/KTp 609 ! dphiQ = (dd_KQs/KQs) phiQ + (beta dd_AQ - dqsatsrf )*KxKwQs/KQs 610 ! dphiU = (dd_KUp/KUp) phiU + ( dd_AU )*KxKwUp/KUp 611 ! dphiV = (dd_KVp/KVp) phiV + ( dd_AV )*KxKwVp/KVp 612 613 !! 614 sigx(:) = 1. - sigw(:) 615 616 ! print *,' AAAA wx_pbl_split, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:) 617 618 IF (iflag_split == 2 .AND. nsrf /= is_oce) THEN 619 620 ! Delta_tsurf and Delta_qsurf computation 621 ! ----------------------------------------- 622 IF (prt_level >=10) THEN 623 print *, ' wx_pbl_split, dTs_ins, dTs0 , Gamma_phiT, g_T ', dTs_ins, dTs0, Gamma_phiT, g_T 624 print *, ' wx_pbl_split, dqsatsrf_ins, Gamma_phiQ, g_q ', dqsatsrf_ins, Gamma_phiQ, g_q 630 625 ENDIF 631 626 632 DO j = 1, knon633 dTs(j) = dTs0(j) + (dTs_ins(j) - dTs0(j) + Gamma_phiT(j)*(phiT(j)-phiT0_b(j)) )/(1 - g_T(j))634 dqsatsrf(j) = dqsatsrf0(j) + (dqsatsrf_ins(j) - (beta(j) /QQ_b(j))*dTs0(j) + &635 Gamma_phiQ(j)*(phiQ(j)-phiQ0_b(j)) )/(1 - g_Q(j))627 DO j = 1, knon 628 dTs(j) = dTs0(j) + (dTs_ins(j) - dTs0(j) + Gamma_phiT(j) * (phiT(j) - phiT0_b(j))) / (1 - g_T(j)) 629 dqsatsrf(j) = dqsatsrf0(j) + (dqsatsrf_ins(j) - (beta(j) / QQ_b(j)) * dTs0(j) + & 630 Gamma_phiQ(j) * (phiQ(j) - phiQ0_b(j))) / (1 - g_Q(j)) 636 631 ENDDO ! j = 1,knon 637 632 638 IF (prt_level >=10 ) THEN 639 print *,' wx_pbl_split, dqsatsrf0, QQ_b ', dqsatsrf0, QQ_b 640 print *,' wx_pbl_split, phiT0_b, phiT, dTs ', phiT0_b, phiT, dTs 641 print *,' wx_pbl_split, phiQ0_b, phiQ, dqsatsrf ', phiQ0_b, phiQ, dqsatsrf 642 ENDIF 643 ELSE 644 dTs(:) = 0. 645 dqsatsrf(:) = 0. 646 ENDIF ! (iflag_split .EQ. 2 .AND. nsrf .NE. is_oce) 647 648 DO j = 1,knon 649 dphiT(j) = (phiT(j)*dd_KTp(j) + ( dd_AT(j) - C_p(j)*dTs(j))*KxKwTp(j))/Kech_Tp(j) 650 dphiQ(j) = (phiQ(j)*dd_KQs(j) + (beta(j)*dd_AQ(j) - dqsatsrf(j))*KxKwQs(j))/Kech_Qs(j) 651 dphiU(j) = (phiU(j)*dd_KUp(j) + dd_AU(j) *KxKwUp(j))/Kech_Up(j) 652 dphiV(j) = (phiV(j)*dd_KVp(j) + dd_AV(j) *KxKwVp(j))/Kech_Vp(j) 653 654 phiT_x(j)=phiT(j) - sigw(j)*dphiT(j) 655 phiT_w(j)=phiT(j) + sigx(j)*dphiT(j) 656 phiQ_x(j)=phiQ(j) - sigw(j)*dphiQ(j) 657 phiQ_w(j)=phiQ(j) + sigx(j)*dphiQ(j) 658 phiU_x(j)=phiU(j) - sigw(j)*dphiU(j) 659 phiU_w(j)=phiU(j) + sigx(j)*dphiU(j) 660 phiV_x(j)=phiV(j) - sigw(j)*dphiV(j) 661 phiV_w(j)=phiV(j) + sigx(j)*dphiV(j) 662 663 philat_x(j)=phiQ_x(j)*RLVTT 664 philat_w(j)=phiQ_w(j)*RLVTT 665 ENDDO ! j = 1,knon 666 667 DO j = 1,knon 668 q1_x(j) = AQ_x(j) + BQ_x(j)*phiQ_x(j)*dtime 669 q1_w(j) = AQ_w(j) + BQ_w(j)*phiQ_w(j)*dtime 670 ENDDO ! j = 1,knon 671 DO j = 1,knon 672 delta_qsurf(j) = (1.-beta(j))*(q1_w(j) - q1_x(j)) + dqsatsrf(j) 673 ENDDO ! j = 1,knon 674 675 !! Do j = 1,knon 676 !! print *,'XXXsplit : j, q1_x(j), AQ_x(j), BQ_x(j), phiQ_x(j) ', j, q1_x(j), AQ_x(j), BQ_x(j), phiQ_x(j) 677 !! print *,'XXXsplit : j, q1_w(j), AQ_w(j), BQ_w(j), phiQ_w(j) ', j, q1_w(j), AQ_w(j), BQ_w(j), phiQ_w(j) 678 !! ENDDO 679 680 IF (prt_level >=10 ) THEN 681 print *,' wx_pbl_split, phiT, dphiT, dTs ', phiT, dphiT, dTs 682 print *,' wx_pbl_split, phiQ, dphiQ, dqsatsrf ', phiQ, dphiQ, dqsatsrf 683 ENDIF 684 685 IF (prt_level >=10 ) THEN 686 !! print *,' wx_pbl_split, verif dqsatsrf = beta dqsatdT0 dTs ' 687 !! print *,' wx_pbl_split, dqsatsrf, dqsatdT0*dTs ', dqsatsrf, dqsatdT0*dTs 688 ENDIF 689 690 !! IF (knon .NE. 0) THEN 691 !! call iophys_ecrit('sigw', 1,'sigw', '.',sigw) 692 !! call iophys_ecrit('phit', 1,'phit', 'W/m2',phit) 693 !! call iophys_ecrit('phit_w', 1,'phit_w', 'W/m2',phit_w) 694 !! call iophys_ecrit('phit_x', 1,'phit_x', 'W/m2',phit_x) 695 !! call iophys_ecrit('phiq', 1,'phiq', 'kg/m2/s',phiq) 696 !! call iophys_ecrit('phiq_w', 1,'phiq_w', 'kg/m2/s',phiq_w) 697 !! call iophys_ecrit('phiq_x', 1,'phiq_x', 'kg/m2/s',phiq_x) 698 !! call iophys_ecrit('q1_w', 1,'q1_w', '.',q1_w) 699 !! call iophys_ecrit('q1_x', 1,'q1_x', '.',q1_x) 700 !! ENDIF ! (knon .NE. 0) 701 702 703 704 END SUBROUTINE wx_pbl_split 705 706 SUBROUTINE wx_pbl_check( knon, dtime, ypplay, ypaprs, & 707 sigw, beta, iflag_split, & 708 Ts0_b9, dTs09, & 709 qs_b9, Ts_b9, & ! yqsurf, Tsurf_new 710 dTs9, dqsatsrf9, & 711 AcoefT_x, AcoefT_w, & 712 BcoefT_x, BcoefT_w, & 713 AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, & 714 AcoefT, AcoefQ, BcoefT, BcoefQ, & 715 phiT_b9, phiQ_b9, & 716 phiT_x9, phiT_w9, & 717 phiQ_x9, phiQ_w9 & 718 ) 633 IF (prt_level >=10) THEN 634 print *, ' wx_pbl_split, dqsatsrf0, QQ_b ', dqsatsrf0, QQ_b 635 print *, ' wx_pbl_split, phiT0_b, phiT, dTs ', phiT0_b, phiT, dTs 636 print *, ' wx_pbl_split, phiQ0_b, phiQ, dqsatsrf ', phiQ0_b, phiQ, dqsatsrf 637 ENDIF 638 ELSE 639 dTs(:) = 0. 640 dqsatsrf(:) = 0. 641 ENDIF ! (iflag_split .EQ. 2 .AND. nsrf .NE. is_oce) 642 643 DO j = 1, knon 644 dphiT(j) = (phiT(j) * dd_KTp(j) + (dd_AT(j) - C_p(j) * dTs(j)) * KxKwTp(j)) / Kech_Tp(j) 645 dphiQ(j) = (phiQ(j) * dd_KQs(j) + (beta(j) * dd_AQ(j) - dqsatsrf(j)) * KxKwQs(j)) / Kech_Qs(j) 646 dphiU(j) = (phiU(j) * dd_KUp(j) + dd_AU(j) * KxKwUp(j)) / Kech_Up(j) 647 dphiV(j) = (phiV(j) * dd_KVp(j) + dd_AV(j) * KxKwVp(j)) / Kech_Vp(j) 648 649 phiT_x(j) = phiT(j) - sigw(j) * dphiT(j) 650 phiT_w(j) = phiT(j) + sigx(j) * dphiT(j) 651 phiQ_x(j) = phiQ(j) - sigw(j) * dphiQ(j) 652 phiQ_w(j) = phiQ(j) + sigx(j) * dphiQ(j) 653 phiU_x(j) = phiU(j) - sigw(j) * dphiU(j) 654 phiU_w(j) = phiU(j) + sigx(j) * dphiU(j) 655 phiV_x(j) = phiV(j) - sigw(j) * dphiV(j) 656 phiV_w(j) = phiV(j) + sigx(j) * dphiV(j) 657 658 philat_x(j) = phiQ_x(j) * RLVTT 659 philat_w(j) = phiQ_w(j) * RLVTT 660 ENDDO ! j = 1,knon 661 662 DO j = 1, knon 663 q1_x(j) = AQ_x(j) + BQ_x(j) * phiQ_x(j) * dtime 664 q1_w(j) = AQ_w(j) + BQ_w(j) * phiQ_w(j) * dtime 665 ENDDO ! j = 1,knon 666 DO j = 1, knon 667 delta_qsurf(j) = (1. - beta(j)) * (q1_w(j) - q1_x(j)) + dqsatsrf(j) 668 ENDDO ! j = 1,knon 669 670 !! Do j = 1,knon 671 !! print *,'XXXsplit : j, q1_x(j), AQ_x(j), BQ_x(j), phiQ_x(j) ', j, q1_x(j), AQ_x(j), BQ_x(j), phiQ_x(j) 672 !! print *,'XXXsplit : j, q1_w(j), AQ_w(j), BQ_w(j), phiQ_w(j) ', j, q1_w(j), AQ_w(j), BQ_w(j), phiQ_w(j) 673 !! ENDDO 674 675 IF (prt_level >=10) THEN 676 print *, ' wx_pbl_split, phiT, dphiT, dTs ', phiT, dphiT, dTs 677 print *, ' wx_pbl_split, phiQ, dphiQ, dqsatsrf ', phiQ, dphiQ, dqsatsrf 678 ENDIF 679 680 IF (prt_level >=10) THEN 681 !! print *,' wx_pbl_split, verif dqsatsrf = beta dqsatdT0 dTs ' 682 !! print *,' wx_pbl_split, dqsatsrf, dqsatdT0*dTs ', dqsatsrf, dqsatdT0*dTs 683 ENDIF 684 685 !! IF (knon .NE. 0) THEN 686 !! call iophys_ecrit('sigw', 1,'sigw', '.',sigw) 687 !! call iophys_ecrit('phit', 1,'phit', 'W/m2',phit) 688 !! call iophys_ecrit('phit_w', 1,'phit_w', 'W/m2',phit_w) 689 !! call iophys_ecrit('phit_x', 1,'phit_x', 'W/m2',phit_x) 690 !! call iophys_ecrit('phiq', 1,'phiq', 'kg/m2/s',phiq) 691 !! call iophys_ecrit('phiq_w', 1,'phiq_w', 'kg/m2/s',phiq_w) 692 !! call iophys_ecrit('phiq_x', 1,'phiq_x', 'kg/m2/s',phiq_x) 693 !! call iophys_ecrit('q1_w', 1,'q1_w', '.',q1_w) 694 !! call iophys_ecrit('q1_x', 1,'q1_x', '.',q1_x) 695 !! ENDIF ! (knon .NE. 0) 696 697 END SUBROUTINE wx_pbl_split 698 699 SUBROUTINE wx_pbl_check(knon, dtime, ypplay, ypaprs, & 700 sigw, beta, iflag_split, & 701 Ts0_b9, dTs09, & 702 qs_b9, Ts_b9, & ! yqsurf, Tsurf_new 703 dTs9, dqsatsrf9, & 704 AcoefT_x, AcoefT_w, & 705 BcoefT_x, BcoefT_w, & 706 AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, & 707 AcoefT, AcoefQ, BcoefT, BcoefQ, & 708 phiT_b9, phiQ_b9, & 709 phiT_x9, phiT_w9, & 710 phiQ_x9, phiQ_w9 & 711 ) 719 712 720 713 USE wx_pbl_var_mod 721 714 722 USE lmdz_print_control, ONLY: prt_level, lunout723 USE lmdz_ YOETHF715 USE lmdz_print_control, ONLY: prt_level, lunout 716 USE lmdz_yoethf 724 717 USE lmdz_fcttre, ONLY: foeew, foede, qsats, qsatl, dqsats, dqsatl, thermcep 725 726 INCLUDE "YOMCST.h" 727 728 INTEGER, INTENT(IN) :: knon ! number of grid cells 729 REAL, INTENT(IN) :: dtime ! time step size (s) 730 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) 731 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) 732 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area 733 REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor 734 INTEGER, INTENT(IN) :: iflag_split 735 REAL, DIMENSION(knon), INTENT(IN) :: Ts0_b9, dTs09 736 REAL, DIMENSION(knon), INTENT(IN) :: qs_b9, Ts_b9 ! yqsurf, Tsurf_new 737 REAL, DIMENSION(knon), INTENT(IN) :: dTs9, dqsatsrf9 738 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w 739 REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w 740 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0 741 742 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT, AcoefQ, BcoefT, BcoefQ 743 REAL, DIMENSION(knon), INTENT(IN) :: phiT_b9, phiQ_b9 744 REAL, DIMENSION(knon), INTENT(IN) :: phiT_x9, phiT_w9 745 REAL, DIMENSION(knon), INTENT(IN) :: phiQ_x9, phiQ_w9 746 747 !! Local variables 748 INTEGER :: j 749 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region 750 REAL, DIMENSION(knon) :: AcoefT_b, AcoefQ_b ! mean values of AcoefT and AcoefQ 751 REAL :: zzt, zzq, zzqsat 752 REAL :: zdelta, zcvm5, zcor, qsat 753 REAL, DIMENSION(knon) :: qsat_w, qsat_x 754 REAL, DIMENSION(knon) :: dqsatdT_w, dqsatdT_x 755 REAL, DIMENSION(knon) :: qsat_bs ! qsat(Ts_b) 756 REAL, DIMENSION(knon) :: qsat01, dqsatdT01 757 REAL, DIMENSION(knon) :: Ts_x, Ts_w, qs_x, qs_w 758 REAL, DIMENSION(knon) :: T1_x, T1_w, q1_x, q1_w 759 REAL, DIMENSION(knon) :: Rn_x, Rn_w 760 REAL, DIMENSION(knon) :: phiQ0_x, phiQ0_w 761 REAL, DIMENSION(knon) :: Ta, qa 762 REAL, DIMENSION(knon) :: qsatsrf_w, qsatsrf_x, qsatsrf_b 763 REAL, DIMENSION(knon) :: qsurf_w, qsurf_x 764 REAL :: dphiT, dphiQ 765 REAL :: dqsatsrf1 766 REAL :: phiT_w1, phiT_w2 767 REAL :: phiT_x1, phiT_x2 768 REAL :: phiQ_w1, phiQ_w2, phiQ_w3 769 REAL :: phiQ_x1, phiQ_x2, phiQ_x3 770 REAL :: phiT_b1, phiQ_b1 771 REAL :: Kech_Q_sw1, Kech_Q_sx1 772 REAL :: evap_pot 773 774 !---------------------------------------------------------------------------- 775 ! Equations to be checked: 776 ! ----------------------- 777 ! Input : Ts0_b, dTs0, Ts_b, dTs, qsatsrf_b, dqsatsrf, 778 ! phiT_b, phiQ_b, phiT_w, phiT_x, phiQ_w, phiQ_x, 779 780 ! AcoefT, AcoefQ, AcoefT_w, AcoefQ_w, AcoefT_x, AcoefQ_x, 781 ! BcoefT, BcoefQ, BcoefT_w, BcoefQ_w, BcoefT_x, BcoefQ_x 782 783 ! C_p T1_w = AcoefT_w + BcoefT_w phiT_w Delta t C_p T1_x = AcoefT_x + BcoefT_x phiT_x Delta t 784 ! q1_w = AQ_w + BQ_w phiQ_w Delta t q1_x = AQ_x + BQ_x phiQ_x Delta t 785 ! qsatsrf_w = beta qsat(Ts_w) qsatsrf_x = beta qsat(Ts_x) 786 ! qsurf_w = (1-beta) q1_w + qsatsrf_w qsurf_x = (1-beta) q1_x + qsatsrf_x 787 ! phiT_w = Kech_h_w C_p ( T1_w - Ts_w) phiT_x = Kech_h_x C_p ( T1_x - Ts_x) 788 ! phiT_w = Kech_T_pw ( AcoefT_w - C_p Ts_w) phiT_x = Kech_T_px ( AcoefT_x - C_p Ts_x) 789 ! phiq_w = Kech_h_w ( beta q1_w - qsatsrf_w) phiq_x = Kech_h_x ( beta q1_x - qsatsrf_x)) 790 ! phiq_w = Kech_Q_sw (beta AQ_w -qsatsrf_w) phiq_x = Kech_Q_sx (beta AQ_x -qsatsrf_x) 791 ! phiq_w = Kech_h_w (q1_w - qsurf_w) phiq_x = Kech_h_x (q1_x - qsurf_x) 792 ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x 793 ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x 794 ! Ts_b = sigw Ts_w + sigx Ts_x dTs = Ts_w - Ts_x 795 ! qsatsrf_b = sigw qsatsrf_w + sigx qsatsrf_x 796 ! C_p Ta = AcoefT + BcoefT phiT_b Delta t 797 ! qa = AcoefQ + BcoefQ phiQ_b Delta t 798 ! phiT_b = Kech_h C_p (Ta - Ts_b) 799 ! phiQ_b = beta Kech_h (qa - qsatsrf_b) 800 ! dTs = sqrt(tau)/I (dphit + L_v dphiq + dR) 801 802 !---------------------------------------------------------------------------- 803 804 !! 805 sigx(:) = 1.-sigw(:) 806 AcoefT_b(1:knon) = AcoefT_x(1:knon) + sigw(1:knon)*dd_AT(1:knon) 807 AcoefQ_b(1:knon) = AQ_x(1:knon) + sigw(1:knon)*dd_AQ(1:knon) 808 809 ! Compute the three qsat and dqsatdTs 810 ! --------------------------------------------- 811 !! C_p(1:knon) = RCpd 812 !! L_v(1:knon) = RLvtt 813 IF (prt_level >=10 ) THEN 814 print *,' AAAA wx_pbl_check, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:) 718 USE lmdz_yomcst 719 720 IMPLICIT NONE 721 722 INTEGER, INTENT(IN) :: knon ! number of grid cells 723 REAL, INTENT(IN) :: dtime ! time step size (s) 724 REAL, DIMENSION(knon, klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) 725 REAL, DIMENSION(knon, klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) 726 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area 727 REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor 728 INTEGER, INTENT(IN) :: iflag_split 729 REAL, DIMENSION(knon), INTENT(IN) :: Ts0_b9, dTs09 730 REAL, DIMENSION(knon), INTENT(IN) :: qs_b9, Ts_b9 ! yqsurf, Tsurf_new 731 REAL, DIMENSION(knon), INTENT(IN) :: dTs9, dqsatsrf9 732 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w 733 REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w 734 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0 735 736 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT, AcoefQ, BcoefT, BcoefQ 737 REAL, DIMENSION(knon), INTENT(IN) :: phiT_b9, phiQ_b9 738 REAL, DIMENSION(knon), INTENT(IN) :: phiT_x9, phiT_w9 739 REAL, DIMENSION(knon), INTENT(IN) :: phiQ_x9, phiQ_w9 740 741 !! Local variables 742 INTEGER :: j 743 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region 744 REAL, DIMENSION(knon) :: AcoefT_b, AcoefQ_b ! mean values of AcoefT and AcoefQ 745 REAL :: zzt, zzq, zzqsat 746 REAL :: zdelta, zcvm5, zcor, qsat 747 REAL, DIMENSION(knon) :: qsat_w, qsat_x 748 REAL, DIMENSION(knon) :: dqsatdT_w, dqsatdT_x 749 REAL, DIMENSION(knon) :: qsat_bs ! qsat(Ts_b) 750 REAL, DIMENSION(knon) :: qsat01, dqsatdT01 751 REAL, DIMENSION(knon) :: Ts_x, Ts_w, qs_x, qs_w 752 REAL, DIMENSION(knon) :: T1_x, T1_w, q1_x, q1_w 753 REAL, DIMENSION(knon) :: Rn_x, Rn_w 754 REAL, DIMENSION(knon) :: phiQ0_x, phiQ0_w 755 REAL, DIMENSION(knon) :: Ta, qa 756 REAL, DIMENSION(knon) :: qsatsrf_w, qsatsrf_x, qsatsrf_b 757 REAL, DIMENSION(knon) :: qsurf_w, qsurf_x 758 REAL :: dphiT, dphiQ 759 REAL :: dqsatsrf1 760 REAL :: phiT_w1, phiT_w2 761 REAL :: phiT_x1, phiT_x2 762 REAL :: phiQ_w1, phiQ_w2, phiQ_w3 763 REAL :: phiQ_x1, phiQ_x2, phiQ_x3 764 REAL :: phiT_b1, phiQ_b1 765 REAL :: Kech_Q_sw1, Kech_Q_sx1 766 REAL :: evap_pot 767 768 !---------------------------------------------------------------------------- 769 ! Equations to be checked: 770 ! ----------------------- 771 ! Input : Ts0_b, dTs0, Ts_b, dTs, qsatsrf_b, dqsatsrf, 772 ! phiT_b, phiQ_b, phiT_w, phiT_x, phiQ_w, phiQ_x, 773 774 ! AcoefT, AcoefQ, AcoefT_w, AcoefQ_w, AcoefT_x, AcoefQ_x, 775 ! BcoefT, BcoefQ, BcoefT_w, BcoefQ_w, BcoefT_x, BcoefQ_x 776 777 ! C_p T1_w = AcoefT_w + BcoefT_w phiT_w Delta t C_p T1_x = AcoefT_x + BcoefT_x phiT_x Delta t 778 ! q1_w = AQ_w + BQ_w phiQ_w Delta t q1_x = AQ_x + BQ_x phiQ_x Delta t 779 ! qsatsrf_w = beta qsat(Ts_w) qsatsrf_x = beta qsat(Ts_x) 780 ! qsurf_w = (1-beta) q1_w + qsatsrf_w qsurf_x = (1-beta) q1_x + qsatsrf_x 781 ! phiT_w = Kech_h_w C_p ( T1_w - Ts_w) phiT_x = Kech_h_x C_p ( T1_x - Ts_x) 782 ! phiT_w = Kech_T_pw ( AcoefT_w - C_p Ts_w) phiT_x = Kech_T_px ( AcoefT_x - C_p Ts_x) 783 ! phiq_w = Kech_h_w ( beta q1_w - qsatsrf_w) phiq_x = Kech_h_x ( beta q1_x - qsatsrf_x)) 784 ! phiq_w = Kech_Q_sw (beta AQ_w -qsatsrf_w) phiq_x = Kech_Q_sx (beta AQ_x -qsatsrf_x) 785 ! phiq_w = Kech_h_w (q1_w - qsurf_w) phiq_x = Kech_h_x (q1_x - qsurf_x) 786 ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x 787 ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x 788 ! Ts_b = sigw Ts_w + sigx Ts_x dTs = Ts_w - Ts_x 789 ! qsatsrf_b = sigw qsatsrf_w + sigx qsatsrf_x 790 ! C_p Ta = AcoefT + BcoefT phiT_b Delta t 791 ! qa = AcoefQ + BcoefQ phiQ_b Delta t 792 ! phiT_b = Kech_h C_p (Ta - Ts_b) 793 ! phiQ_b = beta Kech_h (qa - qsatsrf_b) 794 ! dTs = sqrt(tau)/I (dphit + L_v dphiq + dR) 795 796 !---------------------------------------------------------------------------- 797 798 !! 799 sigx(:) = 1. - sigw(:) 800 AcoefT_b(1:knon) = AcoefT_x(1:knon) + sigw(1:knon) * dd_AT(1:knon) 801 AcoefQ_b(1:knon) = AQ_x(1:knon) + sigw(1:knon) * dd_AQ(1:knon) 802 803 ! Compute the three qsat and dqsatdTs 804 ! --------------------------------------------- 805 !! C_p(1:knon) = RCpd 806 !! L_v(1:knon) = RLvtt 807 IF (prt_level >=10) THEN 808 print *, ' AAAA wx_pbl_check, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:) 815 809 ENDIF ! (prt_level >=10 ) 816 810 817 DO j = 1, knon 818 zdelta = MAX(0.,SIGN(1.,RTT-Ts0_b9(j))) 819 zcvm5 = R5LES*(1.-zdelta) + R5IES*zdelta 820 qsat = R2ES*FOEEW(Ts0_b9(j),zdelta)/ypaprs(j,1) 821 qsat = MIN(0.5,qsat) 822 zcor = 1./(1.-RETV*qsat) 823 qsat01(j) = fqsat*qsat*zcor 824 !! dqsatdT0(j) = FOEDE(Ts0_b(j),zdelta,zcvm5,qsat0(j),zcor)/RLVTT ! jyg 20210116 825 !! dqsatdT0(j) = (RLvtt*(1.-zdelta)+RLSTT*zdelta)*qsat0(j)/(Rv*Ts0_b(j)*Ts0_b(j)) 826 dqsatdT01(j) = fqsat*FOEDE(Ts0_b9(j),zdelta,zcvm5,qsat01(j),zcor) 827 ENDDO 828 829 !-------------------------------------------------------------------------------------------------- 830 IF (prt_level >=10 ) THEN 831 832 DO j = 1, knon 833 834 print *,'wx_pbl_check: Kech_h, Kech_q ', Kech_h(j), Kech_q(j) 835 836 Ta(j) = (AcoefT(j) + BcoefT(j)*phiT_b9(j)*dtime)/C_p(j) 837 qa(j) = AcoefQ(j) + BcoefQ(j)*phiQ_b9(j)*dtime 838 print *, 'wx_pbl_check: j, Ta, qa ', Ta(j), qa(j) 839 840 qsat_bs(j) = qsat01(j) + dqsatdT01(j)*(Ts_b9(j)-Ts0_b9(j)) 841 842 print *,'wx_pbl_check: qsat01, qsat_bs ', j,qsat01(j), qsat_bs(j) 843 844 Ts_x(j) = Ts_b9(j) - sigw(j)*dTs9(j) 845 Ts_w(j) = Ts_b9(j) + sigx(j)*dTs9(j) 846 print *, 'wx_pbl_check: j, Ts_b9, Ts_w, Ts_x ', j, Ts_b9(j), Ts_w(j), Ts_x(j) 847 848 qsat_x(j) = qsat0_x(j) + dqsatdT0_x(j)*(Ts_x(j)-Ts0_x(j)) 849 qsat_w(j) = qsat0_w(j) + dqsatdT0_w(j)*(Ts_w(j)-Ts0_w(j)) 850 851 print *,'wx_pbl_check: qsat0_w, qsat0_x, qsat_w, qsat_x ', qsat0_w(j), qsat0_x(j), qsat_w(j), qsat_x(j) 852 853 T1_x(j) = (AcoefT_x(j) + BcoefT_x(j)*phiT_x9(j)*dtime) / C_p(j) 854 T1_w(j) = (AcoefT_w(j) + BcoefT_w(j)*phiT_w9(j)*dtime) / C_p(j) 855 print *, 'wx_pbl_check: j, T1_w, T1_x ', j, T1_w(j), T1_x(j) 856 857 q1_x(j) = AQ_x(j) + BQ_x(j)*phiQ_x9(j)*dtime 858 q1_w(j) = AQ_w(j) + BQ_w(j)*phiQ_w9(j)*dtime 859 print *, 'wx_pbl_check: j, q1_w, q1_x ', j, q1_w(j), q1_x(j) 860 861 qsatsrf_x(j) = beta(j)*qsat_x(j) 862 qsatsrf_w(j) = beta(j)*qsat_w(j) 863 qsatsrf_b(j) = sigw(j)*qsatsrf_w(j) + sigx(j)*qsatsrf_x(j) 864 865 dqsatsrf1 = qsatsrf_w(j) - qsatsrf_x(j) 866 print *, 'wx_pbl_check: j, qsatsrf_w, qsatsrf_x, dqsatsrf1, dqsatsrf9 ', & 867 qsatsrf_w(j), qsatsrf_x(j), dqsatsrf1, dqsatsrf9(j) 868 869 qsurf_x(j) = (1-beta(j))*q1_x(j) + qsatsrf_x(j) 870 qsurf_w(j) = (1-beta(j))*q1_w(j) + qsatsrf_w(j) 871 print *, 'wx_pbl_check: j, qsurf_w, qsurf_x ', j, qsurf_w(j), qsurf_x(j) 872 873 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 874 ! Test qsat01 = qsat0 et dqsatdT01 = dqsatdT0 875 !------------------------------------------------------------------------------------------------------ 876 print *, 'wx_pbl_check: j, qsat01(j), qsat0(j) ', j, qsat01(j), qsat0(j) 877 print *, 'wx_pbl_check: j, dqsatdT01(j), dqsatdT0(j) ', j, dqsatdT01(j), dqsatdT0(j) 878 879 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 880 ! Test Kexh_Q_sw = Kech_q_w/(1.-beta*Kech_q_w*BcoefQ) Kexh_Q_sx = Kech_q_x/(1.-beta*Kech_q_x*BcoefQ) 881 !------------------------------------------------------------------------------------------------------ 882 Kech_Q_sx1 = Kech_q_x(j)/(1.-beta(j)*Kech_q_x(j)*BQ_x(j)*dtime) 883 Kech_Q_sw1 = Kech_q_w(j)/(1.-beta(j)*Kech_q_w(j)*BQ_w(j)*dtime) 884 print *, 'wx_pbl_check: j, Kech_Q_sx1, Kech_Q_sx(j)', j, Kech_Q_sx1, Kech_Q_sx(j) 885 print *, 'wx_pbl_check: j, Kech_Q_sw1, Kech_Q_sw(j)', j, Kech_Q_sw1, Kech_Q_sw(j) 886 887 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 888 ! Test phiT_w = Kech_h_w*C_p(j)*(T1_w(j)-Ts_w(j)) phiT_x = Kech_h_x*C_p(j)*(T1_x(j)-Ts_x(j)) 889 !----------------------------------------------------- 890 phiT_x1 = Kech_h_x(j)*C_p(j)*(T1_x(j)-Ts_x(j)) 891 phiT_w1 = Kech_h_w(j)*C_p(j)*(T1_w(j)-Ts_w(j)) 892 893 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 894 ! Test phiT_w = Kech_T_pw*(AcoefT_w(j)-C_p(j)*Ts_w(j)) phiT_x = Kech_T_px*(AcoefT_x(j)-C_p(j)*Ts_x(j)) 895 !----------------------------------------------------- 896 phiT_x2 = Kech_T_px(j)*(AcoefT_x(j)-C_p(j)*Ts_x(j)) 897 phiT_w2 = Kech_T_pw(j)*(AcoefT_w(j)-C_p(j)*Ts_w(j)) 898 print *, 'wx_pbl_check: j, phiT_w1, phiT_w2, phiT_w9 ', j, phiT_w1, phiT_w2, phiT_w9(j) 899 print *, 'wx_pbl_check: j, phiT_x1, phiT_x2, phiT_x9 ', j, phiT_x1, phiT_x2, phiT_x9(j) 900 901 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 902 ! Test phiq_w = Kech_q_w ( beta q1_w - qsatsrf_w) phiq_x = Kech_q_x ( beta q1_x - qsatsrf_x)) 903 !-------------------------------------------------------------- 904 phiq_x1 = Kech_q_x(j)*( beta(j)*q1_x(j) - qsatsrf_x(j)) 905 phiq_w1 = Kech_q_w(j)*( beta(j)*q1_w(j) - qsatsrf_w(j)) 906 907 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 908 ! Test phiq_w = Kech_Q_sw (beta AQ_w -qsatsrf_w) phiq_x = Kech_Q_sx (beta AQ_x -qsatsrf_x) 909 !-------------------------------------------------------------- 910 phiq_x2 = Kech_Q_sx(j)*(beta(j)*AQ_x(j) -qsatsrf_x(j)) 911 phiq_w2 = Kech_Q_sw(j)*(beta(j)*AQ_w(j) -qsatsrf_w(j)) 912 913 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 914 ! Test phiq_w = Kech_q_w ( q1_w - qsurf_w) phiq_x = Kech_q_x ( q1_x - qsurf_x)) 915 !-------------------------------------------------------------- 916 phiq_x3 = Kech_q_x(j)*( q1_x(j) - qsurf_x(j)) 917 phiq_w3 = Kech_q_w(j)*( q1_w(j) - qsurf_w(j)) 918 print *, 'wx_pbl_check: j, phiQ_w1, phiQ_w2, phiQ_w3, phiQ_w9 ', j, phiQ_w1, phiQ_w2, phiQ_w3, phiQ_w9(j) 919 print *, 'wx_pbl_check: j, phiQ_x1, phiQ_x2, phiQ_x3, phiQ_x9 ', j, phiQ_x1, phiQ_x2, phiQ_x3, phiQ_x9(j) 920 921 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 922 ! Test phiT_b = Kech_h C_p (Ta - Ts_b) 923 !-------------------------------------------------------------- 924 phiT_b1 = Kech_h(j)*C_p(j)*(Ta(j) - Ts_b9(j)) 925 print *, 'wx_pbl_check: j, phiT_b1, PhiT_b9 ', j, phiT_b1, PhiT_b9(j) 926 927 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 928 ! Test phiQ_b = beta Kech_q (qa - qsat_bs) 929 !-------------------------------------------------------------- 930 evap_pot = Kech_q(j)*(qa(j) - qsat_bs(j)) 931 phiQ_b1 = beta(j)*Kech_q(j)*(qa(j) - qsat_bs(j)) 932 print *, 'wx_pbl_check: j, beta, evap_pot, phiQ_b1, PhiQ_b9 ', j, beta(j), evap_pot, phiQ_b1, PhiQ_b9(j) 933 934 935 ENDDO ! j = 1, knon 936 937 ENDIF ! (prt_level >=10 ) 938 !-------------------------------------------------------------------------------------------------- 939 940 941 942 END SUBROUTINE wx_pbl_check 943 944 SUBROUTINE wx_pbl_dts_check( knon, dtime, ypplay, ypaprs, & 945 sigw, beta, iflag_split, & 946 Ts0_b9, dTs09, & 947 qs_b9, Ts_b9, & ! yqsurf, Tsurf_new 948 dqsatsrf9, dTs9, delta_qsurf9, & 949 AcoefT_x, AcoefT_w, & 950 BcoefT_x, BcoefT_w, & 951 AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, & 952 AcoefT, AcoefQ, BcoefT, BcoefQ, & 953 HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, & 954 phiT0_b9, dphiT09, phiQ0_b9, dphiQ09, Rn0_b9, dRn09, & 955 g_T, g_Q, & 956 Gamma_phiT, Gamma_phiQ, & 957 dTs_ins, dqsatsrf_ins, & 958 phiT_b9, phiQ_b9, & 959 phiT_x9, phiT_w9, & 960 phiQ_x9, phiQ_w9 & 961 ) 811 DO j = 1, knon 812 zdelta = MAX(0., SIGN(1., RTT - Ts0_b9(j))) 813 zcvm5 = R5LES * (1. - zdelta) + R5IES * zdelta 814 qsat = R2ES * FOEEW(Ts0_b9(j), zdelta) / ypaprs(j, 1) 815 qsat = MIN(0.5, qsat) 816 zcor = 1. / (1. - RETV * qsat) 817 qsat01(j) = fqsat * qsat * zcor 818 !! dqsatdT0(j) = FOEDE(Ts0_b(j),zdelta,zcvm5,qsat0(j),zcor)/RLVTT ! jyg 20210116 819 !! dqsatdT0(j) = (RLvtt*(1.-zdelta)+RLSTT*zdelta)*qsat0(j)/(Rv*Ts0_b(j)*Ts0_b(j)) 820 dqsatdT01(j) = fqsat * FOEDE(Ts0_b9(j), zdelta, zcvm5, qsat01(j), zcor) 821 ENDDO 822 823 !-------------------------------------------------------------------------------------------------- 824 IF (prt_level >=10) THEN 825 826 DO j = 1, knon 827 828 print *, 'wx_pbl_check: Kech_h, Kech_q ', Kech_h(j), Kech_q(j) 829 830 Ta(j) = (AcoefT(j) + BcoefT(j) * phiT_b9(j) * dtime) / C_p(j) 831 qa(j) = AcoefQ(j) + BcoefQ(j) * phiQ_b9(j) * dtime 832 print *, 'wx_pbl_check: j, Ta, qa ', Ta(j), qa(j) 833 834 qsat_bs(j) = qsat01(j) + dqsatdT01(j) * (Ts_b9(j) - Ts0_b9(j)) 835 836 print *, 'wx_pbl_check: qsat01, qsat_bs ', j, qsat01(j), qsat_bs(j) 837 838 Ts_x(j) = Ts_b9(j) - sigw(j) * dTs9(j) 839 Ts_w(j) = Ts_b9(j) + sigx(j) * dTs9(j) 840 print *, 'wx_pbl_check: j, Ts_b9, Ts_w, Ts_x ', j, Ts_b9(j), Ts_w(j), Ts_x(j) 841 842 qsat_x(j) = qsat0_x(j) + dqsatdT0_x(j) * (Ts_x(j) - Ts0_x(j)) 843 qsat_w(j) = qsat0_w(j) + dqsatdT0_w(j) * (Ts_w(j) - Ts0_w(j)) 844 845 print *, 'wx_pbl_check: qsat0_w, qsat0_x, qsat_w, qsat_x ', qsat0_w(j), qsat0_x(j), qsat_w(j), qsat_x(j) 846 847 T1_x(j) = (AcoefT_x(j) + BcoefT_x(j) * phiT_x9(j) * dtime) / C_p(j) 848 T1_w(j) = (AcoefT_w(j) + BcoefT_w(j) * phiT_w9(j) * dtime) / C_p(j) 849 print *, 'wx_pbl_check: j, T1_w, T1_x ', j, T1_w(j), T1_x(j) 850 851 q1_x(j) = AQ_x(j) + BQ_x(j) * phiQ_x9(j) * dtime 852 q1_w(j) = AQ_w(j) + BQ_w(j) * phiQ_w9(j) * dtime 853 print *, 'wx_pbl_check: j, q1_w, q1_x ', j, q1_w(j), q1_x(j) 854 855 qsatsrf_x(j) = beta(j) * qsat_x(j) 856 qsatsrf_w(j) = beta(j) * qsat_w(j) 857 qsatsrf_b(j) = sigw(j) * qsatsrf_w(j) + sigx(j) * qsatsrf_x(j) 858 859 dqsatsrf1 = qsatsrf_w(j) - qsatsrf_x(j) 860 print *, 'wx_pbl_check: j, qsatsrf_w, qsatsrf_x, dqsatsrf1, dqsatsrf9 ', & 861 qsatsrf_w(j), qsatsrf_x(j), dqsatsrf1, dqsatsrf9(j) 862 863 qsurf_x(j) = (1 - beta(j)) * q1_x(j) + qsatsrf_x(j) 864 qsurf_w(j) = (1 - beta(j)) * q1_w(j) + qsatsrf_w(j) 865 print *, 'wx_pbl_check: j, qsurf_w, qsurf_x ', j, qsurf_w(j), qsurf_x(j) 866 867 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 868 ! Test qsat01 = qsat0 et dqsatdT01 = dqsatdT0 869 !------------------------------------------------------------------------------------------------------ 870 print *, 'wx_pbl_check: j, qsat01(j), qsat0(j) ', j, qsat01(j), qsat0(j) 871 print *, 'wx_pbl_check: j, dqsatdT01(j), dqsatdT0(j) ', j, dqsatdT01(j), dqsatdT0(j) 872 873 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 874 ! Test Kexh_Q_sw = Kech_q_w/(1.-beta*Kech_q_w*BcoefQ) Kexh_Q_sx = Kech_q_x/(1.-beta*Kech_q_x*BcoefQ) 875 !------------------------------------------------------------------------------------------------------ 876 Kech_Q_sx1 = Kech_q_x(j) / (1. - beta(j) * Kech_q_x(j) * BQ_x(j) * dtime) 877 Kech_Q_sw1 = Kech_q_w(j) / (1. - beta(j) * Kech_q_w(j) * BQ_w(j) * dtime) 878 print *, 'wx_pbl_check: j, Kech_Q_sx1, Kech_Q_sx(j)', j, Kech_Q_sx1, Kech_Q_sx(j) 879 print *, 'wx_pbl_check: j, Kech_Q_sw1, Kech_Q_sw(j)', j, Kech_Q_sw1, Kech_Q_sw(j) 880 881 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 882 ! Test phiT_w = Kech_h_w*C_p(j)*(T1_w(j)-Ts_w(j)) phiT_x = Kech_h_x*C_p(j)*(T1_x(j)-Ts_x(j)) 883 !----------------------------------------------------- 884 phiT_x1 = Kech_h_x(j) * C_p(j) * (T1_x(j) - Ts_x(j)) 885 phiT_w1 = Kech_h_w(j) * C_p(j) * (T1_w(j) - Ts_w(j)) 886 887 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 888 ! Test phiT_w = Kech_T_pw*(AcoefT_w(j)-C_p(j)*Ts_w(j)) phiT_x = Kech_T_px*(AcoefT_x(j)-C_p(j)*Ts_x(j)) 889 !----------------------------------------------------- 890 phiT_x2 = Kech_T_px(j) * (AcoefT_x(j) - C_p(j) * Ts_x(j)) 891 phiT_w2 = Kech_T_pw(j) * (AcoefT_w(j) - C_p(j) * Ts_w(j)) 892 print *, 'wx_pbl_check: j, phiT_w1, phiT_w2, phiT_w9 ', j, phiT_w1, phiT_w2, phiT_w9(j) 893 print *, 'wx_pbl_check: j, phiT_x1, phiT_x2, phiT_x9 ', j, phiT_x1, phiT_x2, phiT_x9(j) 894 895 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 896 ! Test phiq_w = Kech_q_w ( beta q1_w - qsatsrf_w) phiq_x = Kech_q_x ( beta q1_x - qsatsrf_x)) 897 !-------------------------------------------------------------- 898 phiq_x1 = Kech_q_x(j) * (beta(j) * q1_x(j) - qsatsrf_x(j)) 899 phiq_w1 = Kech_q_w(j) * (beta(j) * q1_w(j) - qsatsrf_w(j)) 900 901 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 902 ! Test phiq_w = Kech_Q_sw (beta AQ_w -qsatsrf_w) phiq_x = Kech_Q_sx (beta AQ_x -qsatsrf_x) 903 !-------------------------------------------------------------- 904 phiq_x2 = Kech_Q_sx(j) * (beta(j) * AQ_x(j) - qsatsrf_x(j)) 905 phiq_w2 = Kech_Q_sw(j) * (beta(j) * AQ_w(j) - qsatsrf_w(j)) 906 907 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 908 ! Test phiq_w = Kech_q_w ( q1_w - qsurf_w) phiq_x = Kech_q_x ( q1_x - qsurf_x)) 909 !-------------------------------------------------------------- 910 phiq_x3 = Kech_q_x(j) * (q1_x(j) - qsurf_x(j)) 911 phiq_w3 = Kech_q_w(j) * (q1_w(j) - qsurf_w(j)) 912 print *, 'wx_pbl_check: j, phiQ_w1, phiQ_w2, phiQ_w3, phiQ_w9 ', j, phiQ_w1, phiQ_w2, phiQ_w3, phiQ_w9(j) 913 print *, 'wx_pbl_check: j, phiQ_x1, phiQ_x2, phiQ_x3, phiQ_x9 ', j, phiQ_x1, phiQ_x2, phiQ_x3, phiQ_x9(j) 914 915 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 916 ! Test phiT_b = Kech_h C_p (Ta - Ts_b) 917 !-------------------------------------------------------------- 918 phiT_b1 = Kech_h(j) * C_p(j) * (Ta(j) - Ts_b9(j)) 919 print *, 'wx_pbl_check: j, phiT_b1, PhiT_b9 ', j, phiT_b1, PhiT_b9(j) 920 921 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 922 ! Test phiQ_b = beta Kech_q (qa - qsat_bs) 923 !-------------------------------------------------------------- 924 evap_pot = Kech_q(j) * (qa(j) - qsat_bs(j)) 925 phiQ_b1 = beta(j) * Kech_q(j) * (qa(j) - qsat_bs(j)) 926 print *, 'wx_pbl_check: j, beta, evap_pot, phiQ_b1, PhiQ_b9 ', j, beta(j), evap_pot, phiQ_b1, PhiQ_b9(j) 927 928 ENDDO ! j = 1, knon 929 930 ENDIF ! (prt_level >=10 ) 931 !-------------------------------------------------------------------------------------------------- 932 933 END SUBROUTINE wx_pbl_check 934 935 SUBROUTINE wx_pbl_dts_check(knon, dtime, ypplay, ypaprs, & 936 sigw, beta, iflag_split, & 937 Ts0_b9, dTs09, & 938 qs_b9, Ts_b9, & ! yqsurf, Tsurf_new 939 dqsatsrf9, dTs9, delta_qsurf9, & 940 AcoefT_x, AcoefT_w, & 941 BcoefT_x, BcoefT_w, & 942 AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, & 943 AcoefT, AcoefQ, BcoefT, BcoefQ, & 944 HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, & 945 phiT0_b9, dphiT09, phiQ0_b9, dphiQ09, Rn0_b9, dRn09, & 946 g_T, g_Q, & 947 Gamma_phiT, Gamma_phiQ, & 948 dTs_ins, dqsatsrf_ins, & 949 phiT_b9, phiQ_b9, & 950 phiT_x9, phiT_w9, & 951 phiQ_x9, phiQ_w9 & 952 ) 962 953 963 954 USE wx_pbl_var_mod 964 955 965 USE lmdz_print_control, ONLY: prt_level, lunout966 USE lmdz_ YOETHF956 USE lmdz_print_control, ONLY: prt_level, lunout 957 USE lmdz_yoethf 967 958 USE lmdz_fcttre, ONLY: foeew, foede, qsats, qsatl, dqsats, dqsatl, thermcep 968 969 INCLUDE "YOMCST.h" 970 971 INTEGER, INTENT(IN) :: knon ! number of grid cells 972 REAL, INTENT(IN) :: dtime ! time step size (s)973 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa)974 REAL, DIMENSION(knon, klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa)975 REAL, DIMENSION(knon ), INTENT(IN) :: sigw ! cold pools fractional area976 REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor977 INTEGER, INTENT(IN) :: iflag_split978 REAL, DIMENSION(knon), INTENT(IN) :: Ts0_b9, dTs09979 REAL, DIMENSION(knon), INTENT(IN) :: qs_b9, Ts_b9 ! yqsurf, Tsurf_new980 REAL, DIMENSION(knon), INTENT(IN) :: dTs9, dqsatsrf9981 REAL, DIMENSION(knon), INTENT(IN) :: delta_qsurf9982 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w983 REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w984 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0985 986 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT, AcoefQ, BcoefT, BcoefQ 987 REAL, DIMENSION(knon), INTENT(IN) :: HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn988 REAL, DIMENSION(knon), INTENT(IN) :: phiT0_b9, dphiT09, phiQ0_b9, dphiQ09, Rn0_b9, dRn09989 REAL, DIMENSION(knon), INTENT(IN) :: g_T, g_Q990 REAL, DIMENSION(knon), INTENT(IN) :: Gamma_phiT, Gamma_phiQ991 REAL, DIMENSION(knon), INTENT(IN) :: dTs_ins, dqsatsrf_ins992 REAL, DIMENSION(knon), INTENT(IN) :: phiT_b9, phiQ_b9993 REAL, DIMENSION(knon), INTENT(IN) :: phiT_x9, phiT_w9994 REAL, DIMENSION(knon), INTENT(IN) :: phiQ_x9, phiQ_w9995 996 !! Local variables 997 INTEGER :: j998 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region999 REAL, DIMENSION(knon) :: AcoefT_b, AcoefQ_b ! mean values of AcoefT and AcoefQ1000 REAL :: zzt, zzq, zzqsat1001 REAL :: zdelta, zcvm5, zcor,qsat1002 REAL , DIMENSION(knon) :: qsat_w, qsat_x1003 REAL, DIMENSION(knon) :: Ts_x, Ts_w, qs_x, qs_w1004 REAL, DIMENSION(knon) :: T1_x, T1_w, q1_x, q1_w1005 REAL, DIMENSION(knon) :: Rn_x, Rn_w1006 REAL, DIMENSION(knon) :: Rn_b, dRn1007 REAL, DIMENSION(knon) :: phiQ0_x, phiQ0_w1008 REAL, DIMENSION(knon) :: Ta, qa1009 REAL, DIMENSION(knon) :: err_phiT_w, err_phiT_x1010 REAL, DIMENSION(knon) :: err_phiq_w, err_phiq_x1011 REAL, DIMENSION(knon) :: err_phiT_b1012 REAL, DIMENSION(knon) :: err_phiQ_b1013 REAL, DIMENSION(knon) :: err2_phiT_b1014 REAL :: T1A_x, T1A_w, q1A_x, q1A_w1015 REAL :: qsatsrf_w, qsatsrf_x, qsatsrfb, qsbA1016 REAL :: dphiT, dphiQ1017 REAL :: dphiT_H, dphiQ_H1018 REAL :: phiQ_pot1019 REAL :: phiQ_w_m_phiQ0_w1020 REAL :: phiQ_x_m_phiQ0_x1021 REAL :: dphiQ_m_dphiQ01022 REAL :: dphiT_m_dphiT01023 REAL :: dRN_m_dRn01024 REAL :: phiTb_m_phiT0b1025 1026 !---------------------------------------------------------------------------- 1027 ! Equations to be checked: 1028 ! ----------------------- 1029 ! Input : Ts0_b, dTs0, Ts_b, dTs, qsatsrf_b, dqsatsrf, 1030 ! phiT_b, phiQ_b, phiT_w, phiT_x, phiQ_w, phiQ_x, 1031 1032 ! AcoefT, AcoefQ, AcoefT_w, AcoefQ_w, AcoefT_x, AcoefQ_x, 1033 ! BcoefT, BcoefQ, BcoefT_w, BcoefQ_w, BcoefT_x, BcoefQ_x 1034 1035 ! Ts_w = Ts_b + sigx dTs Ts_x = Ts_b - sigw dTs 1036 ! T1_w = AcoefT_w + BcoefT_w phiT_w Delta t T1_x = AcoefT_x + BcoefT_x phiT_x Delta t 1037 ! q1_w = AcoefQ_w + BcoefQ_w phiQ_w Delta t q1_x = AcoefQ_x + BcoefQ_x phiQ_x Delta t1038 ! phiT_w = Kech_h_w ( T1_w - Ts_w) phiT_x = Kech_h_x ( T1_x - Ts_x) 1039 ! phiq_w = beta Kech_h_w ( q1_w - qsat(Ts_w)) phiq_x = beta Kech_h_x ( q1_x - qsat(Ts_x))1040 ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x 1041 ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x1042 ! Ts_b = sigw Ts_w + sigx Ts_x dTs = Ts_w - Ts_x1043 ! Ta = AcoefT + BcoefT phiT_b Delta t 1044 ! qa = AcoefQ + BcoefQ phiQ_b Delta t1045 ! phiT_b = Kech_h (Ta - Ts_b) 1046 ! phiQ_b = beta Kech_h (qa - qsat(Ts_b))1047 ! dTs = sqrt(tau)/I (dphit + L_v dphiq + dR)1048 1049 !---------------------------------------------------------------------------- 1050 1051 !! 1052 sigx(:) = 1.-sigw(:)1053 AcoefT_b(1:knon) = AcoefT_x(1:knon) + sigw(1:knon)*dd_AT(1:knon)1054 AcoefQ_b(1:knon) = AQ_x(1:knon) + sigw(1:knon)*dd_AQ(1:knon)1055 1056 IF (prt_level >=10 ) THEN 1057 print *,'->wx_pbl_dts_check, HTphiT_b, HTphiQ_b, HTRn_b ', &1058 HTphiT_b, HTphiQ_b, HTRn_b1059 print *,'->wx_pbl_dts_check, dd_HTphiT, dd_HTphiQ, dd_HTRn ', &1060 dd_HTphiT, dd_HTphiQ, dd_HTRn1061 ENDIF ! (prt_level >=10 )1062 1063 ! Compute the three qsat and dqsatdTs 1064 ! --------------------------------------------- 1065 !! print *,' AAAA wx_pbl_dts_check, C_p(j), qsat0(j), Ts0(j) : ', & 1066 !! (C_p(j), qsat0(j), Ts0(j), j = 1,knon) 1067 1068 1069 !-------------------------------------------------------------------------------------------------- 1070 IF (prt_level >=10 ) THEN1071 1072 DO j = 1, knon 1073 Ts_x(j) = Ts_b9(j) - sigw(j)*dTs9(j)1074 Ts_w(j) = Ts_b9(j) + sigx(j)*dTs9(j)1075 print *, 'wx_pbl_dts_check: j, Ts_b9, Ts_w, Ts_x ', j, Ts_b9(j), Ts_w(j), Ts_x(j)1076 1077 qsat_x(j) = qsat0_x(j) + dqsatdT0_x(j)*(Ts_x(j)-Ts0_x(j)) 1078 qsat_w(j) = qsat0_w(j) + dqsatdT0_w(j)*(Ts_w(j)-Ts0_w(j))1079 1080 T1_x(j) = (AcoefT_x(j) + BcoefT_x(j)*phiT_x9(j)*dtime) / C_p(j) 1081 T1_w(j) = (AcoefT_w(j) + BcoefT_w(j)*phiT_w9(j)*dtime) / C_p(j)1082 print *, 'wx_pbl_dts_check: j, T1_w, T1_x ', j, T1_w(j), T1_x(j)1083 1084 q1_x(j) = AQ_x(j) + BQ_x(j)*phiQ_x9(j)*dtime 1085 q1_w(j) = AQ_w(j) + BQ_w(j)*phiQ_w9(j)*dtime1086 print *, 'wx_pbl_dts_check: j, q1_w, q1_x ', j, q1_w(j), q1_x(j)1087 1088 Rn_x(j) = eps_1*Rsigma*T1_x(j)**4 - Rsigma*Ts_x(j)**4 1089 Rn_w(j) = eps_1*Rsigma*T1_w(j)**4 - Rsigma*Ts_w(j)**41090 Rn_b(j) = sigw(j)*Rn_w(j) + sigx(j)*Rn_x(j)1091 dRn(j) = dRn09(j) - ( HTRn_b(j) &1092 +(sigx(j)-sigw(j))*dd_HTRn(j) &1093 -sigw(j)*sigx(j)*dd_HTRn(j)*dd_HTphiT(j)/HTphiT_b(j) &1094 )*(dTs9(j)-dTs09(j)) &1095 + dd_HTRn(j)/HTphiT_b(j)*(phiT_b9(j)-phiT0_b9(j))1096 1097 print *,'wx_pbl_dts_check, dphiT, L_v*dphiQ, dRn, dTs ', & 1098 phiT_w9(j)-phiT_x9(j), L_v(j)*(phiQ_w9(j)-phiQ_x9(j)), dRn(j), dTs9(j)1099 1100 phiQ0_x(j) = PhiQ0_b9(j) - sigw(j)*dphiQ09(j) 1101 phiQ0_w(j) = PhiQ0_b9(j) + sigx(j)*dphiQ09(j)1102 1103 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1104 ! Test phiQ_w-phiQ0_w = -beta*Kech_Q_sw*dqsatdT_w*(Ts_w-Ts0_w) 1105 !-------------------------------------------------------------- 1106 print *,'wx_pbl_dts_check: beta(j), Kech_Q_sw(j), dqsatdT0_w(j), Ts_w(j), Ts0_w(j) ', &1107 beta(j), Kech_Q_sw(j), dqsatdT0_w(j), Ts_w(j), Ts0_w(j)1108 phiQ_w_m_phiQ0_w = -beta(j)*Kech_Q_sw(j)*dqsatdT0_w(j)*(Ts_w(j)-Ts0_w(j))1109 print *,'wx_pbl_dts_check: j, phiQ_w9-phiQ0_w, phiQ_w_m_phiQ0_w ', &1110 j, phiQ_w9(j)-phiQ0_w(j), phiQ_w_m_phiQ0_w1111 1112 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1113 ! Test phiQ_x-phiQ0_x = -beta*Kech_Q_sx*dqsatdT_x*(Ts_x-Ts0_x) 1114 !-------------------------------------------------------------- 1115 phiQ_x_m_phiQ0_x = -beta(j)*Kech_Q_sx(j)*dqsatdT0_x(j)*(Ts_x(j)-Ts0_x(j))1116 print *,'wx_pbl_dts_check: j, phiQ_x9-phiQ0_x, phiQ_x_m_phiQ0_x ', &1117 j, phiQ_x9(j)-phiQ0_x(j), phiQ_x_m_phiQ0_x1118 1119 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1120 ! Test dphiT-dphiT0 = -(HTphiT_b+(sigx-sigw)*dd_HTphiT)*(dTs-dTs0) - dd_HTphiT*(Ts_b-Ts0_b) 1121 !------------------------------------------------------------------------------------------- 1122 dphiT = phiT_w9(j) - phiT_x9(j)1123 dphiT_m_dphiT0 = -(HTphiT_b(j)+(sigx(j)-sigw(j))*dd_HTphiT(j))*(dTs9(j)-dTs09(j)) &1124 - dd_HTphiT(j)*(Ts_b9(j)-Ts0_b9(j))1125 print *,'wx_pbl_dts_check: j, dphiT-dphiT09, dphiT_m_dphiT0 ',j, dphiT-dphiT09(j), dphiT_m_dphiT01126 1127 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1128 ! Test dphiQ-dphiQ0 = -(HTphiQ_b+(sigx-sigw)*dd_HTphiQ)*(dTs-dTs0) - dd_HTphiQ*(Ts_b-Ts0_b) 1129 !------------------------------------------------------------------------------------------- 1130 dphiQ = phiQ_w9(j) - phiQ_x9(j)1131 dphiQ_m_dphiQ0 = -(HTphiQ_b(j)+(sigx(j)-sigw(j))*dd_HTphiQ(j))*(dTs9(j)-dTs09(j)) &1132 - dd_HTphiQ(j)*(Ts_b9(j)-Ts0_b9(j))1133 print *,'wx_pbl_dts_check: j, dphiQ-dphiQ09, dphiQ_m_dphiQ0 ',j, dphiQ-dphiQ09(j), dphiQ_m_dphiQ01134 1135 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1136 ! Test dRn-dRn0 = -(HTRn_b+(sigx-sigw)*dd_HTRn)*(dTs-dTs0) - dd_HTRn*(Ts_b-Ts0_b) 1137 !------------------------------------------------------------------------------------------- 1138 dRn_m_dRn0 = -(HTRn_b(j)+(sigx(j)-sigw(j))*dd_HTRn(j))*(dTs9(j)-dTs09(j)) &1139 - dd_HTRn(j)*(Ts_b9(j)-Ts0_b9(j))1140 print *,'wx_pbl_dts_check: j, dRn-dRn09, dRn_m_dRn0 ',j, dRn-dRn09(j), dRn_m_dRn01141 1142 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1143 ! Test phiT_b-phiT0_b = -sigx*sigw*dd_HTphiT*(dTs-dTs0) - HTphiT_b*(Ts_b-Ts0_b) 1144 !------------------------------------------------------------------------------- 1145 phiTb_m_phiT0b = -sigx(j)*sigw(j)*dd_HTphiT(j)*(dTs9(j)-dTs09(j)) - HTphiT_b(j)*(Ts_b9(j)-Ts0_b9(j))1146 print *,'wx_pbl_dts_check: j, phiT_b9-phiT0_b9, phiTb_m_phiT0b ',j ,phiT_b9(j)-phiT0_b9(j), phiTb_m_phiT0b1147 1148 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1149 ! Test phiT_w, phiT_x, dphiT from HTphiT 1150 !------------------------------------------ 1151 ! phiT_w = Kech_h_w C_p ( T1_w - Ts_w) phiT_x = Kech_h_x C_p ( T1_x - Ts_x) 1152 err_phiT_x(j) = Kech_h_x(j)*C_p(j)*(T1_x(j) - Ts_x(j)) - phiT_x9(j)1153 err_phiT_w(j) = Kech_h_w(j)*C_p(j)*(T1_w(j) - Ts_w(j)) - phiT_w9(j)1154 print *, 'wx_pbl_dts_check: j, phiT_w9, phiT_x9, err_phiT_w, err_phiT_x ', &1155 j, phiT_w9(j), phiT_x9(j), err_phiT_w(j), err_phiT_x(j)1156 dphiT = phiT_w9(j) - phiT_x9(j)1157 dphiT_H = dphiT09(j) - ( HTphiT_b(j) &1158 +(sigx(j)-sigw(j))*dd_HTphiT(j) &1159 -sigw(j)*sigx(j)*dd_HTphiT(j)*dd_HTphiT(j)/HTphiT_b(j) &1160 )*(dTs9(j)-dTs09(j)) &1161 + dd_HTphiT(j)/HTphiT_b(j)*(phiT_b9(j)-phiT0_b9(j))1162 print *,'wx_pbl_dts_check: j, dphiT, dphiT_H ', j, dphiT, dphiT_H1163 1164 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1165 ! Test phiq_w, phiq_x, dphiq from HTphiq 1166 !------------------------------------------ 1167 1168 ! phiq_w = beta Kech_q_w ( q1_w - qsat(Ts_w)) phiq_x = beta Kech_q_x ( q1_x - qsat(Ts_x)) 1169 err_phiq_x(j) = beta(j)*Kech_q_x(j)*( q1_x(j) - qsat_x(j)) - phiq_x9(j)1170 err_phiq_w(j) = beta(j)*Kech_q_w(j)*( q1_w(j) - qsat_w(j)) - phiq_w9(j)1171 dphiQ = phiQ_w9(j) - phiQ_x9(j)1172 dphiQ_H = dphiQ09(j) - ( HTphiQ_b(j) &1173 +(sigx(j)-sigw(j))*dd_HTphiQ(j) &1174 -sigw(j)*sigx(j)*dd_HTphiQ(j)*dd_HTphiT(j)/HTphiT_b(j) &1175 )*(dTs9(j)-dTs09(j)) &1176 + dd_HTphiQ(j)/HTphiT_b(j)*(phiT_b9(j)-phiT0_b9(j))1177 print *,'wx_pbl_dts_check: j, dphiQ, dphiQ_H ', j, dphiQ, dphiQ_H1178 1179 ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x 1180 err_phiT_b(j) = sigw(j)*phiT_w9(j) + sigx(j)*phiT_x9(j) - phiT_b9(j)1181 1182 ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x 1183 err_phiQ_b(j) = sigw(j)*phiQ_w9(j) + sigx(j)*phiQ_x9(j) - phiQ_b9(j)1184 1185 ! Ta = AcoefT + BcoefT phiT_b Delta t 1186 ! phiT_b = Kech_h C_p (Ta - Ts_b) 1187 Ta(j) = (AcoefT(j) + BcoefT(j)*phiT_b9(j)*dtime) / C_p(j)1188 err2_phiT_b(j) = Kech_h(j)*C_p(j)*(Ta(j) - Ts_b9(j)) - phiT_b9(j)1189 print *, 'wx_pbl_dts_check: j, Ta, phiT_b9, err2_phiT_b ', &1190 j, Ta(j), phiT_b9(j), err2_phiT_b(j)1191 1192 ENDDO ! j = 1, knon 1193 1194 ENDIF ! (prt_level >=10 ) 1195 !-------------------------------------------------------------------------------------------------- 1196 1197 1198 END SUBROUTINE wx_pbl_dts_check1199 1200 SUBROUTINE wx_evappot(knon, q1, Ts, evap_pot)959 USE lmdz_yomcst 960 961 IMPLICIT NONE 962 963 INTEGER, INTENT(IN) :: knon ! number of grid cells 964 REAL, INTENT(IN) :: dtime ! time step size (s) 965 REAL, DIMENSION(knon, klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) 966 REAL, DIMENSION(knon, klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) 967 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area 968 REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor 969 INTEGER, INTENT(IN) :: iflag_split 970 REAL, DIMENSION(knon), INTENT(IN) :: Ts0_b9, dTs09 971 REAL, DIMENSION(knon), INTENT(IN) :: qs_b9, Ts_b9 ! yqsurf, Tsurf_new 972 REAL, DIMENSION(knon), INTENT(IN) :: dTs9, dqsatsrf9 973 REAL, DIMENSION(knon), INTENT(IN) :: delta_qsurf9 974 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w 975 REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w 976 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0 977 978 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT, AcoefQ, BcoefT, BcoefQ 979 REAL, DIMENSION(knon), INTENT(IN) :: HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn 980 REAL, DIMENSION(knon), INTENT(IN) :: phiT0_b9, dphiT09, phiQ0_b9, dphiQ09, Rn0_b9, dRn09 981 REAL, DIMENSION(knon), INTENT(IN) :: g_T, g_Q 982 REAL, DIMENSION(knon), INTENT(IN) :: Gamma_phiT, Gamma_phiQ 983 REAL, DIMENSION(knon), INTENT(IN) :: dTs_ins, dqsatsrf_ins 984 REAL, DIMENSION(knon), INTENT(IN) :: phiT_b9, phiQ_b9 985 REAL, DIMENSION(knon), INTENT(IN) :: phiT_x9, phiT_w9 986 REAL, DIMENSION(knon), INTENT(IN) :: phiQ_x9, phiQ_w9 987 988 !! Local variables 989 INTEGER :: j 990 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region 991 REAL, DIMENSION(knon) :: AcoefT_b, AcoefQ_b ! mean values of AcoefT and AcoefQ 992 REAL :: zzt, zzq, zzqsat 993 REAL :: zdelta, zcvm5, zcor, qsat 994 REAL, DIMENSION(knon) :: qsat_w, qsat_x 995 REAL, DIMENSION(knon) :: Ts_x, Ts_w, qs_x, qs_w 996 REAL, DIMENSION(knon) :: T1_x, T1_w, q1_x, q1_w 997 REAL, DIMENSION(knon) :: Rn_x, Rn_w 998 REAL, DIMENSION(knon) :: Rn_b, dRn 999 REAL, DIMENSION(knon) :: phiQ0_x, phiQ0_w 1000 REAL, DIMENSION(knon) :: Ta, qa 1001 REAL, DIMENSION(knon) :: err_phiT_w, err_phiT_x 1002 REAL, DIMENSION(knon) :: err_phiq_w, err_phiq_x 1003 REAL, DIMENSION(knon) :: err_phiT_b 1004 REAL, DIMENSION(knon) :: err_phiQ_b 1005 REAL, DIMENSION(knon) :: err2_phiT_b 1006 REAL :: T1A_x, T1A_w, q1A_x, q1A_w 1007 REAL :: qsatsrf_w, qsatsrf_x, qsatsrfb, qsbA 1008 REAL :: dphiT, dphiQ 1009 REAL :: dphiT_H, dphiQ_H 1010 REAL :: phiQ_pot 1011 REAL :: phiQ_w_m_phiQ0_w 1012 REAL :: phiQ_x_m_phiQ0_x 1013 REAL :: dphiQ_m_dphiQ0 1014 REAL :: dphiT_m_dphiT0 1015 REAL :: dRN_m_dRn0 1016 REAL :: phiTb_m_phiT0b 1017 1018 !---------------------------------------------------------------------------- 1019 ! Equations to be checked: 1020 ! ----------------------- 1021 ! Input : Ts0_b, dTs0, Ts_b, dTs, qsatsrf_b, dqsatsrf, 1022 ! phiT_b, phiQ_b, phiT_w, phiT_x, phiQ_w, phiQ_x, 1023 1024 ! AcoefT, AcoefQ, AcoefT_w, AcoefQ_w, AcoefT_x, AcoefQ_x, 1025 ! BcoefT, BcoefQ, BcoefT_w, BcoefQ_w, BcoefT_x, BcoefQ_x 1026 1027 ! Ts_w = Ts_b + sigx dTs Ts_x = Ts_b - sigw dTs 1028 ! T1_w = AcoefT_w + BcoefT_w phiT_w Delta t T1_x = AcoefT_x + BcoefT_x phiT_x Delta t 1029 ! q1_w = AcoefQ_w + BcoefQ_w phiQ_w Delta t q1_x = AcoefQ_x + BcoefQ_x phiQ_x Delta t 1030 ! phiT_w = Kech_h_w ( T1_w - Ts_w) phiT_x = Kech_h_x ( T1_x - Ts_x) 1031 ! phiq_w = beta Kech_h_w ( q1_w - qsat(Ts_w)) phiq_x = beta Kech_h_x ( q1_x - qsat(Ts_x)) 1032 ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x 1033 ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x 1034 ! Ts_b = sigw Ts_w + sigx Ts_x dTs = Ts_w - Ts_x 1035 ! Ta = AcoefT + BcoefT phiT_b Delta t 1036 ! qa = AcoefQ + BcoefQ phiQ_b Delta t 1037 ! phiT_b = Kech_h (Ta - Ts_b) 1038 ! phiQ_b = beta Kech_h (qa - qsat(Ts_b)) 1039 ! dTs = sqrt(tau)/I (dphit + L_v dphiq + dR) 1040 1041 !---------------------------------------------------------------------------- 1042 1043 !! 1044 sigx(:) = 1. - sigw(:) 1045 AcoefT_b(1:knon) = AcoefT_x(1:knon) + sigw(1:knon) * dd_AT(1:knon) 1046 AcoefQ_b(1:knon) = AQ_x(1:knon) + sigw(1:knon) * dd_AQ(1:knon) 1047 1048 IF (prt_level >=10) THEN 1049 print *, '->wx_pbl_dts_check, HTphiT_b, HTphiQ_b, HTRn_b ', & 1050 HTphiT_b, HTphiQ_b, HTRn_b 1051 print *, '->wx_pbl_dts_check, dd_HTphiT, dd_HTphiQ, dd_HTRn ', & 1052 dd_HTphiT, dd_HTphiQ, dd_HTRn 1053 ENDIF ! (prt_level >=10 ) 1054 1055 ! Compute the three qsat and dqsatdTs 1056 ! --------------------------------------------- 1057 !! print *,' AAAA wx_pbl_dts_check, C_p(j), qsat0(j), Ts0(j) : ', & 1058 !! (C_p(j), qsat0(j), Ts0(j), j = 1,knon) 1059 1060 1061 !-------------------------------------------------------------------------------------------------- 1062 IF (prt_level >=10) THEN 1063 1064 DO j = 1, knon 1065 Ts_x(j) = Ts_b9(j) - sigw(j) * dTs9(j) 1066 Ts_w(j) = Ts_b9(j) + sigx(j) * dTs9(j) 1067 print *, 'wx_pbl_dts_check: j, Ts_b9, Ts_w, Ts_x ', j, Ts_b9(j), Ts_w(j), Ts_x(j) 1068 1069 qsat_x(j) = qsat0_x(j) + dqsatdT0_x(j) * (Ts_x(j) - Ts0_x(j)) 1070 qsat_w(j) = qsat0_w(j) + dqsatdT0_w(j) * (Ts_w(j) - Ts0_w(j)) 1071 1072 T1_x(j) = (AcoefT_x(j) + BcoefT_x(j) * phiT_x9(j) * dtime) / C_p(j) 1073 T1_w(j) = (AcoefT_w(j) + BcoefT_w(j) * phiT_w9(j) * dtime) / C_p(j) 1074 print *, 'wx_pbl_dts_check: j, T1_w, T1_x ', j, T1_w(j), T1_x(j) 1075 1076 q1_x(j) = AQ_x(j) + BQ_x(j) * phiQ_x9(j) * dtime 1077 q1_w(j) = AQ_w(j) + BQ_w(j) * phiQ_w9(j) * dtime 1078 print *, 'wx_pbl_dts_check: j, q1_w, q1_x ', j, q1_w(j), q1_x(j) 1079 1080 Rn_x(j) = eps_1 * Rsigma * T1_x(j)**4 - Rsigma * Ts_x(j)**4 1081 Rn_w(j) = eps_1 * Rsigma * T1_w(j)**4 - Rsigma * Ts_w(j)**4 1082 Rn_b(j) = sigw(j) * Rn_w(j) + sigx(j) * Rn_x(j) 1083 dRn(j) = dRn09(j) - (HTRn_b(j) & 1084 + (sigx(j) - sigw(j)) * dd_HTRn(j) & 1085 - sigw(j) * sigx(j) * dd_HTRn(j) * dd_HTphiT(j) / HTphiT_b(j) & 1086 ) * (dTs9(j) - dTs09(j)) & 1087 + dd_HTRn(j) / HTphiT_b(j) * (phiT_b9(j) - phiT0_b9(j)) 1088 1089 print *, 'wx_pbl_dts_check, dphiT, L_v*dphiQ, dRn, dTs ', & 1090 phiT_w9(j) - phiT_x9(j), L_v(j) * (phiQ_w9(j) - phiQ_x9(j)), dRn(j), dTs9(j) 1091 1092 phiQ0_x(j) = PhiQ0_b9(j) - sigw(j) * dphiQ09(j) 1093 phiQ0_w(j) = PhiQ0_b9(j) + sigx(j) * dphiQ09(j) 1094 1095 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1096 ! Test phiQ_w-phiQ0_w = -beta*Kech_Q_sw*dqsatdT_w*(Ts_w-Ts0_w) 1097 !-------------------------------------------------------------- 1098 print *, 'wx_pbl_dts_check: beta(j), Kech_Q_sw(j), dqsatdT0_w(j), Ts_w(j), Ts0_w(j) ', & 1099 beta(j), Kech_Q_sw(j), dqsatdT0_w(j), Ts_w(j), Ts0_w(j) 1100 phiQ_w_m_phiQ0_w = -beta(j) * Kech_Q_sw(j) * dqsatdT0_w(j) * (Ts_w(j) - Ts0_w(j)) 1101 print *, 'wx_pbl_dts_check: j, phiQ_w9-phiQ0_w, phiQ_w_m_phiQ0_w ', & 1102 j, phiQ_w9(j) - phiQ0_w(j), phiQ_w_m_phiQ0_w 1103 1104 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1105 ! Test phiQ_x-phiQ0_x = -beta*Kech_Q_sx*dqsatdT_x*(Ts_x-Ts0_x) 1106 !-------------------------------------------------------------- 1107 phiQ_x_m_phiQ0_x = -beta(j) * Kech_Q_sx(j) * dqsatdT0_x(j) * (Ts_x(j) - Ts0_x(j)) 1108 print *, 'wx_pbl_dts_check: j, phiQ_x9-phiQ0_x, phiQ_x_m_phiQ0_x ', & 1109 j, phiQ_x9(j) - phiQ0_x(j), phiQ_x_m_phiQ0_x 1110 1111 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1112 ! Test dphiT-dphiT0 = -(HTphiT_b+(sigx-sigw)*dd_HTphiT)*(dTs-dTs0) - dd_HTphiT*(Ts_b-Ts0_b) 1113 !------------------------------------------------------------------------------------------- 1114 dphiT = phiT_w9(j) - phiT_x9(j) 1115 dphiT_m_dphiT0 = -(HTphiT_b(j) + (sigx(j) - sigw(j)) * dd_HTphiT(j)) * (dTs9(j) - dTs09(j)) & 1116 - dd_HTphiT(j) * (Ts_b9(j) - Ts0_b9(j)) 1117 print *, 'wx_pbl_dts_check: j, dphiT-dphiT09, dphiT_m_dphiT0 ', j, dphiT - dphiT09(j), dphiT_m_dphiT0 1118 1119 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1120 ! Test dphiQ-dphiQ0 = -(HTphiQ_b+(sigx-sigw)*dd_HTphiQ)*(dTs-dTs0) - dd_HTphiQ*(Ts_b-Ts0_b) 1121 !------------------------------------------------------------------------------------------- 1122 dphiQ = phiQ_w9(j) - phiQ_x9(j) 1123 dphiQ_m_dphiQ0 = -(HTphiQ_b(j) + (sigx(j) - sigw(j)) * dd_HTphiQ(j)) * (dTs9(j) - dTs09(j)) & 1124 - dd_HTphiQ(j) * (Ts_b9(j) - Ts0_b9(j)) 1125 print *, 'wx_pbl_dts_check: j, dphiQ-dphiQ09, dphiQ_m_dphiQ0 ', j, dphiQ - dphiQ09(j), dphiQ_m_dphiQ0 1126 1127 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1128 ! Test dRn-dRn0 = -(HTRn_b+(sigx-sigw)*dd_HTRn)*(dTs-dTs0) - dd_HTRn*(Ts_b-Ts0_b) 1129 !------------------------------------------------------------------------------------------- 1130 dRn_m_dRn0 = -(HTRn_b(j) + (sigx(j) - sigw(j)) * dd_HTRn(j)) * (dTs9(j) - dTs09(j)) & 1131 - dd_HTRn(j) * (Ts_b9(j) - Ts0_b9(j)) 1132 print *, 'wx_pbl_dts_check: j, dRn-dRn09, dRn_m_dRn0 ', j, dRn - dRn09(j), dRn_m_dRn0 1133 1134 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1135 ! Test phiT_b-phiT0_b = -sigx*sigw*dd_HTphiT*(dTs-dTs0) - HTphiT_b*(Ts_b-Ts0_b) 1136 !------------------------------------------------------------------------------- 1137 phiTb_m_phiT0b = -sigx(j) * sigw(j) * dd_HTphiT(j) * (dTs9(j) - dTs09(j)) - HTphiT_b(j) * (Ts_b9(j) - Ts0_b9(j)) 1138 print *, 'wx_pbl_dts_check: j, phiT_b9-phiT0_b9, phiTb_m_phiT0b ', j, phiT_b9(j) - phiT0_b9(j), phiTb_m_phiT0b 1139 1140 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1141 ! Test phiT_w, phiT_x, dphiT from HTphiT 1142 !------------------------------------------ 1143 ! phiT_w = Kech_h_w C_p ( T1_w - Ts_w) phiT_x = Kech_h_x C_p ( T1_x - Ts_x) 1144 err_phiT_x(j) = Kech_h_x(j) * C_p(j) * (T1_x(j) - Ts_x(j)) - phiT_x9(j) 1145 err_phiT_w(j) = Kech_h_w(j) * C_p(j) * (T1_w(j) - Ts_w(j)) - phiT_w9(j) 1146 print *, 'wx_pbl_dts_check: j, phiT_w9, phiT_x9, err_phiT_w, err_phiT_x ', & 1147 j, phiT_w9(j), phiT_x9(j), err_phiT_w(j), err_phiT_x(j) 1148 dphiT = phiT_w9(j) - phiT_x9(j) 1149 dphiT_H = dphiT09(j) - (HTphiT_b(j) & 1150 + (sigx(j) - sigw(j)) * dd_HTphiT(j) & 1151 - sigw(j) * sigx(j) * dd_HTphiT(j) * dd_HTphiT(j) / HTphiT_b(j) & 1152 ) * (dTs9(j) - dTs09(j)) & 1153 + dd_HTphiT(j) / HTphiT_b(j) * (phiT_b9(j) - phiT0_b9(j)) 1154 print *, 'wx_pbl_dts_check: j, dphiT, dphiT_H ', j, dphiT, dphiT_H 1155 1156 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1157 ! Test phiq_w, phiq_x, dphiq from HTphiq 1158 !------------------------------------------ 1159 1160 ! phiq_w = beta Kech_q_w ( q1_w - qsat(Ts_w)) phiq_x = beta Kech_q_x ( q1_x - qsat(Ts_x)) 1161 err_phiq_x(j) = beta(j) * Kech_q_x(j) * (q1_x(j) - qsat_x(j)) - phiq_x9(j) 1162 err_phiq_w(j) = beta(j) * Kech_q_w(j) * (q1_w(j) - qsat_w(j)) - phiq_w9(j) 1163 dphiQ = phiQ_w9(j) - phiQ_x9(j) 1164 dphiQ_H = dphiQ09(j) - (HTphiQ_b(j) & 1165 + (sigx(j) - sigw(j)) * dd_HTphiQ(j) & 1166 - sigw(j) * sigx(j) * dd_HTphiQ(j) * dd_HTphiT(j) / HTphiT_b(j) & 1167 ) * (dTs9(j) - dTs09(j)) & 1168 + dd_HTphiQ(j) / HTphiT_b(j) * (phiT_b9(j) - phiT0_b9(j)) 1169 print *, 'wx_pbl_dts_check: j, dphiQ, dphiQ_H ', j, dphiQ, dphiQ_H 1170 1171 ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x 1172 err_phiT_b(j) = sigw(j) * phiT_w9(j) + sigx(j) * phiT_x9(j) - phiT_b9(j) 1173 1174 ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x 1175 err_phiQ_b(j) = sigw(j) * phiQ_w9(j) + sigx(j) * phiQ_x9(j) - phiQ_b9(j) 1176 1177 ! Ta = AcoefT + BcoefT phiT_b Delta t 1178 ! phiT_b = Kech_h C_p (Ta - Ts_b) 1179 Ta(j) = (AcoefT(j) + BcoefT(j) * phiT_b9(j) * dtime) / C_p(j) 1180 err2_phiT_b(j) = Kech_h(j) * C_p(j) * (Ta(j) - Ts_b9(j)) - phiT_b9(j) 1181 print *, 'wx_pbl_dts_check: j, Ta, phiT_b9, err2_phiT_b ', & 1182 j, Ta(j), phiT_b9(j), err2_phiT_b(j) 1183 1184 ENDDO ! j = 1, knon 1185 1186 ENDIF ! (prt_level >=10 ) 1187 !-------------------------------------------------------------------------------------------------- 1188 1189 END SUBROUTINE wx_pbl_dts_check 1190 1191 SUBROUTINE wx_evappot(knon, q1, Ts, evap_pot) 1201 1192 1202 1193 USE wx_pbl_var_mod 1203 1194 1204 INTEGER, INTENT(IN) :: knon ! number of grid cells 1205 REAL, DIMENSION(knon), INTENT(IN) :: q1 ! specific humidity in layer 1 1206 REAL, DIMENSION(knon), INTENT(IN) :: Ts ! surface temperature 1207 1208 REAL, DIMENSION(knon), INTENT(OUT) :: evap_pot ! potential evaporation 1209 1210 INTEGER :: j 1211 REAL :: qsat_bs 1212 1213 DO j = 1,knon 1214 evap_pot(j) = Kech_q(j)*(qsat0(j)+dqsatdT0(j)*(Ts(j)-Ts0(j))-q1(j)) 1215 1216 qsat_bs = qsat0(j)+dqsatdT0(j)*(Ts(j)-Ts0(j)) 1217 !! print *,'wx_evappot : Kech_q, qsat_bs, qa, evap_pot ', Kech_q(j), qsat_bs, q1(j), evap_pot(j) 1218 ENDDO 1219 1220 1221 END SUBROUTINE wx_evappot 1195 INTEGER, INTENT(IN) :: knon ! number of grid cells 1196 REAL, DIMENSION(knon), INTENT(IN) :: q1 ! specific humidity in layer 1 1197 REAL, DIMENSION(knon), INTENT(IN) :: Ts ! surface temperature 1198 1199 REAL, DIMENSION(knon), INTENT(OUT) :: evap_pot ! potential evaporation 1200 1201 INTEGER :: j 1202 REAL :: qsat_bs 1203 1204 DO j = 1, knon 1205 evap_pot(j) = Kech_q(j) * (qsat0(j) + dqsatdT0(j) * (Ts(j) - Ts0(j)) - q1(j)) 1206 1207 qsat_bs = qsat0(j) + dqsatdT0(j) * (Ts(j) - Ts0(j)) 1208 !! print *,'wx_evappot : Kech_q, qsat_bs, qa, evap_pot ', Kech_q(j), qsat_bs, q1(j), evap_pot(j) 1209 ENDDO 1210 1211 END SUBROUTINE wx_evappot 1222 1212 1223 1213 END MODULE wx_pbl_mod
Note: See TracChangeset
for help on using the changeset viewer.
