source: LMDZ5/trunk/libf/phylmd/thermcell_alim.F90 @ 4100

Last change on this file since 4100 was 2406, checked in by fhourdin, 9 years ago

Small corrections

File size: 4.0 KB
Line 
1!
2! $Id: thermcell_plume.F90 2311 2015-06-25 07:45:24Z emillour $
3!
4      SUBROUTINE thermcell_alim(flag,ngrid,klev,ztv,d_temp,zlev,alim_star,lalim)
5IMPLICIT NONE
6
7!--------------------------------------------------------------------------
8! FH : 2015/11/06
9! thermcell_alim: calcule la distribution verticale de l'alimentation
10! laterale a la base des panaches thermiques
11!--------------------------------------------------------------------------
12
13#include "YOMCST.h"
14#include "YOETHF.h"
15#include "FCTTRE.h"
16#include "thermcell.h"
17
18!      fort(10) ptimestep,ztv,zthl,po,zl,rhobarz,zlev,pplev,pphi,zpspsk,f0
19      INTEGER, INTENT(IN) :: ngrid,klev
20      REAL, INTENT(IN) :: ztv(ngrid,klev)
21      REAL, INTENT(IN) :: d_temp(ngrid)
22      REAL, INTENT(IN) :: zlev(ngrid,klev+1)
23      REAL, INTENT(OUT) :: alim_star(ngrid,klev)
24      INTEGER, INTENT(OUT) :: lalim(ngrid)
25      INTEGER, INTENT(IN) :: flag
26
27      REAL :: alim_star_tot(ngrid),zi(ngrid),zh(ngrid)
28      REAL :: zlay(ngrid,klev)
29      REAL ztv_parcel
30
31      INTEGER ig,l
32
33      REAL h,z,falim
34      falim(h,z)=0.2*((z-h)**5+h**5)
35
36
37!===================================================================
38
39   lalim(:)=1
40   alim_star_tot(:)=0.
41
42!-------------------------------------------------------------------------
43! Definition de l'alimentation a l'origine dans thermcell_init
44!-------------------------------------------------------------------------
45   IF (flag==0) THEN ! CMIP5 version
46      do l=1,klev-1
47         do ig=1,ngrid
48            if (ztv(ig,l)> ztv(ig,l+1) .and. ztv(ig,1)>=ztv(ig,l) ) then
49               alim_star(ig,l)=MAX((ztv(ig,l)-ztv(ig,l+1)),0.)  &
50     &                       *sqrt(zlev(ig,l+1))
51               lalim(ig)=l+1
52               alim_star_tot(ig)=alim_star_tot(ig)+alim_star(ig,l)
53            endif
54         enddo
55      enddo
56      do l=1,klev
57         do ig=1,ngrid
58            if (alim_star_tot(ig) > 1.e-10 ) then
59               alim_star(ig,l)=alim_star(ig,l)/alim_star_tot(ig)
60            endif
61         enddo
62      enddo
63      alim_star_tot(:)=1.
64
65!-------------------------------------------------------------------------
66! Nouvelle definition avec possibilite d'introduire un DT en surface
67! On suppose que la forme du profile d'alimentation scale avec la hauteur
68! d'inversion calculée avec une particule partant de la premieere couche
69
70! Fonction  f(z) = z ( h - z ) , avec h = zi/3
71! On utilise l'integralle
72! Int_0^z f(z') dz' = z^2 ( h/2 - z/3 ) = falim(h,z)
73! Pour calculer l'alimentation des couches
74!-------------------------------------------------------------------------
75   ELSE
76! Computing inversion height zi and zh=zi/3.
77      zi(:)=0.
78! Il faut recalculer zlay qui n'est pas dispo dans thermcell_plume
79! A changer eventuellement.
80      do l=1,klev
81         zlay(:,l)=0.5*(zlev(:,l)+zlev(:,l+1))
82      enddo
83
84      do l=klev-1,1,-1
85         do ig=1,ngrid
86            ztv_parcel=ztv(ig,1)+d_temp(ig)
87            if (ztv_parcel<ztv(ig,l+1)) lalim(ig)=l
88         enddo
89      enddo
90
91      do ig=1,ngrid
92         l=lalim(ig)
93         IF (l==1) THEN
94            zi(ig)=0.
95         ELSE
96            ztv_parcel=ztv(ig,1)+d_temp(ig)
97            zi(ig)=zlay(ig,l)+(zlay(ig,l+1)-zlay(ig,l))/(ztv(ig,l+1)-ztv(ig,l))*(ztv_parcel-ztv(ig,l))
98         ENDIF
99      enddo
100
101      zh(:)=zi(:)/2.
102      alim_star_tot(:)=0.
103      alim_star(:,:)=0.
104      lalim(:)=0
105      do l=1,klev-1
106         do ig=1,ngrid
107            IF (zh(ig)==0.) THEN
108               alim_star(ig,l)=0.
109               lalim(ig)=1
110            ELSE IF (zlev(ig,l+1)<=zh(ig)) THEN
111               alim_star(ig,l)=(falim(zh(ig),zlev(ig,l+1))-falim(zh(ig),zlev(ig,l)))/falim(zh(ig),zh(ig))
112               lalim(ig)=l
113            ELSE IF (zlev(ig,l)<=zh(ig)) THEN
114               alim_star(ig,l)=(falim(zh(ig),zh(ig))-falim(zh(ig),zlev(ig,l)))/falim(zh(ig),zh(ig))
115               lalim(ig)=l
116            ELSE
117               alim_star(ig,l)=0.
118            ENDIF
119         ENDDO
120         alim_star_tot(:)=alim_star_tot(:)+alim_star(:,l)
121      ENDDO
122      IF (ngrid==1) print*,'NEW ALIM CALCUL DE ZI ',alim_star_tot,lalim,zi,zh
123      alim_star_tot(:)=1.
124
125   ENDIF
126
127
128RETURN
129END
Note: See TracBrowser for help on using the repository browser.