MODULE wx_pbl_mod ! Split Planetary Boundary Layer ! This module manages the splitting of the boundary layer between two regions; the (w) ! region (inside cold pools) and the (x) region (outside cold pools) USE dimphy IMPLICIT NONE CONTAINS !**************************************************************************************** SUBROUTINE wx_pbl0_merge(knon, ypplay, ypaprs, & sigw, dTs_forcing, dqs_forcing, & yt_x, yt_w, yq_x, yq_w, & yu_x, yu_w, yv_x, yv_w, & ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w, & ycdragm_x, ycdragm_w, & AcoefT_x, AcoefT_w, AcoefQ_x, AcoefQ_w, & AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, & BcoefT_x, BcoefT_w, BcoefQ_x, BcoefQ_w, & BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, & AcoefT, AcoefQ, AcoefU, AcoefV, & BcoefT, BcoefQ, BcoefU, BcoefV, & ycdragh, ycdragq, ycdragm, & yt1, yq1, yu1, yv1 & ) USE wx_pbl_var_mod USE lmdz_print_control, ONLY: prt_level,lunout USE indice_sol_mod, ONLY: is_oce INCLUDE "YOMCST.h" INCLUDE "FCTTRE.h" INCLUDE "YOETHF.h" INCLUDE "clesphys.h" INTEGER, INTENT(IN) :: knon ! number of grid cells REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area REAL, DIMENSION(knon), INTENT(IN) :: dTs_forcing ! forced temperature difference (w)-(x) REAL, DIMENSION(knon), INTENT(IN) :: dqs_forcing ! forced humidity difference (w)-(x) REAL, DIMENSION(knon,klev), INTENT(IN) :: yt_x, yt_w, yq_x, yq_w REAL, DIMENSION(knon,klev), INTENT(IN) :: yu_x, yu_w, yv_x, yv_w REAL, DIMENSION(knon), INTENT(IN) :: ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w REAL, DIMENSION(knon), INTENT(IN) :: ycdragm_x, ycdragm_w REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w, AcoefQ_x, AcoefQ_w REAL, DIMENSION(knon), INTENT(IN) :: AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w, BcoefQ_x, BcoefQ_w REAL, DIMENSION(knon), INTENT(IN) :: BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w REAL, DIMENSION(knon), INTENT(OUT) :: AcoefT, AcoefQ, AcoefU, AcoefV REAL, DIMENSION(knon), INTENT(OUT) :: BcoefT, BcoefQ, BcoefU, BcoefV REAL, DIMENSION(knon), INTENT(OUT) :: ycdragh, ycdragq, ycdragm REAL, DIMENSION(knon), INTENT(OUT) :: yt1, yq1, yu1, yv1 ! Apparent T, q, u, v at first level, as !seen by surface modules ! Local variables INTEGER :: j REAL :: dd_Kh REAL :: dd_Kq REAL :: dd_Km REAL :: dd_u REAL :: dd_v REAL :: dd_t REAL :: dd_q REAL :: LambdaTs, LambdaQs, LambdaUs, LambdaVs REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region sigx(1:knon) = 1.-sigw(1:knon) DO j=1,knon ! Compute w-x differences dd_t = yt_w(j,1) - yt_x(j,1) dd_q = yq_w(j,1) - yq_x(j,1) dd_u = yu_w(j,1) - yu_x(j,1) dd_v = yv_w(j,1) - yv_x(j,1) ! Merged exchange coefficients dd_Kh = Kech_h_w(j) - Kech_h_x(j) dd_Kq = Kech_q_w(j) - Kech_q_x(j) dd_Km = Kech_m_w(j) - Kech_m_x(j) LambdaTs = dd_KTp(j)/Kech_Tp(j) LambdaQs = dd_KQs(j)/Kech_Qs(j) LambdaUs = dd_KUp(j)/Kech_Up(j) LambdaVs = dd_KVp(j)/Kech_Vp(j) ! Calcul des coef A, B \'equivalents dans la couche 1 ! The dTs_forcing and dqs_forcing terms are added for diagnostic purpose ; they should be zero in normal operation. AcoefT(j) = AcoefT_x(j) + sigw(j)*(1.+sigx(j)*LambdaTs)*(dd_AT(j) - C_p(j)*dTs_forcing(j)) AcoefQ(j) = AcoefQ_x(j) + sigw(j)*(1.+sigx(j)*LambdaQs)*(dd_AQ(j) - dqs_forcing(j)) AcoefU(j) = AcoefU_x(j) + sigw(j)*(1.+sigx(j)*LambdaUs)*dd_AU(j) AcoefV(j) = AcoefV_x(j) + sigw(j)*(1.+sigx(j)*LambdaVs)*dd_AV(j) !! BcoefT(j) = (sigw(j)*Kech_h_w(j)*Kech_T_pw(j)*BcoefT_w(j) + & !! sigx(j)*Kech_h_x(j)*Kech_T_px(j)*BcoefT_x(j) )/(Kech_h(j)*Kech_Tp(j)) !! BcoefQ(j) = (sigw(j)*Kech_q_w(j)*Kech_Q_pw(j)*BcoefQ_w(j) + & !! sigx(j)*Kech_q_x(j)*Kech_Q_px(j)*BcoefQ_x(j) )/(Kech_q(j)*Kech_Qp(j)) !! BcoefU(j) = (sigw(j)*Kech_m_w(j)*Kech_U_pw(j)*BcoefU_w(j) + & !! sigx(j)*Kech_m_x(j)*Kech_U_px(j)*BcoefU_x(j) )/(Kech_m(j)*Kech_Up(j)) !! BcoefV(j) = (sigw(j)*Kech_m_w(j)*Kech_V_pw(j)*BcoefV_w(j) + & !! sigx(j)*Kech_m_x(j)*Kech_V_px(j)*BcoefV_x(j) )/(Kech_m(j)*Kech_Vp(j)) !! Print *,'YYYYpbl0: BcoefT_x, sigw, sigx, dd_Kh, dd_KTp, Kech_h_w ', & !! BcoefT_x, sigw, sigx, dd_Kh, dd_KTp, Kech_h_w !! Print *,'YYYYpbl0: Kech_T_pw, dd_BT, Kech_h, Kech_Tp ', & !! Kech_T_pw, dd_BT, Kech_h, Kech_Tp BcoefT(j) = BcoefT_x(j) + sigw(j)*(sigx(j)*dd_Kh*dd_KTp(j)*BcoefT_x(j) + & Kech_h_w(j)*Kech_T_pw(j)*dd_BT(j))/(Kech_h(j)*Kech_Tp(j)) BcoefQ(j) = BcoefQ_x(j) + sigw(j)*(sigx(j)*dd_Kq*dd_KQs(j)*BcoefQ_x(j) + & Kech_q_w(j)*Kech_Q_sw(j)*dd_BQ(j))/(Kech_q(j)*Kech_Qs(j)) BcoefU(j) = BcoefU_x(j) + sigw(j)*(sigx(j)*dd_Km*dd_KUp(j)*BcoefU_x(j) + & Kech_m_w(j)*Kech_U_pw(j)*dd_BU(j))/(Kech_m(j)*Kech_Up(j)) BcoefV(j) = BcoefV_x(j) + sigw(j)*(sigx(j)*dd_Km*dd_KVp(j)*BcoefV_x(j) + & Kech_m_w(j)*Kech_V_pw(j)*dd_BV(j))/(Kech_m(j)*Kech_Vp(j)) !>jyg ! Calcul des cdrag \'equivalents dans la couche ycdragm(j) = ycdragm_x(j) + sigw(j)*dd_Cdragm(j) ycdragh(j) = ycdragh_x(j) + sigw(j)*dd_Cdragh(j) ycdragq(j) = ycdragq_x(j) + sigw(j)*dd_Cdragq(j) ! Calcul de T, q, u et v \'equivalents dans la couche 1 !! yt1(j) = yt_x(j,1) + sigw(j)*dd_t*(1.+sigx(j)*dd_Kh/KCT) !! yq1(j) = yq_x(j,1) + sigw(j)*dd_q*(1.+sigx(j)*dd_Kh/KCQ) !! yu1(j) = yu_x(j,1) + sigw(j)*dd_u*(1.+sigx(j)*dd_Km/KCU) !! yv1(j) = yv_x(j,1) + sigw(j)*dd_v*(1.+sigx(j)*dd_Km/KCV) yt1(j) = yt_x(j,1) + sigw(j)*dd_t yq1(j) = yq_x(j,1) + sigw(j)*dd_q yu1(j) = yu_x(j,1) + sigw(j)*dd_u yv1(j) = yv_x(j,1) + sigw(j)*dd_v ENDDO END SUBROUTINE wx_pbl0_merge SUBROUTINE wx_pbl_dts_merge(knon, dtime, ypplay, ypaprs, & sigw, beta, wcstar, wdens, & AT_x, AT_w, & BT_x, BT_w, & AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, & AcoefT, AcoefQ, BcoefT, BcoefQ, & HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, & phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0, & g_T, g_Q, & Gamma_phiT, Gamma_phiQ, & dTs_ins, dqsatsrf_ins & ) USE wx_pbl_var_mod USE lmdz_print_control, ONLY: prt_level,lunout INCLUDE "YOMCST.h" INCLUDE "FCTTRE.h" INCLUDE "YOETHF.h" INTEGER, INTENT(IN) :: knon ! number of grid cells REAL, INTENT(IN) :: dtime ! time step size (s) REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pool fractional area REAL, DIMENSION(knon), INTENT(IN) :: beta ! evaporation by potential evaporation REAL, DIMENSION(knon), INTENT(IN) :: wcstar ! cold pool gust front speed REAL, DIMENSION(knon), INTENT(IN) :: wdens ! cold pool number density REAL, DIMENSION(knon), INTENT(IN) :: AT_x, AT_w REAL, DIMENSION(knon), INTENT(IN) :: BT_x, BT_w REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0 REAL, DIMENSION(knon), INTENT(OUT) :: AcoefT, AcoefQ, BcoefT, BcoefQ REAL, DIMENSION(knon), INTENT(OUT) :: HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn REAL, DIMENSION(knon), INTENT(OUT) :: phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0 REAL, DIMENSION(knon), INTENT(OUT) :: g_T, g_Q REAL, DIMENSION(knon), INTENT(OUT) :: Gamma_phiT, Gamma_phiQ REAL, DIMENSION(knon), INTENT(OUT) :: dTs_ins, dqsatsrf_ins ! Local variables REAL, DIMENSION(knon) :: qsat_x REAL, DIMENSION(knon) :: qsat_w REAL, DIMENSION(knon) :: dqsatdT_x REAL, DIMENSION(knon) :: dqsatdT_w REAL, DIMENSION(knon) :: T10_x REAL, DIMENSION(knon) :: T10_w REAL, DIMENSION(knon) :: phiT0_x REAL, DIMENSION(knon) :: phiT0_w REAL, DIMENSION(knon) :: phiQ0_x REAL, DIMENSION(knon) :: phiQ0_w REAL, DIMENSION(knon) :: Rn0_x REAL, DIMENSION(knon) :: Rn0_w REAL, DIMENSION(knon) :: Rp1_x REAL, DIMENSION(knon) :: Rp1_w REAL, DIMENSION(knon) :: Rps_x REAL, DIMENSION(knon) :: Rps_w REAL, DIMENSION(knon) :: HTphiT_x REAL, DIMENSION(knon) :: HTphiT_w REAL, DIMENSION(knon) :: HTphiQ_x REAL, DIMENSION(knon) :: HTphiQ_w REAL, DIMENSION(knon) :: HTRn_x REAL, DIMENSION(knon) :: HTRn_w REAL, DIMENSION(knon) :: HQphiT_x REAL, DIMENSION(knon) :: HQphiT_w REAL, DIMENSION(knon) :: HQphiQ_x REAL, DIMENSION(knon) :: HQphiQ_w REAL, DIMENSION(knon) :: HQRn_x REAL, DIMENSION(knon) :: HQRn_w REAL, DIMENSION(knon) :: HQphiT_b REAL, DIMENSION(knon) :: dd_HQphiT REAL, DIMENSION(knon) :: HQphiQ_b REAL, DIMENSION(knon) :: dd_HQphiQ REAL, DIMENSION(knon) :: HQRn_b REAL, DIMENSION(knon) :: dd_HQRn REAL, DIMENSION(knon) :: sigx REAL, DIMENSION(knon) :: Ts, T1 !!! REAL, DIMENSION(knon) :: qsat, dqsat_dT !!! REAL, DIMENSION(knon) :: phiT0 !!! REAL, DIMENSION(knon) :: Cp, Lv REAL, DIMENSION(knon) :: tau, Inert REAL :: dd_Kh REAL :: zdelta, zcvm5, zcor REAL :: qsat INTEGER :: j !---------------------------------------------------------------------------- ! Reference state ! --------------- ! dqsat_dT_w = dqsat_dT(Ts0_w) dqsat_dT_x = dqsat_dT(Ts0_x) ! T10_w = (AT_w/Cp - Kech_T_w BT_w dtime Ts0_w)/(1 - Kech_T_w BT_w dtime) ! T10_x = (AT_x/Cp - Kech_T_x BT_x dtime Ts0_x)/(1 - Kech_T_x BT_x dtime) ! phiT0_w = Kech_T_pw (AT_w - Cp Ts0_w) phiT0_x = Kech_T_px (AT_x - Cp Ts0_x) ! phiQ0_w = Kech_Q_sw (beta AQ_w - qsatsrf0_w) phiQ0_x = Kech_Q_sx (beta AQ_x - qsatsrf0_x) ! Rn0_w = eps_1 Rsigma T10_w^4 - Rsigma Ts0_w^4 Rn0_x = eps_1 Rsigma T10_x^4 - Rsigma Ts0_x^4 ! Rp1_w = 4 eps_1 Rsigma T10_w^3 Rp1_x = 4 eps_1 Rsigma T10_x^3 ! Rps_w = 4 Rsigma Ts0_w^3 Rps_x = 4 Rsigma Ts0_x^3 ! phiT0_b = sigw phiT0_w + sigx phiT0_x ! dphiT0 = phiT0_w - phiT0_x ! phiQ0_b = sigw phiQ0_w + sigx phiQ0_x ! dphiQ0 = phiQ0_w - phiQ0_x ! Rn0_b = sigw Rn0_w + sigx Rn0_x dRn0 = Rn0_w - Rn0_x !---------------------------------------------------------------------------- ! Elementary enthalpy equations ! ----------------------------- ! phiT_w = phiT0_w - HTphiT_w (Ts_w-Ts0_w) phiT_x = phiT0_x - HTphiT_x (Ts_x-Ts0_x) ! phiQ_w = phiQ0_w - HTphiQ_w (Ts_w-Ts0_w) phiQ_x = phiQ0_x - HTphiQ_x (Ts_x-Ts0_x) ! Rn_w = Rn0_w - HTRn_w (Ts_w-Ts0_w) Rn_x = Rn0_x - HTRn_x (Ts_x-Ts0_x) ! DFlux_DT coefficients ! --------------------- ! Heat flux equation ! HTphiT_w = Cp Kech_T_pw HTphiT_x = Cp Kech_T_px ! Moisture flux equation ! HTphiQ_w = beta Kech_Q_sw dqsat_dT_w HTphiQ_x = beta Kech_Q_sx dqsat_dT_x ! Radiation equation ! HTRn_w = Rp1_w Kech_T_pw BcoefT_w dtime + Rps_w HTRn_x = Rp1_x Kech_T_px BcoefT_x dtime + Rps_x !---------------------------------------------------------------------------- ! Elementary moisture equations ! ----------------------------- ! beta Ts_w = beta Ts0_w + QQ_w (qsatsrf_w-qsatsrf0_w) beta Ts_x = beta Ts0_x + QQ_x (qsatsrf_x-qsatsrf0_x) ! beta phiT_w = beta phiT0_w - HQphiT_w (qsatsrf_w-qsatsrf0_w) beta phiQ_x = beta phiQ0_x - HTphiQ_x (qsatsrf_x-qsatsrf0_x) ! beta phiQ_w = beta phiQ0_w - HQphiQ_w (qsatsrf_w-qsatsrf0_w) beta phiQ_x = beta phiQ0_x - HTphiQ_x (qsatsrf_x-qsatsrf0_x) ! beta Rn_w = beta Rn0_w - HQRn_w (qsatsrf_w-qsatsrf0_w) beta Rn_x = beta Rn0_x - HTRn_x (qsatsrf_x-qsatsrf0_x) ! DFluxDQ coefficients ! --------------------- ! dqsat_dT equation ! QQ_w = 1. / dqsat_dT_w QQ_x = 1. / dqsat_dT_x ! Heat flux equation ! HQphiT_w = Cp Kech_T_pw QQ_w HQphiT_x = Cp Kech_T_px QQ_x ! Moisture flux equation ! HQphiQ_w = beta Kech_Q_sw HQphiQ_x = beta Kech_Q_sx ! Radiation equation ! HQRn_w = (Rp1_w Kech_T_pw BcoefT_w dtime + Rps_w) QQ_w ! HQRn_x = (Rp1_x Kech_T_px BcoefT_x dtime + Rps_x) QQ_x !---------------------------------------------------------------------------- ! Mean values and w-x differences ! ------------------------------- ! HTphiT_b = sigw HTphiT_w + sigx HTphiT_x dd_HTphiT = HTphiT_w - HTphiT_x ! HTphiQ_b = sigw HTphiQ_w + sigx HTphiQ_x dd_HTphiQ = HTphiQ_w - HTphiQ_x ! HTRn_b = sigw HTRn_w + sigx HTRn_x dd_HTRn = HTRn_w - HTRn_x ! QQ_b = sigw QQ_w + sigx QQ_x dd_QQ = QQ_w - QQ_x ! HQphiT_b = sigw HQphiT_w + sigx HQphiT_x dd_HQphiT = HQphiT_w - HQphiT_x ! HQphiQ_b = sigw HQphiQ_w + sigx HQphiQ_x dd_HQphiQ = HQphiQ_w - HQphiQ_x ! HQRn_b = sigw HQRn_w + sigx HQRn_x dd_HQRn = HQRn_w - HQRn_x !---------------------------------------------------------------------------- ! Equations ! --------- ! (1 - g_T) dTs = dTs_ins + Gamma_phiT phiT ! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ ! Feedback Gains ! -------------- ! g_T = - (sqrt(tau)/I) [ HTphiT_b + Lv HTphiQ_b + HTRn_b + & ! (dd_HTphiT + Lv dd_HTphiQ + dd_HTRn) (sigx - sigw - sigw sigx dd_HTphiT/HTphiT_b) ] ! g_Q = - (sqrt(tau)/(I QQ_b)) ( HQphiT_b + Lv HQphiQ_b + HQRn_b ) - & ! (sigx - sigw - sigw sigx dd_HQphiQ/HQphiQ_b) & ! [ dd_QQ/QQ_b + (sqrt(tau)/(I QQ_b))(dd_HQphiT + Lv dd_HQphiQ + dd_HQRn) ] ! Ts, qs Coupling coefficients / ! ---------------------------- ! Gamma_phiT = (sqrt(tau)/(I HTphiT_b)) (dd_HTphiT + Lv dd_HTphiQ + dd_HTRn) ! Gamma_phiQ = (1/(HQphiQ_b QQ_b)) [ dd_QQ + (sqrt(tau)/(I )) (dd_HQphiT + Lv dd_HQphiQ + dd_HQRn) ] ! Insensitive changes ! ------------------- ! dTs_ins = (1 - g_T) dTs0 - Gamma_phiT phiT0_b ! dqsatsrf_ins = (1 - g_Q) dqsatsrf0 - Gamma_phiQ phiQ0_b !---------------------------------------------------------------------------- ! Effective coefficients Acoef and Bcoef ! -------------------------------------- ! Equations ! --------- ! Cp Ta = AcoefT + BcoefT phiT dtime ! qa = AcoefQ + BcoefQ phiQ dtime ! Coefficients ! ------------ ! AcoefT = AcoefT0 - sigw sigx (dd_KTp/Kech_Tp) Cp dTs_ins/(1 - g_T) ! BcoefT = BcoefT0 - sigw sigx (dd_KTp/Kech_Tp) Cp Gamma_phiT/(1 - g_T)/dtime ! AcoefQ = AcoefQ0 - sigw sigx (dd_KQp/Kech_Qp) dqs_ins/(1 - g_Q) ! BcoefQ = BcoefQ0 - sigw sigx (dd_KQp/Kech_Qp) Gamma_phiq/(1 - g_Q)/dtime !============================================================================== ! Parameters ! ---------- Inert(1:knon) = 2000. tau(1:knon) = sqrt(sigw(1:knon)/max(rpi*wdens(1:knon)*wcstar(1:knon)**2 , & sigw(1:knon)*1.e-12,smallestreal)) sigx(1:knon) = 1.-sigw(1:knon) !! Compute Cp, Lv, qsat, dqsat_dT. ! C_p(1:knon) = RCpd ! L_v(1:knon) = RLvtt ! print *,' AAAA wx_pbl_dTs, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:) 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))/ & (1 - Kech_h_x(1:knon)*BT_x(1:knon)*dtime) 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))/ & (1 - Kech_h_w(1:knon)*BT_w(1:knon)*dtime) phiT0_x(1:knon) = Kech_T_px(1:knon)*(AT_x(1:knon) - C_p(1:knon)*Ts0_x(1:knon)) phiT0_w(1:knon) = Kech_T_pw(1:knon)*(AT_w(1:knon) - C_p(1:knon)*Ts0_w(1:knon)) phiQ0_x(1:knon) = Kech_Q_sx(1:knon)*(beta(1:knon)*AQ_x(1:knon) - qsatsrf0_x(1:knon)) phiQ0_w(1:knon) = Kech_Q_sw(1:knon)*(beta(1:knon)*AQ_w(1:knon) - qsatsrf0_w(1:knon)) Rn0_x(1:knon) = eps_1*Rsigma*T10_x(1:knon)**4 - Rsigma*Ts0_x(1:knon)**4 Rn0_w(1:knon) = eps_1*Rsigma*T10_w(1:knon)**4 - Rsigma*Ts0_w(1:knon)**4 Rp1_x(1:knon) = 4*eps_1*Rsigma*T10_x(1:knon)**3 Rp1_w(1:knon) = 4*eps_1*Rsigma*T10_w(1:knon)**3 Rps_x(1:knon) = 4*Rsigma*Ts0_x(1:knon)**3 Rps_w(1:knon) = 4*Rsigma*Ts0_w(1:knon)**3 ! DFlux_DT coefficients ! --------------------- ! Heat flux equation HTphiT_x(1:knon) = C_p(1:knon)*Kech_T_px(1:knon) HTphiT_w(1:knon) = C_p(1:knon)*Kech_T_pw(1:knon) ! Moisture flux equation HTphiQ_x(1:knon) = beta(1:knon)*Kech_Q_sx(1:knon)*dqsatdT0_x(1:knon) HTphiQ_w(1:knon) = beta(1:knon)*Kech_Q_sw(1:knon)*dqsatdT0_w(1:knon) ! Radiation equation HTRn_x(1:knon) = Rp1_x(1:knon)*Kech_T_px(1:knon)*BT_x(1:knon)*dtime + Rps_x(1:knon) HTRn_w(1:knon) = Rp1_w(1:knon)*Kech_T_pw(1:knon)*BT_w(1:knon)*dtime + Rps_w(1:knon) ! DFluxDQ coefficients ! --------------------- ! Heat flux equation HQphiT_x(1:knon) = C_p(1:knon)*Kech_T_px(1:knon)*QQ_x(1:knon) HQphiT_w(1:knon) = C_p(1:knon)*Kech_T_pw(1:knon)*QQ_w(1:knon) ! Moisture flux equation HQphiQ_x(1:knon) = beta(1:knon)*Kech_Q_sx(1:knon) HQphiQ_w(1:knon) = beta(1:knon)*Kech_Q_sw(1:knon) ! Radiation equation 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) 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) ! Mean values and w-x differences ! ------------------------------- phiT0_b(1:knon) = sigw(1:knon)*phiT0_w(1:knon) + sigx(1:knon)*phiT0_x(1:knon) phiQ0_b(1:knon) = sigw(1:knon)*phiQ0_w(1:knon) + sigx(1:knon)*phiQ0_x(1:knon) Rn0_b(1:knon) = sigw(1:knon)*Rn0_w(1:knon) + sigx(1:knon)*Rn0_x(1:knon) dphiT0(1:knon) = phiT0_w(1:knon) - phiT0_x(1:knon) dphiQ0(1:knon) = phiQ0_w(1:knon) - phiQ0_x(1:knon) dRn0(1:knon) = Rn0_w(1:knon) - Rn0_x(1:knon) HTphiT_b(1:knon) = sigw(1:knon)*HTphiT_w(1:knon) + sigx(1:knon)*HTphiT_x(1:knon) dd_HTphiT(1:knon) = HTphiT_w(1:knon) - HTphiT_x(1:knon) HTphiQ_b(1:knon) = sigw(1:knon)*HTphiQ_w(1:knon) + sigx(1:knon)*HTphiQ_x(1:knon) dd_HTphiQ(1:knon) = HTphiQ_w(1:knon) - HTphiQ_x(1:knon) HTRn_b(1:knon) = sigw(1:knon)*HTRn_w(1:knon) + sigx(1:knon)*HTRn_x(1:knon) dd_HTRn(1:knon) = HTRn_w(1:knon) - HTRn_x(1:knon) HQphiT_b(1:knon) = sigw(1:knon)*HQphiT_w(1:knon) + sigx(1:knon)*HQphiT_x(1:knon) dd_HQphiT(1:knon) = HQphiT_w(1:knon) - HQphiT_x(1:knon) HQphiQ_b(1:knon) = sigw(1:knon)*HQphiQ_w(1:knon) + sigx(1:knon)*HQphiQ_x(1:knon) dd_HQphiQ(1:knon) = HQphiQ_w - HQphiQ_x(1:knon) HQRn_b(1:knon) = sigw(1:knon)*HQRn_w(1:knon) + sigx(1:knon)*HQRn_x(1:knon) dd_HQRn(1:knon) = HQRn_w(1:knon) - HQRn_x(1:knon) ! Feedback Gains ! -------------- g_T(1:knon) = - (sqrt(tau(1:knon))/Inert(1:knon)) & * (HTphiT_b(1:knon) + L_v(1:knon)*HTphiQ_b(1:knon) + HTRn_b(1:knon) & + (dd_HTphiT(1:knon) + L_v(1:knon)*dd_HTphiQ(1:knon) + dd_HTRn(1:knon)) & * (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HTphiT(1:knon)/HTphiT_b(1:knon)) ) !!!! DO j = 1,knon !!!! IF (mod(j,20) .EQ.0) THEN !!!! print *, ' j dd_QQ QQ_b dd_HQphiQ dd_HQphiT dd_HQRn HQphiQ_b HQphiT_b HQRn_b ' !!!! ENDIF !!!! 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) !!!! 1789 FORMAT( I4, 10(1X,E10.2)) !!!! ENDDO g_Q(1:knon) = - (dd_QQ(1:knon)/QQ_b(1:knon)) * & (sigx(1:knon)-sigw(1:knon)-sigw(1:knon)*sigx(1:knon)*dd_KQs(1:knon)/Kech_Qs(1:knon)) & - sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon)) * & ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) + & (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_KQs(1:knon)/Kech_Qs(1:knon)) * & (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) !! g_Q(1:knon) = - (dd_QQ(1:knon)/QQ_b(1:knon)) * & !! (sigx(1:knon)-sigw(1:knon)-sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) & !! - sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon)) * & !! ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) + & !! (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) * & !! (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) !! g_Q(1:knon) = - (sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon))) * & !! ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) ) & !! - (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) * & !! ( dd_QQ(1:knon)/QQ_b(1:knon) & !! + (sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon))) & !! * (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) ! Ts, qs Coupling coefficients / ! ---------------------------- Gamma_phiT(1:knon) = (sqrt(tau(1:knon))/(Inert(1:knon)*HTphiT_b(1:knon))) & * (dd_HTphiT(1:knon) + L_v(1:knon)*dd_HTphiQ(1:knon) + dd_HTRn(1:knon)) Gamma_phiQ(1:knon) = (1./(Kech_Qs(1:knon)*QQ_b(1:knon))) * & ( dd_QQ(1:knon) & + (sqrt(tau(1:knon))/(Inert(1:knon))) * & (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) !! Gamma_phiQ(1:knon) = (beta(1:knon)/(HQphiQ_b(1:knon)*QQ_b(1:knon))) * & !! ( dd_QQ(1:knon) & !! + (sqrt(tau(1:knon))/(Inert(1:knon))) * & !! (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) !! Gamma_phiQ(1:knon) = (1/(HQphiQ_b(1:knon)*QQ_b(1:knon))) & !! * ( dd_QQ(1:knon) & !! + (sqrt(tau(1:knon))/(Inert(1:knon))) & !! * (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) ) ! Insensitive changes ! ------------------- dTs_ins(1:knon) = (sqrt(tau(1:knon))/Inert(1:knon))* & (dphiT0(1:knon) + L_v(1:knon)*dphiQ0(1:knon) + dRn0(1:knon)) dqsatsrf_ins(1:knon) = (beta(1:knon)/QQ_b(1:knon))*dTs_ins(1:knon) IF (prt_level >= 10) THEN print *,'wx_pbl_merge, tau ', tau print *,'wx_pbl_merge, AcoefT0 ', AcoefT0 print *,'wx_pbl_merge, AcoefQ0 ', AcoefQ0 print *,'wx_pbl_merge, BcoefT0 ', BcoefT0 print *,'wx_pbl_merge, BcoefQ0 ', BcoefQ0 print *,'wx_pbl_merge, qsat0_w, qsat0_x ', (qsat0_w(j), qsat0_x(j),j=1,knon) print *,'wx_pbl_merge, dqsatdT0_w, dqsatdT0_x ', (dqsatdT0_w(j), dqsatdT0_x(j),j=1,knon) ENDIF !---------------------------------------------------------------------------- !------------------------------------------------------------------------------ ! Effective coefficients Acoef and Bcoef ! -------------------------------------- DO j = 1,knon AcoefT(j) = AcoefT0(j) - sigw(j)*sigx(j)*(dd_KTp(j)/Kech_Tp(j))*C_p(j)* & (dTs0(j) + (dTs_ins(j)-dTs0(j)-Gamma_phiT(j)*phiT0_b(j))/(1. - g_T(j))) 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 AcoefQ(j) = AcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))* & (dqsatsrf0(j) + (dqsatsrf_ins(j)-(beta(j)/QQ_b(j))*dTs0(j)-Gamma_phiQ(j)*phiQ0_b(j))/(1 - g_Q(j)))/ & max(beta(j),1.e-4) BcoefQ(j) = BcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))*Gamma_phiQ(j)/(1 - g_Q(j))/ & (max(beta(j),1.e-4)*dtime) !! AcoefQ(j) = AcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))* & !! (dqsatsrf0(j) + (dqsatsrf_ins(j)-(beta(j)/QQ_b(j))*dTs0(j)-Gamma_phiQ(j)*phiQ0_b(j))/(1 - g_Q(j)))/ & !! beta(j) !! BcoefQ(j) = BcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))*Gamma_phiQ(j)/(1 - g_Q(j))/(beta(j)*dtime) ENDDO ! j = 1,knon IF (prt_level >= 10) THEN print *,'wx_pbl_dts AAAA BcoefQ, BcoefQ0, sigw ', & BcoefQ, BcoefQ0, sigw print *,'wx_pbl_dts_merge, dTs_ins ', dTs_ins print *,'wx_pbl_dts_merge, dqs_ins ', dqsatsrf_ins ENDIF END SUBROUTINE wx_pbl_dts_merge SUBROUTINE wx_pbl_split(knon, nsrf, dtime, sigw, beta, iflag_split, & g_T, g_Q, & Gamma_phiT, Gamma_phiQ, & dTs_ins, dqsatsrf_ins, & phiT, phiQ, phiU, phiV, & !!!! HTRn_b, dd_HTRn, HTphiT_b, dd_HTphiT, & phiQ0_b, phiT0_b, & phiT_x, phiT_w, & phiQ_x, phiQ_w, & phiU_x, phiU_w, & phiV_x, phiV_w, & philat_x, philat_w, & !!!! Rn_b, dRn, & dqsatsrf, & dTs, delta_qsurf & ) USE wx_pbl_var_mod USE lmdz_print_control, ONLY: prt_level,lunout USE indice_sol_mod, ONLY: is_oce INCLUDE "YOMCST.h" INTEGER, INTENT(IN) :: knon ! number of grid cells INTEGER, INTENT(IN) :: nsrf ! surface type REAL, INTENT(IN) :: dtime ! time step size (s) REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor INTEGER, INTENT(IN) :: iflag_split REAL, DIMENSION(knon), INTENT(IN) :: g_T, g_Q REAL, DIMENSION(knon), INTENT(IN) :: Gamma_phiT, Gamma_phiQ REAL, DIMENSION(knon), INTENT(IN) :: dTs_ins, dqsatsrf_ins REAL, DIMENSION(knon), INTENT(IN) :: phiT, phiQ, phiU, phiV REAL, DIMENSION(knon), INTENT(IN) :: phiQ0_b, phiT0_b REAL, DIMENSION(knon), INTENT(OUT) :: phiT_x, phiT_w REAL, DIMENSION(knon), INTENT(OUT) :: phiQ_x, phiQ_w REAL, DIMENSION(knon), INTENT(OUT) :: phiU_x, phiU_w REAL, DIMENSION(knon), INTENT(OUT) :: phiV_x, phiV_w REAL, DIMENSION(knon), INTENT(OUT) :: philat_x, philat_w REAL, DIMENSION(knon), INTENT(OUT) :: dqsatsrf ! beta delta(qsat(Ts)) REAL, DIMENSION(knon), INTENT(OUT) :: dTs ! Temperature difference at surface REAL, DIMENSION(knon), INTENT(OUT) :: delta_qsurf !! Local variables INTEGER :: j REAL, DIMENSION(knon) :: dphiT, dphiQ, dphiU, dphiV REAL, DIMENSION(knon) :: q1_x, q1_w REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region !---------------------------------------------------------------------------- ! Equations ! --------- !!!!!! (1 - g_T) dTs = dTs_ins + Gamma_phiT phiT !!!!!! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ !!!!!! dphiT = (dd_KTp/KTp) phiT + ( dd_AT - C_p dTs)*KxKwTp/KTp !!!!!! dphiQ = (dd_KQs/KQs) phiQ + (beta dd_AQ - dqsatsrf )*KxKwQs/KQs !!!!!! dphiU = (dd_KUp/KUp) phiU + ( dd_AU )*KxKwUp/KUp !!!!!! dphiV = (dd_KVp/KVp) phiV + ( dd_AV )*KxKwVp/KVp ! (1 - g_T) (dTs-dTs0) = dTs_ins-dTs0 + Gamma_phiT (phiT-phiT0) ! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ ! dphiT = (dd_KTp/KTp) phiT + ( dd_AT - C_p dTs)*KxKwTp/KTp ! dphiQ = (dd_KQs/KQs) phiQ + (beta dd_AQ - dqsatsrf )*KxKwQs/KQs ! dphiU = (dd_KUp/KUp) phiU + ( dd_AU )*KxKwUp/KUp ! dphiV = (dd_KVp/KVp) phiV + ( dd_AV )*KxKwVp/KVp !! sigx(:) = 1.-sigw(:) ! print *,' AAAA wx_pbl_split, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:) IF (iflag_split == 2 .AND. nsrf /= is_oce) THEN ! Delta_tsurf and Delta_qsurf computation ! ----------------------------------------- IF (prt_level >=10 ) THEN print *,' wx_pbl_split, dTs_ins, dTs0 , Gamma_phiT, g_T ', dTs_ins, dTs0, Gamma_phiT, g_T print *,' wx_pbl_split, dqsatsrf_ins, Gamma_phiQ, g_q ', dqsatsrf_ins, Gamma_phiQ, g_q ENDIF DO j = 1,knon dTs(j) = dTs0(j) + (dTs_ins(j) - dTs0(j) + Gamma_phiT(j)*(phiT(j)-phiT0_b(j)) )/(1 - g_T(j)) dqsatsrf(j) = dqsatsrf0(j) + (dqsatsrf_ins(j) - (beta(j)/QQ_b(j))*dTs0(j) + & Gamma_phiQ(j)*(phiQ(j)-phiQ0_b(j)) )/(1 - g_Q(j)) ENDDO ! j = 1,knon IF (prt_level >=10 ) THEN print *,' wx_pbl_split, dqsatsrf0, QQ_b ', dqsatsrf0, QQ_b print *,' wx_pbl_split, phiT0_b, phiT, dTs ', phiT0_b, phiT, dTs print *,' wx_pbl_split, phiQ0_b, phiQ, dqsatsrf ', phiQ0_b, phiQ, dqsatsrf ENDIF ELSE dTs(:) = 0. dqsatsrf(:) = 0. ENDIF ! (iflag_split .EQ. 2 .AND. nsrf .NE. is_oce) DO j = 1,knon dphiT(j) = (phiT(j)*dd_KTp(j) + ( dd_AT(j) - C_p(j)*dTs(j))*KxKwTp(j))/Kech_Tp(j) dphiQ(j) = (phiQ(j)*dd_KQs(j) + (beta(j)*dd_AQ(j) - dqsatsrf(j))*KxKwQs(j))/Kech_Qs(j) dphiU(j) = (phiU(j)*dd_KUp(j) + dd_AU(j) *KxKwUp(j))/Kech_Up(j) dphiV(j) = (phiV(j)*dd_KVp(j) + dd_AV(j) *KxKwVp(j))/Kech_Vp(j) phiT_x(j)=phiT(j) - sigw(j)*dphiT(j) phiT_w(j)=phiT(j) + sigx(j)*dphiT(j) phiQ_x(j)=phiQ(j) - sigw(j)*dphiQ(j) phiQ_w(j)=phiQ(j) + sigx(j)*dphiQ(j) phiU_x(j)=phiU(j) - sigw(j)*dphiU(j) phiU_w(j)=phiU(j) + sigx(j)*dphiU(j) phiV_x(j)=phiV(j) - sigw(j)*dphiV(j) phiV_w(j)=phiV(j) + sigx(j)*dphiV(j) philat_x(j)=phiQ_x(j)*RLVTT philat_w(j)=phiQ_w(j)*RLVTT ENDDO ! j = 1,knon DO j = 1,knon q1_x(j) = AQ_x(j) + BQ_x(j)*phiQ_x(j)*dtime q1_w(j) = AQ_w(j) + BQ_w(j)*phiQ_w(j)*dtime ENDDO ! j = 1,knon DO j = 1,knon delta_qsurf(j) = (1.-beta(j))*(q1_w(j) - q1_x(j)) + dqsatsrf(j) ENDDO ! j = 1,knon !! Do j = 1,knon !! 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) !! 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) !! ENDDO IF (prt_level >=10 ) THEN print *,' wx_pbl_split, phiT, dphiT, dTs ', phiT, dphiT, dTs print *,' wx_pbl_split, phiQ, dphiQ, dqsatsrf ', phiQ, dphiQ, dqsatsrf ENDIF IF (prt_level >=10 ) THEN !! print *,' wx_pbl_split, verif dqsatsrf = beta dqsatdT0 dTs ' !! print *,' wx_pbl_split, dqsatsrf, dqsatdT0*dTs ', dqsatsrf, dqsatdT0*dTs ENDIF !! IF (knon .NE. 0) THEN !! call iophys_ecrit('sigw', 1,'sigw', '.',sigw) !! call iophys_ecrit('phit', 1,'phit', 'W/m2',phit) !! call iophys_ecrit('phit_w', 1,'phit_w', 'W/m2',phit_w) !! call iophys_ecrit('phit_x', 1,'phit_x', 'W/m2',phit_x) !! call iophys_ecrit('phiq', 1,'phiq', 'kg/m2/s',phiq) !! call iophys_ecrit('phiq_w', 1,'phiq_w', 'kg/m2/s',phiq_w) !! call iophys_ecrit('phiq_x', 1,'phiq_x', 'kg/m2/s',phiq_x) !! call iophys_ecrit('q1_w', 1,'q1_w', '.',q1_w) !! call iophys_ecrit('q1_x', 1,'q1_x', '.',q1_x) !! ENDIF ! (knon .NE. 0) END SUBROUTINE wx_pbl_split SUBROUTINE wx_pbl_check( knon, dtime, ypplay, ypaprs, & sigw, beta, iflag_split, & Ts0_b9, dTs09, & qs_b9, Ts_b9, & ! yqsurf, Tsurf_new dTs9, dqsatsrf9, & AcoefT_x, AcoefT_w, & BcoefT_x, BcoefT_w, & AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, & AcoefT, AcoefQ, BcoefT, BcoefQ, & phiT_b9, phiQ_b9, & phiT_x9, phiT_w9, & phiQ_x9, phiQ_w9 & ) USE wx_pbl_var_mod USE lmdz_print_control, ONLY: prt_level,lunout INCLUDE "YOMCST.h" INCLUDE "FCTTRE.h" INCLUDE "YOETHF.h" INTEGER, INTENT(IN) :: knon ! number of grid cells REAL, INTENT(IN) :: dtime ! time step size (s) REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor INTEGER, INTENT(IN) :: iflag_split REAL, DIMENSION(knon), INTENT(IN) :: Ts0_b9, dTs09 REAL, DIMENSION(knon), INTENT(IN) :: qs_b9, Ts_b9 ! yqsurf, Tsurf_new REAL, DIMENSION(knon), INTENT(IN) :: dTs9, dqsatsrf9 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT, AcoefQ, BcoefT, BcoefQ REAL, DIMENSION(knon), INTENT(IN) :: phiT_b9, phiQ_b9 REAL, DIMENSION(knon), INTENT(IN) :: phiT_x9, phiT_w9 REAL, DIMENSION(knon), INTENT(IN) :: phiQ_x9, phiQ_w9 !! Local variables INTEGER :: j REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region REAL, DIMENSION(knon) :: AcoefT_b, AcoefQ_b ! mean values of AcoefT and AcoefQ REAL :: zzt, zzq, zzqsat REAL :: zdelta, zcvm5, zcor, qsat REAL, DIMENSION(knon) :: qsat_w, qsat_x REAL, DIMENSION(knon) :: dqsatdT_w, dqsatdT_x REAL, DIMENSION(knon) :: qsat_bs ! qsat(Ts_b) REAL, DIMENSION(knon) :: qsat01, dqsatdT01 REAL, DIMENSION(knon) :: Ts_x, Ts_w, qs_x, qs_w REAL, DIMENSION(knon) :: T1_x, T1_w, q1_x, q1_w REAL, DIMENSION(knon) :: Rn_x, Rn_w REAL, DIMENSION(knon) :: phiQ0_x, phiQ0_w REAL, DIMENSION(knon) :: Ta, qa REAL, DIMENSION(knon) :: qsatsrf_w, qsatsrf_x, qsatsrf_b REAL, DIMENSION(knon) :: qsurf_w, qsurf_x REAL :: dphiT, dphiQ REAL :: dqsatsrf1 REAL :: phiT_w1, phiT_w2 REAL :: phiT_x1, phiT_x2 REAL :: phiQ_w1, phiQ_w2, phiQ_w3 REAL :: phiQ_x1, phiQ_x2, phiQ_x3 REAL :: phiT_b1, phiQ_b1 REAL :: Kech_Q_sw1, Kech_Q_sx1 REAL :: evap_pot !---------------------------------------------------------------------------- ! Equations to be checked: ! ----------------------- ! Input : Ts0_b, dTs0, Ts_b, dTs, qsatsrf_b, dqsatsrf, ! phiT_b, phiQ_b, phiT_w, phiT_x, phiQ_w, phiQ_x, ! AcoefT, AcoefQ, AcoefT_w, AcoefQ_w, AcoefT_x, AcoefQ_x, ! BcoefT, BcoefQ, BcoefT_w, BcoefQ_w, BcoefT_x, BcoefQ_x ! C_p T1_w = AcoefT_w + BcoefT_w phiT_w Delta t C_p T1_x = AcoefT_x + BcoefT_x phiT_x Delta t ! q1_w = AQ_w + BQ_w phiQ_w Delta t q1_x = AQ_x + BQ_x phiQ_x Delta t ! qsatsrf_w = beta qsat(Ts_w) qsatsrf_x = beta qsat(Ts_x) ! qsurf_w = (1-beta) q1_w + qsatsrf_w qsurf_x = (1-beta) q1_x + qsatsrf_x ! phiT_w = Kech_h_w C_p ( T1_w - Ts_w) phiT_x = Kech_h_x C_p ( T1_x - Ts_x) ! phiT_w = Kech_T_pw ( AcoefT_w - C_p Ts_w) phiT_x = Kech_T_px ( AcoefT_x - C_p Ts_x) ! phiq_w = Kech_h_w ( beta q1_w - qsatsrf_w) phiq_x = Kech_h_x ( beta q1_x - qsatsrf_x)) ! phiq_w = Kech_Q_sw (beta AQ_w -qsatsrf_w) phiq_x = Kech_Q_sx (beta AQ_x -qsatsrf_x) ! phiq_w = Kech_h_w (q1_w - qsurf_w) phiq_x = Kech_h_x (q1_x - qsurf_x) ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x ! Ts_b = sigw Ts_w + sigx Ts_x dTs = Ts_w - Ts_x ! qsatsrf_b = sigw qsatsrf_w + sigx qsatsrf_x ! C_p Ta = AcoefT + BcoefT phiT_b Delta t ! qa = AcoefQ + BcoefQ phiQ_b Delta t ! phiT_b = Kech_h C_p (Ta - Ts_b) ! phiQ_b = beta Kech_h (qa - qsatsrf_b) ! dTs = sqrt(tau)/I (dphit + L_v dphiq + dR) !---------------------------------------------------------------------------- !! sigx(:) = 1.-sigw(:) AcoefT_b(1:knon) = AcoefT_x(1:knon) + sigw(1:knon)*dd_AT(1:knon) AcoefQ_b(1:knon) = AQ_x(1:knon) + sigw(1:knon)*dd_AQ(1:knon) ! Compute the three qsat and dqsatdTs ! --------------------------------------------- !! C_p(1:knon) = RCpd !! L_v(1:knon) = RLvtt IF (prt_level >=10 ) THEN print *,' AAAA wx_pbl_check, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:) ENDIF ! (prt_level >=10 ) DO j = 1, knon zdelta = MAX(0.,SIGN(1.,RTT-Ts0_b9(j))) zcvm5 = R5LES*(1.-zdelta) + R5IES*zdelta qsat = R2ES*FOEEW(Ts0_b9(j),zdelta)/ypaprs(j,1) qsat = MIN(0.5,qsat) zcor = 1./(1.-RETV*qsat) qsat01(j) = fqsat*qsat*zcor !! dqsatdT0(j) = FOEDE(Ts0_b(j),zdelta,zcvm5,qsat0(j),zcor)/RLVTT ! jyg 20210116 !! dqsatdT0(j) = (RLvtt*(1.-zdelta)+RLSTT*zdelta)*qsat0(j)/(Rv*Ts0_b(j)*Ts0_b(j)) dqsatdT01(j) = fqsat*FOEDE(Ts0_b9(j),zdelta,zcvm5,qsat01(j),zcor) ENDDO !-------------------------------------------------------------------------------------------------- IF (prt_level >=10 ) THEN DO j = 1, knon print *,'wx_pbl_check: Kech_h, Kech_q ', Kech_h(j), Kech_q(j) Ta(j) = (AcoefT(j) + BcoefT(j)*phiT_b9(j)*dtime)/C_p(j) qa(j) = AcoefQ(j) + BcoefQ(j)*phiQ_b9(j)*dtime print *, 'wx_pbl_check: j, Ta, qa ', Ta(j), qa(j) qsat_bs(j) = qsat01(j) + dqsatdT01(j)*(Ts_b9(j)-Ts0_b9(j)) print *,'wx_pbl_check: qsat01, qsat_bs ', j,qsat01(j), qsat_bs(j) Ts_x(j) = Ts_b9(j) - sigw(j)*dTs9(j) Ts_w(j) = Ts_b9(j) + sigx(j)*dTs9(j) print *, 'wx_pbl_check: j, Ts_b9, Ts_w, Ts_x ', j, Ts_b9(j), Ts_w(j), Ts_x(j) qsat_x(j) = qsat0_x(j) + dqsatdT0_x(j)*(Ts_x(j)-Ts0_x(j)) qsat_w(j) = qsat0_w(j) + dqsatdT0_w(j)*(Ts_w(j)-Ts0_w(j)) print *,'wx_pbl_check: qsat0_w, qsat0_x, qsat_w, qsat_x ', qsat0_w(j), qsat0_x(j), qsat_w(j), qsat_x(j) T1_x(j) = (AcoefT_x(j) + BcoefT_x(j)*phiT_x9(j)*dtime) / C_p(j) T1_w(j) = (AcoefT_w(j) + BcoefT_w(j)*phiT_w9(j)*dtime) / C_p(j) print *, 'wx_pbl_check: j, T1_w, T1_x ', j, T1_w(j), T1_x(j) q1_x(j) = AQ_x(j) + BQ_x(j)*phiQ_x9(j)*dtime q1_w(j) = AQ_w(j) + BQ_w(j)*phiQ_w9(j)*dtime print *, 'wx_pbl_check: j, q1_w, q1_x ', j, q1_w(j), q1_x(j) qsatsrf_x(j) = beta(j)*qsat_x(j) qsatsrf_w(j) = beta(j)*qsat_w(j) qsatsrf_b(j) = sigw(j)*qsatsrf_w(j) + sigx(j)*qsatsrf_x(j) dqsatsrf1 = qsatsrf_w(j) - qsatsrf_x(j) print *, 'wx_pbl_check: j, qsatsrf_w, qsatsrf_x, dqsatsrf1, dqsatsrf9 ', & qsatsrf_w(j), qsatsrf_x(j), dqsatsrf1, dqsatsrf9(j) qsurf_x(j) = (1-beta(j))*q1_x(j) + qsatsrf_x(j) qsurf_w(j) = (1-beta(j))*q1_w(j) + qsatsrf_w(j) print *, 'wx_pbl_check: j, qsurf_w, qsurf_x ', j, qsurf_w(j), qsurf_x(j) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test qsat01 = qsat0 et dqsatdT01 = dqsatdT0 !------------------------------------------------------------------------------------------------------ print *, 'wx_pbl_check: j, qsat01(j), qsat0(j) ', j, qsat01(j), qsat0(j) print *, 'wx_pbl_check: j, dqsatdT01(j), dqsatdT0(j) ', j, dqsatdT01(j), dqsatdT0(j) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! 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) !------------------------------------------------------------------------------------------------------ Kech_Q_sx1 = Kech_q_x(j)/(1.-beta(j)*Kech_q_x(j)*BQ_x(j)*dtime) Kech_Q_sw1 = Kech_q_w(j)/(1.-beta(j)*Kech_q_w(j)*BQ_w(j)*dtime) print *, 'wx_pbl_check: j, Kech_Q_sx1, Kech_Q_sx(j)', j, Kech_Q_sx1, Kech_Q_sx(j) print *, 'wx_pbl_check: j, Kech_Q_sw1, Kech_Q_sw(j)', j, Kech_Q_sw1, Kech_Q_sw(j) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! 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)) !----------------------------------------------------- phiT_x1 = Kech_h_x(j)*C_p(j)*(T1_x(j)-Ts_x(j)) phiT_w1 = Kech_h_w(j)*C_p(j)*(T1_w(j)-Ts_w(j)) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! 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)) !----------------------------------------------------- phiT_x2 = Kech_T_px(j)*(AcoefT_x(j)-C_p(j)*Ts_x(j)) phiT_w2 = Kech_T_pw(j)*(AcoefT_w(j)-C_p(j)*Ts_w(j)) print *, 'wx_pbl_check: j, phiT_w1, phiT_w2, phiT_w9 ', j, phiT_w1, phiT_w2, phiT_w9(j) print *, 'wx_pbl_check: j, phiT_x1, phiT_x2, phiT_x9 ', j, phiT_x1, phiT_x2, phiT_x9(j) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test phiq_w = Kech_q_w ( beta q1_w - qsatsrf_w) phiq_x = Kech_q_x ( beta q1_x - qsatsrf_x)) !-------------------------------------------------------------- phiq_x1 = Kech_q_x(j)*( beta(j)*q1_x(j) - qsatsrf_x(j)) phiq_w1 = Kech_q_w(j)*( beta(j)*q1_w(j) - qsatsrf_w(j)) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test phiq_w = Kech_Q_sw (beta AQ_w -qsatsrf_w) phiq_x = Kech_Q_sx (beta AQ_x -qsatsrf_x) !-------------------------------------------------------------- phiq_x2 = Kech_Q_sx(j)*(beta(j)*AQ_x(j) -qsatsrf_x(j)) phiq_w2 = Kech_Q_sw(j)*(beta(j)*AQ_w(j) -qsatsrf_w(j)) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test phiq_w = Kech_q_w ( q1_w - qsurf_w) phiq_x = Kech_q_x ( q1_x - qsurf_x)) !-------------------------------------------------------------- phiq_x3 = Kech_q_x(j)*( q1_x(j) - qsurf_x(j)) phiq_w3 = Kech_q_w(j)*( q1_w(j) - qsurf_w(j)) print *, 'wx_pbl_check: j, phiQ_w1, phiQ_w2, phiQ_w3, phiQ_w9 ', j, phiQ_w1, phiQ_w2, phiQ_w3, phiQ_w9(j) print *, 'wx_pbl_check: j, phiQ_x1, phiQ_x2, phiQ_x3, phiQ_x9 ', j, phiQ_x1, phiQ_x2, phiQ_x3, phiQ_x9(j) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test phiT_b = Kech_h C_p (Ta - Ts_b) !-------------------------------------------------------------- phiT_b1 = Kech_h(j)*C_p(j)*(Ta(j) - Ts_b9(j)) print *, 'wx_pbl_check: j, phiT_b1, PhiT_b9 ', j, phiT_b1, PhiT_b9(j) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test phiQ_b = beta Kech_q (qa - qsat_bs) !-------------------------------------------------------------- evap_pot = Kech_q(j)*(qa(j) - qsat_bs(j)) phiQ_b1 = beta(j)*Kech_q(j)*(qa(j) - qsat_bs(j)) print *, 'wx_pbl_check: j, beta, evap_pot, phiQ_b1, PhiQ_b9 ', j, beta(j), evap_pot, phiQ_b1, PhiQ_b9(j) ENDDO ! j = 1, knon ENDIF ! (prt_level >=10 ) !-------------------------------------------------------------------------------------------------- END SUBROUTINE wx_pbl_check SUBROUTINE wx_pbl_dts_check( knon, dtime, ypplay, ypaprs, & sigw, beta, iflag_split, & Ts0_b9, dTs09, & qs_b9, Ts_b9, & ! yqsurf, Tsurf_new dqsatsrf9, dTs9, delta_qsurf9, & AcoefT_x, AcoefT_w, & BcoefT_x, BcoefT_w, & AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, & AcoefT, AcoefQ, BcoefT, BcoefQ, & HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, & phiT0_b9, dphiT09, phiQ0_b9, dphiQ09, Rn0_b9, dRn09, & g_T, g_Q, & Gamma_phiT, Gamma_phiQ, & dTs_ins, dqsatsrf_ins, & phiT_b9, phiQ_b9, & phiT_x9, phiT_w9, & phiQ_x9, phiQ_w9 & ) USE wx_pbl_var_mod USE lmdz_print_control, ONLY: prt_level,lunout INCLUDE "YOMCST.h" INCLUDE "FCTTRE.h" INCLUDE "YOETHF.h" INTEGER, INTENT(IN) :: knon ! number of grid cells REAL, INTENT(IN) :: dtime ! time step size (s) REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa) REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa) REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor INTEGER, INTENT(IN) :: iflag_split REAL, DIMENSION(knon), INTENT(IN) :: Ts0_b9, dTs09 REAL, DIMENSION(knon), INTENT(IN) :: qs_b9, Ts_b9 ! yqsurf, Tsurf_new REAL, DIMENSION(knon), INTENT(IN) :: dTs9, dqsatsrf9 REAL, DIMENSION(knon), INTENT(IN) :: delta_qsurf9 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT, AcoefQ, BcoefT, BcoefQ REAL, DIMENSION(knon), INTENT(IN) :: HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn REAL, DIMENSION(knon), INTENT(IN) :: phiT0_b9, dphiT09, phiQ0_b9, dphiQ09, Rn0_b9, dRn09 REAL, DIMENSION(knon), INTENT(IN) :: g_T, g_Q REAL, DIMENSION(knon), INTENT(IN) :: Gamma_phiT, Gamma_phiQ REAL, DIMENSION(knon), INTENT(IN) :: dTs_ins, dqsatsrf_ins REAL, DIMENSION(knon), INTENT(IN) :: phiT_b9, phiQ_b9 REAL, DIMENSION(knon), INTENT(IN) :: phiT_x9, phiT_w9 REAL, DIMENSION(knon), INTENT(IN) :: phiQ_x9, phiQ_w9 !! Local variables INTEGER :: j REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region REAL, DIMENSION(knon) :: AcoefT_b, AcoefQ_b ! mean values of AcoefT and AcoefQ REAL :: zzt, zzq, zzqsat REAL :: zdelta, zcvm5, zcor, qsat REAL, DIMENSION(knon) :: qsat_w, qsat_x REAL, DIMENSION(knon) :: Ts_x, Ts_w, qs_x, qs_w REAL, DIMENSION(knon) :: T1_x, T1_w, q1_x, q1_w REAL, DIMENSION(knon) :: Rn_x, Rn_w REAL, DIMENSION(knon) :: Rn_b, dRn REAL, DIMENSION(knon) :: phiQ0_x, phiQ0_w REAL, DIMENSION(knon) :: Ta, qa REAL, DIMENSION(knon) :: err_phiT_w, err_phiT_x REAL, DIMENSION(knon) :: err_phiq_w, err_phiq_x REAL, DIMENSION(knon) :: err_phiT_b REAL, DIMENSION(knon) :: err_phiQ_b REAL, DIMENSION(knon) :: err2_phiT_b REAL :: T1A_x, T1A_w, q1A_x, q1A_w REAL :: qsatsrf_w, qsatsrf_x, qsatsrfb, qsbA REAL :: dphiT, dphiQ REAL :: dphiT_H, dphiQ_H REAL :: phiQ_pot REAL :: phiQ_w_m_phiQ0_w REAL :: phiQ_x_m_phiQ0_x REAL :: dphiQ_m_dphiQ0 REAL :: dphiT_m_dphiT0 REAL :: dRN_m_dRn0 REAL :: phiTb_m_phiT0b !---------------------------------------------------------------------------- ! Equations to be checked: ! ----------------------- ! Input : Ts0_b, dTs0, Ts_b, dTs, qsatsrf_b, dqsatsrf, ! phiT_b, phiQ_b, phiT_w, phiT_x, phiQ_w, phiQ_x, ! AcoefT, AcoefQ, AcoefT_w, AcoefQ_w, AcoefT_x, AcoefQ_x, ! BcoefT, BcoefQ, BcoefT_w, BcoefQ_w, BcoefT_x, BcoefQ_x ! Ts_w = Ts_b + sigx dTs Ts_x = Ts_b - sigw dTs ! T1_w = AcoefT_w + BcoefT_w phiT_w Delta t T1_x = AcoefT_x + BcoefT_x phiT_x Delta t ! q1_w = AcoefQ_w + BcoefQ_w phiQ_w Delta t q1_x = AcoefQ_x + BcoefQ_x phiQ_x Delta t ! phiT_w = Kech_h_w ( T1_w - Ts_w) phiT_x = Kech_h_x ( T1_x - Ts_x) ! phiq_w = beta Kech_h_w ( q1_w - qsat(Ts_w)) phiq_x = beta Kech_h_x ( q1_x - qsat(Ts_x)) ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x ! Ts_b = sigw Ts_w + sigx Ts_x dTs = Ts_w - Ts_x ! Ta = AcoefT + BcoefT phiT_b Delta t ! qa = AcoefQ + BcoefQ phiQ_b Delta t ! phiT_b = Kech_h (Ta - Ts_b) ! phiQ_b = beta Kech_h (qa - qsat(Ts_b)) ! dTs = sqrt(tau)/I (dphit + L_v dphiq + dR) !---------------------------------------------------------------------------- !! sigx(:) = 1.-sigw(:) AcoefT_b(1:knon) = AcoefT_x(1:knon) + sigw(1:knon)*dd_AT(1:knon) AcoefQ_b(1:knon) = AQ_x(1:knon) + sigw(1:knon)*dd_AQ(1:knon) IF (prt_level >=10 ) THEN print *,'->wx_pbl_dts_check, HTphiT_b, HTphiQ_b, HTRn_b ', & HTphiT_b, HTphiQ_b, HTRn_b print *,'->wx_pbl_dts_check, dd_HTphiT, dd_HTphiQ, dd_HTRn ', & dd_HTphiT, dd_HTphiQ, dd_HTRn ENDIF ! (prt_level >=10 ) ! Compute the three qsat and dqsatdTs ! --------------------------------------------- !! print *,' AAAA wx_pbl_dts_check, C_p(j), qsat0(j), Ts0(j) : ', & !! (C_p(j), qsat0(j), Ts0(j), j = 1,knon) !-------------------------------------------------------------------------------------------------- IF (prt_level >=10 ) THEN DO j = 1, knon Ts_x(j) = Ts_b9(j) - sigw(j)*dTs9(j) Ts_w(j) = Ts_b9(j) + sigx(j)*dTs9(j) print *, 'wx_pbl_dts_check: j, Ts_b9, Ts_w, Ts_x ', j, Ts_b9(j), Ts_w(j), Ts_x(j) qsat_x(j) = qsat0_x(j) + dqsatdT0_x(j)*(Ts_x(j)-Ts0_x(j)) qsat_w(j) = qsat0_w(j) + dqsatdT0_w(j)*(Ts_w(j)-Ts0_w(j)) T1_x(j) = (AcoefT_x(j) + BcoefT_x(j)*phiT_x9(j)*dtime) / C_p(j) T1_w(j) = (AcoefT_w(j) + BcoefT_w(j)*phiT_w9(j)*dtime) / C_p(j) print *, 'wx_pbl_dts_check: j, T1_w, T1_x ', j, T1_w(j), T1_x(j) q1_x(j) = AQ_x(j) + BQ_x(j)*phiQ_x9(j)*dtime q1_w(j) = AQ_w(j) + BQ_w(j)*phiQ_w9(j)*dtime print *, 'wx_pbl_dts_check: j, q1_w, q1_x ', j, q1_w(j), q1_x(j) Rn_x(j) = eps_1*Rsigma*T1_x(j)**4 - Rsigma*Ts_x(j)**4 Rn_w(j) = eps_1*Rsigma*T1_w(j)**4 - Rsigma*Ts_w(j)**4 Rn_b(j) = sigw(j)*Rn_w(j) + sigx(j)*Rn_x(j) dRn(j) = dRn09(j) - ( HTRn_b(j) & +(sigx(j)-sigw(j))*dd_HTRn(j) & -sigw(j)*sigx(j)*dd_HTRn(j)*dd_HTphiT(j)/HTphiT_b(j) & )*(dTs9(j)-dTs09(j)) & + dd_HTRn(j)/HTphiT_b(j)*(phiT_b9(j)-phiT0_b9(j)) print *,'wx_pbl_dts_check, dphiT, L_v*dphiQ, dRn, dTs ', & phiT_w9(j)-phiT_x9(j), L_v(j)*(phiQ_w9(j)-phiQ_x9(j)), dRn(j), dTs9(j) phiQ0_x(j) = PhiQ0_b9(j) - sigw(j)*dphiQ09(j) phiQ0_w(j) = PhiQ0_b9(j) + sigx(j)*dphiQ09(j) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test phiQ_w-phiQ0_w = -beta*Kech_Q_sw*dqsatdT_w*(Ts_w-Ts0_w) !-------------------------------------------------------------- print *,'wx_pbl_dts_check: beta(j), Kech_Q_sw(j), dqsatdT0_w(j), Ts_w(j), Ts0_w(j) ', & beta(j), Kech_Q_sw(j), dqsatdT0_w(j), Ts_w(j), Ts0_w(j) phiQ_w_m_phiQ0_w = -beta(j)*Kech_Q_sw(j)*dqsatdT0_w(j)*(Ts_w(j)-Ts0_w(j)) print *,'wx_pbl_dts_check: j, phiQ_w9-phiQ0_w, phiQ_w_m_phiQ0_w ', & j, phiQ_w9(j)-phiQ0_w(j), phiQ_w_m_phiQ0_w !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test phiQ_x-phiQ0_x = -beta*Kech_Q_sx*dqsatdT_x*(Ts_x-Ts0_x) !-------------------------------------------------------------- phiQ_x_m_phiQ0_x = -beta(j)*Kech_Q_sx(j)*dqsatdT0_x(j)*(Ts_x(j)-Ts0_x(j)) print *,'wx_pbl_dts_check: j, phiQ_x9-phiQ0_x, phiQ_x_m_phiQ0_x ', & j, phiQ_x9(j)-phiQ0_x(j), phiQ_x_m_phiQ0_x !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test dphiT-dphiT0 = -(HTphiT_b+(sigx-sigw)*dd_HTphiT)*(dTs-dTs0) - dd_HTphiT*(Ts_b-Ts0_b) !------------------------------------------------------------------------------------------- dphiT = phiT_w9(j) - phiT_x9(j) dphiT_m_dphiT0 = -(HTphiT_b(j)+(sigx(j)-sigw(j))*dd_HTphiT(j))*(dTs9(j)-dTs09(j)) & - dd_HTphiT(j)*(Ts_b9(j)-Ts0_b9(j)) print *,'wx_pbl_dts_check: j, dphiT-dphiT09, dphiT_m_dphiT0 ',j, dphiT-dphiT09(j), dphiT_m_dphiT0 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test dphiQ-dphiQ0 = -(HTphiQ_b+(sigx-sigw)*dd_HTphiQ)*(dTs-dTs0) - dd_HTphiQ*(Ts_b-Ts0_b) !------------------------------------------------------------------------------------------- dphiQ = phiQ_w9(j) - phiQ_x9(j) dphiQ_m_dphiQ0 = -(HTphiQ_b(j)+(sigx(j)-sigw(j))*dd_HTphiQ(j))*(dTs9(j)-dTs09(j)) & - dd_HTphiQ(j)*(Ts_b9(j)-Ts0_b9(j)) print *,'wx_pbl_dts_check: j, dphiQ-dphiQ09, dphiQ_m_dphiQ0 ',j, dphiQ-dphiQ09(j), dphiQ_m_dphiQ0 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test dRn-dRn0 = -(HTRn_b+(sigx-sigw)*dd_HTRn)*(dTs-dTs0) - dd_HTRn*(Ts_b-Ts0_b) !------------------------------------------------------------------------------------------- dRn_m_dRn0 = -(HTRn_b(j)+(sigx(j)-sigw(j))*dd_HTRn(j))*(dTs9(j)-dTs09(j)) & - dd_HTRn(j)*(Ts_b9(j)-Ts0_b9(j)) print *,'wx_pbl_dts_check: j, dRn-dRn09, dRn_m_dRn0 ',j, dRn-dRn09(j), dRn_m_dRn0 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test phiT_b-phiT0_b = -sigx*sigw*dd_HTphiT*(dTs-dTs0) - HTphiT_b*(Ts_b-Ts0_b) !------------------------------------------------------------------------------- phiTb_m_phiT0b = -sigx(j)*sigw(j)*dd_HTphiT(j)*(dTs9(j)-dTs09(j)) - HTphiT_b(j)*(Ts_b9(j)-Ts0_b9(j)) print *,'wx_pbl_dts_check: j, phiT_b9-phiT0_b9, phiTb_m_phiT0b ',j ,phiT_b9(j)-phiT0_b9(j), phiTb_m_phiT0b !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test phiT_w, phiT_x, dphiT from HTphiT !------------------------------------------ ! phiT_w = Kech_h_w C_p ( T1_w - Ts_w) phiT_x = Kech_h_x C_p ( T1_x - Ts_x) err_phiT_x(j) = Kech_h_x(j)*C_p(j)*(T1_x(j) - Ts_x(j)) - phiT_x9(j) err_phiT_w(j) = Kech_h_w(j)*C_p(j)*(T1_w(j) - Ts_w(j)) - phiT_w9(j) print *, 'wx_pbl_dts_check: j, phiT_w9, phiT_x9, err_phiT_w, err_phiT_x ', & j, phiT_w9(j), phiT_x9(j), err_phiT_w(j), err_phiT_x(j) dphiT = phiT_w9(j) - phiT_x9(j) dphiT_H = dphiT09(j) - ( HTphiT_b(j) & +(sigx(j)-sigw(j))*dd_HTphiT(j) & -sigw(j)*sigx(j)*dd_HTphiT(j)*dd_HTphiT(j)/HTphiT_b(j) & )*(dTs9(j)-dTs09(j)) & + dd_HTphiT(j)/HTphiT_b(j)*(phiT_b9(j)-phiT0_b9(j)) print *,'wx_pbl_dts_check: j, dphiT, dphiT_H ', j, dphiT, dphiT_H !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Test phiq_w, phiq_x, dphiq from HTphiq !------------------------------------------ ! phiq_w = beta Kech_q_w ( q1_w - qsat(Ts_w)) phiq_x = beta Kech_q_x ( q1_x - qsat(Ts_x)) err_phiq_x(j) = beta(j)*Kech_q_x(j)*( q1_x(j) - qsat_x(j)) - phiq_x9(j) err_phiq_w(j) = beta(j)*Kech_q_w(j)*( q1_w(j) - qsat_w(j)) - phiq_w9(j) dphiQ = phiQ_w9(j) - phiQ_x9(j) dphiQ_H = dphiQ09(j) - ( HTphiQ_b(j) & +(sigx(j)-sigw(j))*dd_HTphiQ(j) & -sigw(j)*sigx(j)*dd_HTphiQ(j)*dd_HTphiT(j)/HTphiT_b(j) & )*(dTs9(j)-dTs09(j)) & + dd_HTphiQ(j)/HTphiT_b(j)*(phiT_b9(j)-phiT0_b9(j)) print *,'wx_pbl_dts_check: j, dphiQ, dphiQ_H ', j, dphiQ, dphiQ_H ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x err_phiT_b(j) = sigw(j)*phiT_w9(j) + sigx(j)*phiT_x9(j) - phiT_b9(j) ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x err_phiQ_b(j) = sigw(j)*phiQ_w9(j) + sigx(j)*phiQ_x9(j) - phiQ_b9(j) ! Ta = AcoefT + BcoefT phiT_b Delta t ! phiT_b = Kech_h C_p (Ta - Ts_b) Ta(j) = (AcoefT(j) + BcoefT(j)*phiT_b9(j)*dtime) / C_p(j) err2_phiT_b(j) = Kech_h(j)*C_p(j)*(Ta(j) - Ts_b9(j)) - phiT_b9(j) print *, 'wx_pbl_dts_check: j, Ta, phiT_b9, err2_phiT_b ', & j, Ta(j), phiT_b9(j), err2_phiT_b(j) ENDDO ! j = 1, knon ENDIF ! (prt_level >=10 ) !-------------------------------------------------------------------------------------------------- END SUBROUTINE wx_pbl_dts_check SUBROUTINE wx_evappot(knon, q1, Ts, evap_pot) USE wx_pbl_var_mod INTEGER, INTENT(IN) :: knon ! number of grid cells REAL, DIMENSION(knon), INTENT(IN) :: q1 ! specific humidity in layer 1 REAL, DIMENSION(knon), INTENT(IN) :: Ts ! surface temperature REAL, DIMENSION(knon), INTENT(OUT) :: evap_pot ! potential evaporation INTEGER :: j REAL :: qsat_bs DO j = 1,knon evap_pot(j) = Kech_q(j)*(qsat0(j)+dqsatdT0(j)*(Ts(j)-Ts0(j))-q1(j)) qsat_bs = qsat0(j)+dqsatdT0(j)*(Ts(j)-Ts0(j)) !! print *,'wx_evappot : Kech_q, qsat_bs, qa, evap_pot ', Kech_q(j), qsat_bs, q1(j), evap_pot(j) ENDDO END SUBROUTINE wx_evappot END MODULE wx_pbl_mod