source: LMDZ6/branches/Ocean_skin/libf/phylmd/calcratqs.F90 @ 4741

Last change on this file since 4741 was 4013, checked in by lguez, 3 years ago

Sync latest trunk changes to Ocean_skin
  • 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: 9.1 KB
RevLine 
[1689]1SUBROUTINE calcratqs(klon,klev,prt_level,lunout,       &
[2236]2           iflag_ratqs,iflag_con,iflag_cld_th,pdtphys, &
[2534]3           ratqsbas,ratqshaut,ratqsp0,ratqsdp, &
[4013]4           tau_ratqs,fact_cldcon,wake_s, wake_deltaq,   &
5           ptconv,ptconvth,clwcon0th, rnebcon0th,       &
6           paprs,pplay,t_seri,q_seri,                   &
7           qtc_cv, sigt_cv, zqsat,             &
8           tke,tke_dissip,lmix,wprime, &
9           t2m,q2m,fm_therm, &
10           ratqs,ratqsc,ratqs_inter)
[1689]11
[4013]12
13USE indice_sol_mod
14USE phys_state_var_mod, ONLY: pctsrf
15USE calcratqs_multi_mod, ONLY: calcratqs_inter, calcratqs_oro, calcratqs_hetero, calcratqs_tke
16
[1689]17implicit none
18
19!========================================================================
20! Computation of ratqs, the width of the subrid scale water distribution
21! (normalized by the mean value)
22! Various options controled by flags iflag_con and iflag_ratqs
23! F Hourdin 2012/12/06
24!========================================================================
25
26! Declarations
27
28! Input
29integer,intent(in) :: klon,klev,prt_level,lunout
[2236]30integer,intent(in) :: iflag_con,iflag_cld_th,iflag_ratqs
[1689]31real,intent(in) :: pdtphys,ratqsbas,ratqshaut,fact_cldcon,tau_ratqs
[2534]32real,intent(in) :: ratqsp0, ratqsdp
[4013]33real, dimension(klon,klev+1),intent(in) :: paprs,tke,tke_dissip,lmix,wprime
34real, dimension(klon,klev),intent(in) :: pplay,t_seri,q_seri,zqsat,fm_therm, qtc_cv, sigt_cv
[1689]35logical, dimension(klon,klev),intent(in) :: ptconv
36real, dimension(klon,klev),intent(in) :: rnebcon0th,clwcon0th
[4013]37real, dimension(klon,klev),intent(in) :: wake_deltaq,wake_s
38real, dimension(klon,nbsrf),intent(in) :: t2m,q2m
39! Output
40real, dimension(klon,klev),intent(inout) :: ratqs,ratqsc,ratqs_inter
[1689]41
42logical, dimension(klon,klev),intent(inout) :: ptconvth
43
44! local
45integer i,k
46real, dimension(klon,klev) :: ratqss
47real facteur,zfratqs1,zfratqs2
[4013]48real, dimension(klon,klev) :: ratqs_hetero,ratqs_oro,ratqs_tke
[1689]49
[4013]50
[1689]51!-------------------------------------------------------------------------
52!  Caclul des ratqs
53!-------------------------------------------------------------------------
54
55!      print*,'calcul des ratqs'
56!   ratqs convectifs a l'ancienne en fonction de q(z=0)-q / q
57!   ----------------
58!   on ecrase le tableau ratqsc calcule par clouds_gno
[2236]59      if (iflag_cld_th.eq.1) then
[1689]60         do k=1,klev
61         do i=1,klon
62            if(ptconv(i,k)) then
63              ratqsc(i,k)=ratqsbas &
64              +fact_cldcon*(q_seri(i,1)-q_seri(i,k))/q_seri(i,k)
65            else
66               ratqsc(i,k)=0.
67            endif
68         enddo
69         enddo
70
71!-----------------------------------------------------------------------
72!  par nversion de la fonction log normale
73!-----------------------------------------------------------------------
[2236]74      else if (iflag_cld_th.eq.4) then
[1689]75         ptconvth(:,:)=.false.
76         ratqsc(:,:)=0.
77         if(prt_level.ge.9) print*,'avant clouds_gno thermique'
78         call clouds_gno &
79         (klon,klev,q_seri,zqsat,clwcon0th,ptconvth,ratqsc,rnebcon0th)
80         if(prt_level.ge.9) print*,' CLOUDS_GNO OK'
81       
82       endif
83
84!   ratqs stables
85!   -------------
86
87      if (iflag_ratqs.eq.0) then
88
89! Le cas iflag_ratqs=0 correspond a la version IPCC 2005 du modele.
90         do k=1,klev
91            do i=1, klon
92               ratqss(i,k)=ratqsbas+(ratqshaut-ratqsbas)* &
93               min((paprs(i,1)-pplay(i,k))/(paprs(i,1)-30000.),1.)
94            enddo
95         enddo
96
97! Pour iflag_ratqs=1 ou 2, le ratqs est constant au dessus de
98! 300 hPa (ratqshaut), varie lineariement en fonction de la pression
99! entre 600 et 300 hPa et est soit constant (ratqsbas) pour iflag_ratqs=1
100! soit lineaire (entre 0 a la surface et ratqsbas) pour iflag_ratqs=2
101! Il s'agit de differents tests dans la phase de reglage du modele
102! avec thermiques.
103
104      else if (iflag_ratqs.eq.1) then
105
106         do k=1,klev
107            do i=1, klon
108               if (pplay(i,k).ge.60000.) then
109                  ratqss(i,k)=ratqsbas
110               else if ((pplay(i,k).ge.30000.).and.(pplay(i,k).lt.60000.)) then
111                  ratqss(i,k)=ratqsbas+(ratqshaut-ratqsbas)*(60000.-pplay(i,k))/(60000.-30000.)
112               else
113                  ratqss(i,k)=ratqshaut
114               endif
115            enddo
116         enddo
117
118      else if (iflag_ratqs.eq.2) then
119
120         do k=1,klev
121            do i=1, klon
122               if (pplay(i,k).ge.60000.) then
123                  ratqss(i,k)=ratqsbas*(paprs(i,1)-pplay(i,k))/(paprs(i,1)-60000.)
124               else if ((pplay(i,k).ge.30000.).and.(pplay(i,k).lt.60000.)) then
125                    ratqss(i,k)=ratqsbas+(ratqshaut-ratqsbas)*(60000.-pplay(i,k))/(60000.-30000.)
126               else
127                    ratqss(i,k)=ratqshaut
128               endif
129            enddo
130         enddo
131
132      else if (iflag_ratqs==3) then
133         do k=1,klev
134           ratqss(:,k)=ratqsbas+(ratqshaut-ratqsbas) &
135           *min( ((paprs(:,1)-pplay(:,k))/70000.)**2 , 1. )
136         enddo
137
[4013]138      else if (iflag_ratqs==4) then
[1689]139         do k=1,klev
140           ratqss(:,k)=ratqsbas+0.5*(ratqshaut-ratqsbas) &
[2534]141!          *( tanh( (50000.-pplay(:,k))/20000.) + 1.)
142           *( tanh( (ratqsp0-pplay(:,k))/ratqsdp) + 1.)
[1689]143         enddo
144
[4013]145       else if (iflag_ratqs .GT. 9) then
146 
147       ! interactive ratqs calculations that depend on cold pools, orography, surface heterogeneity and small-scale turbulence
148       ! This should help getting a more realistic ratqs in the low and mid troposphere
149       ! We however need a "background" ratqs to account for subgrid distribution of qt (or qt/qs)
150       ! in the high troposphere
151       
152       ! background ratqs and initialisations
153          do k=1,klev
154             do i=1,klon
155              ratqss(i,k)=ratqsbas+0.5*(ratqshaut-ratqsbas) &
156              *( tanh( (ratqsp0-pplay(i,k))/ratqsdp) + 1.)
157              ratqss(i,k)=max(ratqss(i,k),0.0)
158
159              ratqs_hetero(i,k)=0.
160              ratqs_oro(i,k)=0.
161              ratqs_tke(i,k)=0.
162              ratqs_inter(i,k)=0
163             enddo
164          enddo
165     
166          if (iflag_ratqs .EQ. 10) then
167             ! interactive ratqs in presence of cold pools     
168             call calcratqs_inter(klon,klev,iflag_ratqs,pdtphys,ratqsbas,wake_deltaq,wake_s,q_seri,qtc_cv, sigt_cv,ratqs_inter)
169             do k=1,klev
170                do i=1,klon
171                    ratqs_inter(i,k)=ratqs_inter(i,k)-0.5*ratqs_inter(i,k)*(tanh((ratqsp0-pplay(i,k))/ratqsdp)+1.)
172                enddo
173             enddo
174             ratqss=ratqss+ratqs_inter
175          else if (iflag_ratqs .EQ. 11) then
176            ! interactive ratqs with several sources
177            call calcratqs_inter(klon,klev,iflag_ratqs,pdtphys,ratqsbas,wake_deltaq,wake_s,q_seri,qtc_cv, sigt_cv,ratqs_inter)
178             ratqss=ratqss+ratqs_inter
179          else if (iflag_ratqs .EQ. 12) then
180             ! contribution of surface heterogeneities to ratqs
181             call calcratqs_hetero(klon,klev,t2m,q2m,t_seri,q_seri,pplay,paprs,ratqs_hetero)
182             ratqss=ratqss+ratqs_hetero
183          else if (iflag_ratqs .EQ. 13) then
184             ! contribution of ubgrid orography to ratqs
185             call calcratqs_oro(klon,klev,zqsat,t_seri,pplay,paprs,ratqs_oro)
186             ratqss=ratqss+ratqs_oro
187          else if (iflag_ratqs .EQ. 14) then
188             ! effect of subgrid-scale TKE on ratqs (in development)
189             call calcratqs_tke(klon,klev,pdtphys,t_seri,q_seri,zqsat,pplay,paprs,tke,tke_dissip,lmix,wprime,ratqs_tke)     
190             ratqss=ratqss+ratqs_tke
191          endif
192         
193     
[1689]194      endif
195
196
197!  ratqs final
198!  -----------
199
[2236]200      if (iflag_cld_th.eq.1 .or.iflag_cld_th.eq.2.or.iflag_cld_th.eq.4) then
[1689]201
202! On ajoute une constante au ratqsc*2 pour tenir compte de
203! fluctuations turbulentes de petite echelle
204
205         do k=1,klev
206            do i=1,klon
207               if ((fm_therm(i,k).gt.1.e-10)) then
208                  ratqsc(i,k)=sqrt(ratqsc(i,k)**2+0.05**2)
209               endif
210            enddo
211         enddo
212
213!   les ratqs sont une combinaison de ratqss et ratqsc
214       if(prt_level.ge.9) write(lunout,*)'PHYLMD NOUVEAU TAU_RATQS ',tau_ratqs
215
216         if (tau_ratqs>1.e-10) then
217            facteur=exp(-pdtphys/tau_ratqs)
218         else
219            facteur=0.
220         endif
221         ratqs(:,:)=ratqsc(:,:)*(1.-facteur)+ratqs(:,:)*facteur
222!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
223! FH 22/09/2009
224! La ligne ci-dessous faisait osciller le modele et donnait une solution
225! assymptotique bidon et dépendant fortement du pas de temps.
226!        ratqs(:,:)=sqrt(ratqs(:,:)**2+ratqss(:,:)**2)
227!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
228         ratqs(:,:)=max(ratqs(:,:),ratqss(:,:))
[2236]229      else if (iflag_cld_th<=6) then
[1689]230!   on ne prend que le ratqs stable pour fisrtilp
231         ratqs(:,:)=ratqss(:,:)
232      else
233          zfratqs1=exp(-pdtphys/10800.)
234          zfratqs2=exp(-pdtphys/10800.)
235          do k=1,klev
236             do i=1,klon
237                if (ratqsc(i,k).gt.1.e-10) then
[2236]238                   ratqs(i,k)=ratqs(i,k)*zfratqs2+(iflag_cld_th/100.)*ratqsc(i,k)*(1.-zfratqs2)
[1689]239                endif
240                ratqs(i,k)=min(ratqs(i,k)*zfratqs1+ratqss(i,k)*(1.-zfratqs1),0.5)
241             enddo
242          enddo
243      endif
244
245
246return
247end
Note: See TracBrowser for help on using the repository browser.