source: LMDZ6/branches/Amaury_dev/libf/phylmd/lmdz_thermcell_qsat.F90 @ 5143

Last change on this file since 5143 was 5143, checked in by abarral, 3 months ago

Put YOEGWD.h, FCTTRE.h into modules

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
File size: 2.8 KB
Line 
1MODULE lmdz_thermcell_qsat
2CONTAINS
3
4SUBROUTINE thermcell_qsat(klon,active,pplev,ztemp,zqta,zqsat)
5  USE lmdz_YOETHF
6  USE lmdz_fcttre, ONLY: foeew, foede, qsats, qsatl, dqsats, dqsatl, thermcep
7IMPLICIT NONE
8
9  INCLUDE "YOMCST.h"
10
11!====================================================================
12! DECLARATIONS
13!====================================================================
14
15! Arguments
16INTEGER klon
17REAL zpspsk(klon),pplev(klon)
18REAL ztemp(klon),zqta(klon),zqsat(klon)
19LOGICAL active(klon)
20
21! Variables locales
22INTEGER ig,iter
23REAL Tbef(klon),DT(klon)
24REAL tdelta,qsatbef,zcor,qlbef,zdelta,zcvm5,dqsat,num,denom,dqsat_dT
25LOGICAL Zsat
26REAL RLvCp
27
28REAL, SAVE :: DDT0=.01
29!$OMP THREADPRIVATE(DDT0)
30
31LOGICAL afaire(klon),tout_converge
32
33!====================================================================
34! INITIALISATIONS
35!====================================================================
36
37RLvCp = RLVTT/RCPD
38tout_converge=.FALSE.
39afaire(:)=.FALSE.
40DT(:)=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
49DO 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
61END 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
65DO 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
71END DO
72
73DO 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
94END DO
95
96RETURN
97END
98END MODULE lmdz_thermcell_qsat
Note: See TracBrowser for help on using the repository browser.