source: LMDZ4/trunk/libf/phylmd/thermcell_dq.F90 @ 922

Last change on this file since 922 was 878, checked in by Laurent Fairhead, 17 years ago

Bascule de la physique du LMD vers la physique avec thermiques
LF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 2.7 KB
Line 
1      subroutine thermcell_dq(ngrid,nlay,ptimestep,fm,entr,  &
2     &           masse,q,dq,qa,lev_out)
3      implicit none
4
5!=======================================================================
6!
7!   Calcul du transport verticale dans la couche limite en presence
8!   de "thermiques" explicitement representes
9!   calcul du dq/dt une fois qu'on connait les ascendances
10!
11!=======================================================================
12
13      integer ngrid,nlay
14
15      real ptimestep
16      real masse(ngrid,nlay),fm(ngrid,nlay+1)
17      real entr(ngrid,nlay)
18      real q(ngrid,nlay)
19      real dq(ngrid,nlay)
20      integer lev_out                           ! niveau pour les print
21
22      real qa(ngrid,nlay),detr(ngrid,nlay),wqd(ngrid,nlay+1)
23
24      integer ig,k
25
26!   calcul du detrainement
27
28      if (lev_out.ge.1) print*,'Q2 THERMCEL_DQ 0'
29
30      do k=1,nlay
31         do ig=1,ngrid
32            detr(ig,k)=fm(ig,k)-fm(ig,k+1)+entr(ig,k)
33!           print*,'Q2 DQ ',detr(ig,k),fm(ig,k),entr(ig,k)
34!test
35            if (detr(ig,k).lt.0.) then
36               entr(ig,k)=entr(ig,k)-detr(ig,k)
37               detr(ig,k)=0.
38!               print*,'detr2<0!!!','ig=',ig,'k=',k,'f=',fm(ig,k),
39!     s         'f+1=',fm(ig,k+1),'e=',entr(ig,k),'d=',detr(ig,k)
40            endif
41            if (fm(ig,k+1).lt.0.) then
42!               print*,'fm2<0!!!'
43            endif
44            if (entr(ig,k).lt.0.) then
45!               print*,'entr2<0!!!'
46            endif
47         enddo
48      enddo
49
50!   calcul de la valeur dans les ascendances
51      do ig=1,ngrid
52         qa(ig,1)=q(ig,1)
53      enddo
54
55      do k=2,nlay
56         do ig=1,ngrid
57            if ((fm(ig,k+1)+detr(ig,k))*ptimestep.gt.  &
58     &         1.e-5*masse(ig,k)) then
59         qa(ig,k)=(fm(ig,k)*qa(ig,k-1)+entr(ig,k)*q(ig,k))  &
60     &         /(fm(ig,k+1)+detr(ig,k))
61            else
62               qa(ig,k)=q(ig,k)
63            endif
64            if (qa(ig,k).lt.0.) then
65!               print*,'qa<0!!!'
66            endif
67            if (q(ig,k).lt.0.) then
68!               print*,'q<0!!!'
69            endif
70         enddo
71      enddo
72
73      do k=2,nlay
74         do ig=1,ngrid
75!             wqd(ig,k)=fm(ig,k)*0.5*(q(ig,k-1)+q(ig,k))
76            wqd(ig,k)=fm(ig,k)*q(ig,k)
77            if (wqd(ig,k).lt.0.) then
78!               print*,'wqd<0!!!'
79            endif
80         enddo
81      enddo
82      do ig=1,ngrid
83         wqd(ig,1)=0.
84         wqd(ig,nlay+1)=0.
85      enddo
86     
87      do k=1,nlay
88         do ig=1,ngrid
89            dq(ig,k)=(detr(ig,k)*qa(ig,k)-entr(ig,k)*q(ig,k)  &
90     &               -wqd(ig,k)+wqd(ig,k+1))  &
91     &               /masse(ig,k)
92!            if (dq(ig,k).lt.0.) then
93!               print*,'dq<0!!!'
94!            endif
95         enddo
96      enddo
97
98      return
99      end
Note: See TracBrowser for help on using the repository browser.