source: LMDZ6/trunk/libf/phylmd/lmdz_thermcell_qsat.F90 @ 4608

Last change on this file since 4608 was 4593, checked in by yann meurdesoif, 18 months ago

Replace #include (c preprocessor) by INCLUDE (fortran keyword)

in phylmd (except rrtm and ecrad) filtrez, dy3dmem and dyn3dcommon

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