1 | MODULE lmdz_thermcell_qsat |
---|
2 | CONTAINS |
---|
3 | |
---|
4 | SUBROUTINE thermcell_qsat(klon, active, pplev, ztemp, zqta, zqsat) |
---|
5 | USE lmdz_yoethf |
---|
6 | USE lmdz_fcttre, ONLY: foeew, foede, qsats, qsatl, dqsats, dqsatl, thermcep |
---|
7 | USE lmdz_yomcst |
---|
8 | |
---|
9 | IMPLICIT NONE |
---|
10 | |
---|
11 | !==================================================================== |
---|
12 | ! DECLARATIONS |
---|
13 | !==================================================================== |
---|
14 | |
---|
15 | ! Arguments |
---|
16 | INTEGER klon |
---|
17 | REAL zpspsk(klon), pplev(klon) |
---|
18 | REAL ztemp(klon), zqta(klon), zqsat(klon) |
---|
19 | LOGICAL active(klon) |
---|
20 | |
---|
21 | ! Variables locales |
---|
22 | INTEGER ig, iter |
---|
23 | REAL Tbef(klon), DT(klon) |
---|
24 | REAL tdelta, qsatbef, zcor, qlbef, zdelta, zcvm5, dqsat, num, denom, dqsat_dT |
---|
25 | LOGICAL Zsat |
---|
26 | REAL RLvCp |
---|
27 | |
---|
28 | REAL, SAVE :: DDT0 = .01 |
---|
29 | !$OMP THREADPRIVATE(DDT0) |
---|
30 | |
---|
31 | LOGICAL afaire(klon), tout_converge |
---|
32 | |
---|
33 | !==================================================================== |
---|
34 | ! INITIALISATIONS |
---|
35 | !==================================================================== |
---|
36 | |
---|
37 | RLvCp = RLVTT / RCPD |
---|
38 | tout_converge = .FALSE. |
---|
39 | afaire(:) = .FALSE. |
---|
40 | DT(:) = 0. |
---|
41 | |
---|
42 | |
---|
43 | !==================================================================== |
---|
44 | ! Routine a vectoriser en copiant active dans converge et en mettant |
---|
45 | ! la boucle sur les iterations a l'exterieur est en mettant |
---|
46 | ! converge= false des que la convergence est atteinte. |
---|
47 | !==================================================================== |
---|
48 | |
---|
49 | DO ig = 1, klon |
---|
50 | IF (active(ig)) THEN |
---|
51 | Tbef(ig) = ztemp(ig) |
---|
52 | zdelta = MAX(0., SIGN(1., RTT - Tbef(ig))) |
---|
53 | qsatbef = R2ES * FOEEW(Tbef(ig), zdelta) / pplev(ig) |
---|
54 | qsatbef = MIN(0.5, qsatbef) |
---|
55 | zcor = 1. / (1. - retv * qsatbef) |
---|
56 | qsatbef = qsatbef * zcor |
---|
57 | qlbef = max(0., zqta(ig) - qsatbef) |
---|
58 | DT(ig) = 0.5 * RLvCp * qlbef |
---|
59 | zqsat(ig) = qsatbef |
---|
60 | endif |
---|
61 | END DO |
---|
62 | |
---|
63 | ! Traitement du cas ou il y a condensation mais faible |
---|
64 | ! On ne condense pas mais on dit que le qsat est le qta |
---|
65 | DO ig = 1, klon |
---|
66 | IF (active(ig)) THEN |
---|
67 | IF (0.<abs(DT(ig)).AND.abs(DT(ig))<=DDT0) THEN |
---|
68 | zqsat(ig) = zqta(ig) |
---|
69 | endif |
---|
70 | endif |
---|
71 | END DO |
---|
72 | |
---|
73 | DO iter = 1, 10 |
---|
74 | afaire(:) = abs(DT(:))>DDT0 |
---|
75 | do ig = 1, klon |
---|
76 | IF (afaire(ig)) THEN |
---|
77 | Tbef(ig) = Tbef(ig) + DT(ig) |
---|
78 | zdelta = MAX(0., SIGN(1., RTT - Tbef(ig))) |
---|
79 | qsatbef = R2ES * FOEEW(Tbef(ig), zdelta) / pplev(ig) |
---|
80 | qsatbef = MIN(0.5, qsatbef) |
---|
81 | zcor = 1. / (1. - retv * qsatbef) |
---|
82 | qsatbef = qsatbef * zcor |
---|
83 | qlbef = zqta(ig) - qsatbef |
---|
84 | zdelta = MAX(0., SIGN(1., RTT - Tbef(ig))) |
---|
85 | zcvm5 = R5LES * (1. - zdelta) + R5IES * zdelta |
---|
86 | zcor = 1. / (1. - retv * qsatbef) |
---|
87 | dqsat_dT = FOEDE(Tbef(ig), zdelta, zcvm5, qsatbef, zcor) |
---|
88 | num = -Tbef(ig) + ztemp(ig) + RLvCp * qlbef |
---|
89 | denom = 1. + RLvCp * dqsat_dT |
---|
90 | zqsat(ig) = qsatbef |
---|
91 | DT(ig) = num / denom |
---|
92 | endif |
---|
93 | enddo |
---|
94 | END DO |
---|
95 | |
---|
96 | RETURN |
---|
97 | END |
---|
98 | END MODULE lmdz_thermcell_qsat |
---|