da_condens_lin.inc
References to this file elsewhere.
1 subroutine da_condens_lin(DT,SCR31,SCR42,SCR71,DUM31,PRD, &
2 QVT,QCT,QRT,TTT,P_B,T_B,QV_B,QCW_B,QRN_B, &
3 SCR319,SCR429,SCR719,DUM319,PRD9, &
4 QVT9,QCT9,QRT9,TTT9,P_A,T_A,QV_A,QCW_A,QRN_A,kts,kte)
5
6 !-----------------------------------------------------------------------
7 ! Purpose: Condensation
8 !-----------------------------------------------------------------------
9
10 implicit none
11
12 integer, intent(in) :: kts, kte
13 ! real :: SVP1,SVP2,SVP3,SVPT0,TO,gas_constant_v,XLS
14 real, dimension(kts:kte), intent(in) :: DT,SCR31,SCR42,SCR71,PRD,DUM31
15 real, dimension(kts:kte), intent(in) :: P_B,T_B,QV_B,QCW_B,QRN_B
16 real, dimension(kts:kte), intent(in) :: SCR319,SCR429,SCR719,PRD9,DUM319
17 real, dimension(kts:kte), intent(in) :: P_A
18 real, dimension(kts:kte), intent(inout) :: T_A,QV_A,QCW_A,QRN_A
19
20 real, dimension(kts:kte), intent(in) :: QVT,QCT,QRT,TTT
21 real, dimension(kts:kte), intent(in) :: QVT9,QCT9,QRT9,TTT9
22
23 real, dimension(kts:kte) :: TMP,PRC59,DUM2139,PRC5,DUM213,SCR619,SCR89,SCR8
24 real, dimension(kts:kte) :: SCR61
25 real, dimension(kts:kte) :: DUM114,DUM1149,DUM2129,DUM115,DUM212,DUM1159
26
27 integer :: k
28
29 do K=kts, kte
30
31 if (DT(k) <= 0.0) cycle
32
33 DUM114(K)=1.0E3*SVP1*EXP(SVP2*(SCR71(K)-SVPT0)/(SCR71(K)-SVP3))
34 DUM1149(k)=DUM114(K)*SVP2*(SVPT0-SVP3)/(SCR71(K)-SVP3)**2*SCR719(K)
35
36 if(SCR71(K) > TO) then
37 DUM2129(K)=2.0*DUM31(K)/(gas_constant_v*PRD(K))*DUM319(K) &
38 -DUM31(K)*DUM31(K)/(gas_constant_v*PRD(K)*PRD(K))*PRD9(K)
39 DUM212(K)=DUM31(K)*DUM31(K)/(gas_constant_v*PRD(K))
40 else
41 DUM2129(K)=XLS*DUM319(K)/(gas_constant_v*PRD(K)) &
42 -XLS*DUM31(K)/(gas_constant_v*PRD(K)*PRD(K))*PRD9(K)
43 DUM212(K)=XLS*DUM31(K)/(gas_constant_v*PRD(K))
44 end if
45 TMP(K)=0.622/(P_B(K)-DUM114(K))**2
46 PRC59(K)=TMP(K)*P_B(K)*DUM1149(K)-TMP(K)*DUM114(K)*P_A(K)
47 PRC5(K)=0.622*DUM114(K)/(P_B(K)-DUM114(K))
48
49 if(SCR42(K) < PRC5(K) .AND. SCR71(K) < TO) then
50 SCR619(K)=0.0
51 SCR61(K)=0.0
52 else
53 TMP(K)=1./(1.0+DUM212(K)*PRC5(K)/(SCR71(K)*SCR71(K)))
54 SCR89(K)=TMP(K)*SCR429(K) &
55 -TMP(K)*(1.0+(SCR42(K)-PRC5(K))*DUM212(K)/ &
56 (SCR71(K)*SCR71(K))*TMP(K))*PRC59(K) &
57 -TMP(K)*TMP(K)*(SCR42(K)-PRC5(K))*PRC5(K)/ &
58 (SCR71(K)*SCR71(K))*DUM2129(K) &
59 ! WHY?
60 !error -TMP(K)*TMP(K)*2.0*DUM212(K)*PRC5(K) &
61 !error *(SCR42(K)-PRC5(K))/(SCR71(K)*SCR71(K))*SCR719(K)
62 +TMP(K)*TMP(K)*2.0*DUM212(K)*PRC5(K) &
63 *(SCR42(K)-PRC5(K))/(SCR71(K)*SCR71(K)*SCR71(K))*SCR719(K)
64 SCR8(K)=(SCR42(K)-PRC5(K))/(1.0+DUM212(K)*PRC5(K)/ &
65 (SCR71(K)*SCR71(K)))
66
67 DUM1159(K)=SCR319(K)+SCR89(K)
68 DUM115(K)=SCR31(K)+SCR8(K)
69 if(DUM115(K) >= 0.0)then
70 SCR619(K)=SCR89(K)/DT(k)
71 SCR61(K)=SCR8(K)/DT(k)
72 else
73 SCR619(K)=-SCR319(K)/DT(k)
74 SCR61(K)=-SCR31(K)/DT(k)
75 end if
76 end if
77 QV_A(K)=QV_A(K)+(QVT9(K)-SCR619(K))*DT(K)
78 if(QV_B(K) < 1.0E-25) QV_A(K)=0.0
79 QCW_A(K)=QCW_A(K)+(QCT9(K)+SCR619(K))*DT(K)
80 if(QCW_B(K) < 1.0E-25) QCW_A(K)=0.0
81 if(SCR71(K) > TO)then
82 DUM2139(K)=DUM319(K)/PRD(K)-DUM31(K)/(PRD(K)*PRD(K))*PRD9(K)
83 DUM213(K)=DUM31(K)/PRD(K)
84 else
85 DUM2139(K)=-XLS/(PRD(K)*PRD(K))*PRD9(K)
86 DUM213(K)=XLS/PRD(K)
87 end if
88
89 QRN_A(K)=QRN_A(K)+DT(K)*QRT9(K)
90 if(QRN_B(K) < 1.0E-25) QRN_A(K)=0.0
91 T_A(K)=T_A(K)+DT(K)*(TTT9(K)+SCR619(K)*DUM213(K)+SCR61(K)*DUM2139(K))
92
93 end do
94
95 end subroutine da_condens_lin