source: LMDZ4/branches/LMDZ4V5.0-dev/libf/phylmd/calltherm.F90 @ 1306

Last change on this file since 1306 was 1299, checked in by Laurent Fairhead, 15 years ago

Nettoyage general pour se rapprocher des normes et éviter des erreurs a la
compilation:

  • tous les FLOAT() sont remplacés par des REAL()
  • tous les STOP dans phylmd sont remplacés par des appels à abort_gcm
  • le common control défini dans le fichier control.h est remplacé par le module control_mod pour éviter des messages sur l'alignement des variables dans les déclarations
  • des $Header$ remplacés par des $Id$ pour svn

Quelques remplacements à faire ont pu m'échapper


General cleanup of the code to try and adhere to norms and to prevent some
compilation errors:

  • all FLOAT() instructions have been replaced by REAL() instructions
  • all STOP instructions in phylmd have been replaced by calls to abort_gcm
  • the common block control defined in the control.h file has been replaced by the control_mod to prevent compilation warnings on the alignement of declared variables
  • $Header$ replaced by $Id$ for svn

Some changes which should have been made might have escaped me

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 12.7 KB
Line 
1!
2! $Id: calltherm.F90 1299 2010-01-20 14:27:21Z fairhead $
3!
4      subroutine calltherm(dtime  &
5     &      ,pplay,paprs,pphi,weak_inversion  &
6     &      ,u_seri,v_seri,t_seri,q_seri,zqsat,debut  &
7     &      ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs  &
8     &      ,fm_therm,entr_therm,detr_therm,zqasc,clwcon0,lmax,ratqscth,  &
9     &       ratqsdiff,zqsatth,Ale_bl,Alp_bl,lalim_conv,wght_th, &
10     &       zmax0,f0,zw2,fraca)
11
12      USE dimphy
13      implicit none
14#include "dimensions.h"
15!#include "dimphy.h"
16#include "thermcell.h"
17#include "iniprint.h"
18
19!IM 140508
20      INTEGER itap
21      REAL dtime
22      LOGICAL debut
23      LOGICAL logexpr0, logexpr2(klon,klev), logexpr1(klon)
24      REAL fact(klon)
25      INTEGER nbptspb
26
27      REAL u_seri(klon,klev),v_seri(klon,klev)
28      REAL t_seri(klon,klev),q_seri(klon,klev),qmemoire(klon,klev)
29      REAL weak_inversion(klon)
30      REAL paprs(klon,klev+1)
31      REAL pplay(klon,klev)
32      REAL pphi(klon,klev)
33      real zlev(klon,klev+1)
34!test: on sort lentr et a* pour alimenter KE
35      REAL wght_th(klon,klev)
36      INTEGER lalim_conv(klon)
37      REAL zw2(klon,klev+1),fraca(klon,klev+1)
38
39!FH Update Thermiques
40      REAL d_t_ajs(klon,klev), d_q_ajs(klon,klev)
41      REAL d_u_ajs(klon,klev),d_v_ajs(klon,klev)
42      real fm_therm(klon,klev+1)
43      real entr_therm(klon,klev),detr_therm(klon,klev)
44
45!********************************************************
46!     declarations
47      LOGICAL flag_bidouille_stratocu
48      real fmc_therm(klon,klev+1),zqasc(klon,klev)
49      real zqla(klon,klev)
50      real zqta(klon,klev)
51      real wmax_sec(klon)
52      real zmax_sec(klon)
53      real f_sec(klon)
54      real detrc_therm(klon,klev)
55! FH WARNING : il semble que ces save ne servent a rien
56!     save fmc_therm, detrc_therm
57      real clwcon0(klon,klev)
58      real zqsat(klon,klev)
59      real zw_sec(klon,klev+1)
60      integer lmix_sec(klon)
61      integer lmax(klon)
62      real ratqscth(klon,klev)
63      real ratqsdiff(klon,klev)
64      real zqsatth(klon,klev) 
65!nouvelles variables pour la convection
66      real Ale_bl(klon)
67      real Alp_bl(klon)
68      real Ale(klon)
69      real Alp(klon)
70!RC
71      !on garde le zmax du pas de temps precedent
72      real zmax0(klon), f0(klon)
73!********************************************************
74
75
76! variables locales
77      REAL d_t_the(klon,klev), d_q_the(klon,klev)
78      REAL d_u_the(klon,klev),d_v_the(klon,klev)
79!
80      real zfm_therm(klon,klev+1),zdt
81      real zentr_therm(klon,klev),zdetr_therm(klon,klev)
82! FH A VERIFIER : SAVE INUTILES
83!      save zentr_therm,zfm_therm
84
85      character (len=20) :: modname='calltherm'
86      character (len=80) :: abort_message
87
88      integer i,k
89      logical, save :: first=.true.
90!$OMP THREADPRIVATE(first)
91!********************************************************
92      if (first) then
93        itap=0
94        first=.false.
95      endif
96
97! Incrementer le compteur de la physique
98     itap   = itap + 1
99
100!  Modele du thermique
101!  ===================
102!         print*,'thermiques: WARNING on passe t au lieu de t_seri'
103
104
105! On prend comme valeur initiale des thermiques la valeur du pas
106! de temps precedent
107         zfm_therm(:,:)=fm_therm(:,:)
108         zdetr_therm(:,:)=detr_therm(:,:)
109         zentr_therm(:,:)=entr_therm(:,:)
110
111! On reinitialise les flux de masse a zero pour le cumul en
112! cas de splitting
113         fm_therm(:,:)=0.
114         entr_therm(:,:)=0.
115         detr_therm(:,:)=0.
116
117         Ale_bl(:)=0.
118         Alp_bl(:)=0.
119         if (prt_level.ge.10) then
120          print*,'thermV4 nsplit: ',nsplit_thermals,' weak_inversion'
121         endif
122
123!   tests sur les valeurs negatives de l'eau
124         logexpr0=prt_level.ge.10
125         nbptspb=0
126         do k=1,klev
127            do i=1,klon
128! Attention teste abderr 19-03-09
129!               logexpr2(i,k)=.not.q_seri(i,k).ge.0.
130                logexpr2(i,k)=.not.q_seri(i,k).ge.1.e-15
131               if (logexpr2(i,k)) then
132                q_seri(i,k)=1.e-15
133                nbptspb=nbptspb+1
134               endif
135!               if (logexpr0) &
136!    &             print*,'WARN eau<0 avant therm i=',i,'  k=',k  &
137!    &         ,' dq,q',d_q_the(i,k),q_seri(i,k)
138            enddo
139         enddo
140         if(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb   
141
142         zdt=dtime/REAL(nsplit_thermals)
143         do isplit=1,nsplit_thermals
144
145          if (iflag_thermals.eq.1) then
146            CALL thermcell_2002(klon,klev,zdt   &
147     &      ,pplay,paprs,pphi  &
148     &      ,u_seri,v_seri,t_seri,q_seri  &
149     &      ,d_u_the,d_v_the,d_t_the,d_q_the  &
150     &      ,zfm_therm,zentr_therm  &
151     &      ,r_aspect_thermals,30.,w2di_thermals  &
152     &      ,tau_thermals,3)
153          else if (iflag_thermals.eq.2) then
154            CALL thermcell_sec(klon,klev,zdt  &
155     &      ,pplay,paprs,pphi,zlev  &
156     &      ,u_seri,v_seri,t_seri,q_seri  &
157     &      ,d_u_the,d_v_the,d_t_the,d_q_the  &
158     &      ,zfm_therm,zentr_therm  &
159     &      ,r_aspect_thermals,30.,w2di_thermals  &
160     &      ,tau_thermals,3)
161          else if (iflag_thermals.eq.3) then
162            CALL thermcell(klon,klev,zdt  &
163     &      ,pplay,paprs,pphi  &
164     &      ,u_seri,v_seri,t_seri,q_seri  &
165     &      ,d_u_the,d_v_the,d_t_the,d_q_the  &
166     &      ,zfm_therm,zentr_therm  &
167     &      ,r_aspect_thermals,l_mix_thermals,w2di_thermals  &
168     &      ,tau_thermals,3)
169          else if (iflag_thermals.eq.10) then
170            CALL thermcell_eau(klon,klev,zdt  &
171     &      ,pplay,paprs,pphi  &
172     &      ,u_seri,v_seri,t_seri,q_seri  &
173     &      ,d_u_the,d_v_the,d_t_the,d_q_the  &
174     &      ,zfm_therm,zentr_therm  &
175     &      ,r_aspect_thermals,l_mix_thermals,w2di_thermals  &
176     &      ,tau_thermals,3)
177          else if (iflag_thermals.eq.11) then
178              abort_message = 'cas non prevu dans calltherm'
179              CALL abort_gcm (modname,abort_message,1)
180
181!           CALL thermcell_pluie(klon,klev,zdt  &
182!   &      ,pplay,paprs,pphi,zlev  &
183!    &      ,u_seri,v_seri,t_seri,q_seri  &
184!    &      ,d_u_the,d_v_the,d_t_the,d_q_the  &
185!    &      ,zfm_therm,zentr_therm,zqla  &
186!    &      ,r_aspect_thermals,l_mix_thermals,w2di_thermals  &
187!    &      ,tau_thermals,3)
188          else if (iflag_thermals.eq.12) then
189            CALL calcul_sec(klon,klev,zdt  &
190     &      ,pplay,paprs,pphi,zlev  &
191     &      ,u_seri,v_seri,t_seri,q_seri  &
192     &      ,zmax_sec,wmax_sec,zw_sec,lmix_sec  &
193     &      ,r_aspect_thermals,l_mix_thermals,w2di_thermals  &
194     &      ,tau_thermals)
195          else if (iflag_thermals==13.or.iflag_thermals==14) then
196            CALL thermcellV0_main(itap,klon,klev,zdt  &
197     &      ,pplay,paprs,pphi,debut  &
198     &      ,u_seri,v_seri,t_seri,q_seri  &
199     &      ,d_u_the,d_v_the,d_t_the,d_q_the  &
200     &      ,zfm_therm,zentr_therm,zdetr_therm,zqasc,zqla,lmax  &
201     &      ,ratqscth,ratqsdiff,zqsatth  &
202     &      ,r_aspect_thermals,l_mix_thermals  &
203     &      ,tau_thermals,Ale,Alp,lalim_conv,wght_th &
204     &      ,zmax0,f0,zw2,fraca)
205          else if (iflag_thermals==15.or.iflag_thermals==16) then
206            CALL thermcell_main(itap,klon,klev,zdt  &
207     &      ,pplay,paprs,pphi,debut  &
208     &      ,u_seri,v_seri,t_seri,q_seri  &
209     &      ,d_u_the,d_v_the,d_t_the,d_q_the  &
210     &      ,zfm_therm,zentr_therm,zdetr_therm,zqasc,zqla,lmax  &
211     &      ,ratqscth,ratqsdiff,zqsatth  &
212     &      ,r_aspect_thermals,l_mix_thermals  &
213     &      ,tau_thermals,Ale,Alp,lalim_conv,wght_th &
214     &      ,zmax0,f0,zw2,fraca)
215           if (prt_level.gt.10) write(lunout,*)'Apres thermcell_main OK'
216         else
217           abort_message = 'Cas des thermiques non prevu'
218           CALL abort_gcm (modname,abort_message,1)
219         endif
220
221       flag_bidouille_stratocu=iflag_thermals.eq.14.or.iflag_thermals.eq.16
222
223      fact(:)=0.
224      DO i=1,klon
225       logexpr1(i)=flag_bidouille_stratocu.or.weak_inversion(i).gt.0.5
226       IF(logexpr1(i)) fact(i)=1./REAL(nsplit_thermals)
227      ENDDO
228
229     DO k=1,klev
230!  transformation de la derivee en tendance
231            d_t_the(:,k)=d_t_the(:,k)*dtime*fact(:)
232            d_u_the(:,k)=d_u_the(:,k)*dtime*fact(:)
233            d_v_the(:,k)=d_v_the(:,k)*dtime*fact(:)
234            d_q_the(:,k)=d_q_the(:,k)*dtime*fact(:)
235            fm_therm(:,k)=fm_therm(:,k)  &
236     &      +zfm_therm(:,k)*fact(:)
237            entr_therm(:,k)=entr_therm(:,k)  &
238     &       +zentr_therm(:,k)*fact(:)
239            detr_therm(:,k)=detr_therm(:,k)  &
240     &       +zdetr_therm(:,k)*fact(:)
241      ENDDO
242       fm_therm(:,klev+1)=0.
243
244
245
246!  accumulation de la tendance
247            d_t_ajs(:,:)=d_t_ajs(:,:)+d_t_the(:,:)
248            d_u_ajs(:,:)=d_u_ajs(:,:)+d_u_the(:,:)
249            d_v_ajs(:,:)=d_v_ajs(:,:)+d_v_the(:,:)
250            d_q_ajs(:,:)=d_q_ajs(:,:)+d_q_the(:,:)
251
252!  incrementation des variables meteo
253            t_seri(:,:) = t_seri(:,:) + d_t_the(:,:)
254            u_seri(:,:) = u_seri(:,:) + d_u_the(:,:)
255            v_seri(:,:) = v_seri(:,:) + d_v_the(:,:)
256            qmemoire(:,:)=q_seri(:,:)
257            q_seri(:,:) = q_seri(:,:) + d_q_the(:,:)
258           if (prt_level.gt.10) write(lunout,*)'Apres apres thermcell_main OK'
259
260       DO i=1,klon
261        if(prt_level.GE.10) print*,'calltherm i Alp_bl Alp Ale_bl Ale',i,Alp_bl(i),Alp(i),Ale_bl(i),Ale(i)
262            fm_therm(i,klev+1)=0.
263            Ale_bl(i)=Ale_bl(i)+Ale(i)/REAL(nsplit_thermals)
264!            write(22,*)'ALE CALLTHERM',Ale_bl(i),Ale(i)
265            Alp_bl(i)=Alp_bl(i)+Alp(i)/REAL(nsplit_thermals)
266!            write(23,*)'ALP CALLTHERM',Alp_bl(i),Alp(i)
267       ENDDO
268
269!IM 060508 marche pas comme cela !!!        enddo ! isplit
270
271!   tests sur les valeurs negatives de l'eau
272         nbptspb=0
273            DO k = 1, klev
274            DO i = 1, klon
275               logexpr2(i,k)=.not.q_seri(i,k).ge.0.
276               if (logexpr2(i,k)) then
277                q_seri(i,k)=1.e-15
278                nbptspb=nbptspb+1
279!                if (prt_level.ge.10) then
280!                  print*,'WARN eau<0 apres therm i=',i,'  k=',k  &
281!    &         ,' dq,q',d_q_the(i,k),q_seri(i,k),  &
282!    &         'fm=',zfm_therm(i,k),'entr=',entr_therm(i,k)
283                 endif
284            ENDDO
285            ENDDO
286        IF(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb   
287! tests sur les valeurs de la temperature
288        nbptspb=0
289            DO k = 1, klev
290            DO i = 1, klon
291               logexpr2(i,k)=t_seri(i,k).lt.50..or.t_seri(i,k).gt.370.
292               if (logexpr2(i,k)) nbptspb=nbptspb+1
293!              if ((t_seri(i,k).lt.50.) .or.  &
294!    &              (t_seri(i,k).gt.370.)) then
295!                 print*,'WARN temp apres therm i=',i,'  k=',k  &
296!    &         ,' t_seri',t_seri(i,k)
297!              CALL abort
298!              endif
299            ENDDO
300            ENDDO
301        IF(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
302         enddo ! isplit
303
304!
305!***************************************************************
306!     calcul du flux ascencant conservatif
307!            print*,'<<<<calcul flux ascendant conservatif'
308
309      fmc_therm=0.
310               do k=1,klev
311            do i=1,klon
312                  if (entr_therm(i,k).gt.0.) then
313                     fmc_therm(i,k+1)=fmc_therm(i,k)+entr_therm(i,k)
314                  else
315                     fmc_therm(i,k+1)=fmc_therm(i,k)
316                  endif
317                  detrc_therm(i,k)=(fmc_therm(i,k+1)-fm_therm(i,k+1))  &
318     &                 -(fmc_therm(i,k)-fm_therm(i,k))
319               enddo
320            enddo
321     
322     
323!****************************************************************
324!     calcul de l'humidite dans l'ascendance
325!      print*,'<<<<calcul de lhumidite dans thermique'
326!CR:on ne le calcule que pour le cas sec
327      if (iflag_thermals.le.11) then     
328      do i=1,klon
329         zqasc(i,1)=q_seri(i,1)
330         do k=2,klev
331            if (fmc_therm(i,k+1).gt.1.e-6) then
332               zqasc(i,k)=(fmc_therm(i,k)*zqasc(i,k-1)  &
333     &              +entr_therm(i,k)*q_seri(i,k))/fmc_therm(i,k+1)
334!CR:test on asseche le thermique
335!               zqasc(i,k)=zqasc(i,k)/2.
336!            else
337!               zqasc(i,k)=q_seri(i,k)
338            endif
339         enddo
340       enddo
341     
342
343!     calcul de l'eau condensee dans l'ascendance
344!             print*,'<<<<calcul de leau condensee dans thermique'
345             do i=1,klon
346                do k=1,klev
347                   clwcon0(i,k)=zqasc(i,k)-zqsat(i,k)
348                   if (clwcon0(i,k).lt.0. .or.   &
349     &             (fm_therm(i,k+1)+detrc_therm(i,k)).lt.1.e-6) then
350                      clwcon0(i,k)=0.
351                   endif
352                enddo
353             enddo
354       else
355              do i=1,klon
356                do k=1,klev
357                   clwcon0(i,k)=zqla(i,k) 
358                   if (clwcon0(i,k).lt.0. .or.   &
359     &             (fm_therm(i,k+1)+detrc_therm(i,k)).lt.1.e-6) then
360                   clwcon0(i,k)=0.
361                   endif
362                enddo
363             enddo
364       endif
365!*******************************************************************   
366
367
368      return
369
370      end
Note: See TracBrowser for help on using the repository browser.