1 | subroutine thermcell_qsat(klon,active,pplev,ztemp,zqta,zqsat) |
---|
2 | implicit none |
---|
3 | |
---|
4 | #include "YOMCST.h" |
---|
5 | #include "YOETHF.h" |
---|
6 | #include "FCTTRE.h" |
---|
7 | |
---|
8 | |
---|
9 | !==================================================================== |
---|
10 | ! DECLARATIONS |
---|
11 | !==================================================================== |
---|
12 | |
---|
13 | ! Arguments |
---|
14 | INTEGER klon |
---|
15 | REAL zpspsk(klon),pplev(klon) |
---|
16 | REAL ztemp(klon),zqta(klon),zqsat(klon) |
---|
17 | LOGICAL active(klon) |
---|
18 | |
---|
19 | ! Variables locales |
---|
20 | INTEGER ig,iter |
---|
21 | REAL Tbef(klon),DT(klon) |
---|
22 | REAL tdelta,qsatbef,zcor,qlbef,zdelta,zcvm5,dqsat,num,denom,dqsat_dT |
---|
23 | logical Zsat |
---|
24 | REAL RLvCp |
---|
25 | REAL, SAVE :: DDT0=.01 |
---|
26 | LOGICAL afaire(klon),tout_converge |
---|
27 | |
---|
28 | !==================================================================== |
---|
29 | ! INITIALISATIONS |
---|
30 | !==================================================================== |
---|
31 | |
---|
32 | RLvCp = RLVTT/RCPD |
---|
33 | tout_converge=.false. |
---|
34 | afaire(:)=.false. |
---|
35 | DT(:)=0. |
---|
36 | |
---|
37 | |
---|
38 | !==================================================================== |
---|
39 | ! Routine a vectoriser en copiant active dans converge et en mettant |
---|
40 | ! la boucle sur les iterations a l'exterieur est en mettant |
---|
41 | ! converge= false des que la convergence est atteinte. |
---|
42 | !==================================================================== |
---|
43 | |
---|
44 | do ig=1,klon |
---|
45 | if (active(ig)) then |
---|
46 | Tbef(ig)=ztemp(ig) |
---|
47 | zdelta=MAX(0.,SIGN(1.,RTT-Tbef(ig))) |
---|
48 | qsatbef= R2ES * FOEEW(Tbef(ig),zdelta)/pplev(ig) |
---|
49 | qsatbef=MIN(0.5,qsatbef) |
---|
50 | zcor=1./(1.-retv*qsatbef) |
---|
51 | qsatbef=qsatbef*zcor |
---|
52 | qlbef=max(0.,zqta(ig)-qsatbef) |
---|
53 | DT(ig) = 0.5*RLvCp*qlbef |
---|
54 | zqsat(ig)=qsatbef |
---|
55 | endif |
---|
56 | enddo |
---|
57 | |
---|
58 | ! Traitement du cas ou il y a condensation mais faible |
---|
59 | ! On ne condense pas mais on dit que le qsat est le qta |
---|
60 | do ig=1,klon |
---|
61 | if (active(ig)) then |
---|
62 | if (0.<abs(DT(ig)).and.abs(DT(ig))<=DDT0) then |
---|
63 | zqsat(ig)=zqta(ig) |
---|
64 | endif |
---|
65 | endif |
---|
66 | enddo |
---|
67 | |
---|
68 | do iter=1,10 |
---|
69 | afaire(:)=abs(DT(:)).gt.DDT0 |
---|
70 | do ig=1,klon |
---|
71 | if (afaire(ig)) then |
---|
72 | Tbef(ig)=Tbef(ig)+DT(ig) |
---|
73 | zdelta=MAX(0.,SIGN(1.,RTT-Tbef(ig))) |
---|
74 | qsatbef= R2ES * FOEEW(Tbef(ig),zdelta)/pplev(ig) |
---|
75 | qsatbef=MIN(0.5,qsatbef) |
---|
76 | zcor=1./(1.-retv*qsatbef) |
---|
77 | qsatbef=qsatbef*zcor |
---|
78 | qlbef=zqta(ig)-qsatbef |
---|
79 | zdelta=MAX(0.,SIGN(1.,RTT-Tbef(ig))) |
---|
80 | zcvm5=R5LES*(1.-zdelta) + R5IES*zdelta |
---|
81 | zcor=1./(1.-retv*qsatbef) |
---|
82 | dqsat_dT=FOEDE(Tbef(ig),zdelta,zcvm5,qsatbef,zcor) |
---|
83 | num=-Tbef(ig)+ztemp(ig)+RLvCp*qlbef |
---|
84 | denom=1.+RLvCp*dqsat_dT |
---|
85 | zqsat(ig) = qsatbef |
---|
86 | DT(ig)=num/denom |
---|
87 | endif |
---|
88 | enddo |
---|
89 | enddo |
---|
90 | |
---|
91 | return |
---|
92 | end |
---|