/[lmdze]/trunk/phylmd/Thermcell/dvthermcell2.f
ViewVC logotype

Annotation of /trunk/phylmd/Thermcell/dvthermcell2.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 47 - (hide annotations)
Fri Jul 1 15:00:48 2011 UTC (12 years, 11 months ago) by guez
Original Path: trunk/libf/phylmd/Thermcell/dvthermcell2.f90
File size: 4014 byte(s)
Split "thermcell.f" and "cv3_routines.f".
Removed copies of files that are now in "L_util".
Moved "mva9" and "diagetpq" to their own files.
Unified variable names across procedures.

1 guez 47 subroutine dvthermcell2(ngrid,nlay,ptimestep,fm,entr,masse &
2     ,fraca,larga &
3     ,u,v,du,dv,ua,va)
4     use dimens_m
5     use dimphy
6     implicit none
7    
8     !=======================================================================
9     !
10     ! Calcul du transport verticale dans la couche limite en presence
11     ! de "thermiques" explicitement representes
12     ! calcul du dq/dt une fois qu'on connait les ascendances
13     !
14     !=======================================================================
15    
16    
17     integer ngrid,nlay
18    
19     real ptimestep
20     real masse(ngrid,nlay),fm(ngrid,nlay+1)
21     real fraca(ngrid,nlay+1)
22     real larga(ngrid)
23     real entr(ngrid,nlay)
24     real u(ngrid,nlay)
25     real ua(ngrid,nlay)
26     real du(ngrid,nlay)
27     real v(ngrid,nlay)
28     real va(ngrid,nlay)
29     real dv(ngrid,nlay)
30    
31     real qa(klon,klev),detr(klon,klev),zf,zf2
32     real wvd(klon,klev+1),wud(klon,klev+1)
33     real gamma0,gamma(klon,klev+1)
34     real ue(klon,klev),ve(klon,klev)
35     real dua,dva
36     integer iter
37    
38     integer ig,k
39    
40     ! calcul du detrainement
41    
42     do k=1,nlay
43     do ig=1,ngrid
44     detr(ig,k)=fm(ig,k)-fm(ig,k+1)+entr(ig,k)
45     enddo
46     enddo
47    
48     ! calcul de la valeur dans les ascendances
49     do ig=1,ngrid
50     ua(ig,1)=u(ig,1)
51     va(ig,1)=v(ig,1)
52     ue(ig,1)=u(ig,1)
53     ve(ig,1)=v(ig,1)
54     enddo
55    
56     do k=2,nlay
57     do ig=1,ngrid
58     if ((fm(ig,k+1)+detr(ig,k))*ptimestep.gt. &
59     1.e-5*masse(ig,k)) then
60     ! On itère sur la valeur du coeff de freinage.
61     ! gamma0=rho(ig,k)*(zlev(ig,k+1)-zlev(ig,k))
62     gamma0=masse(ig,k) &
63     *sqrt( 0.5*(fraca(ig,k+1)+fraca(ig,k)) ) &
64     *0.5/larga(ig) &
65     *1.
66     ! s *0.5
67     ! gamma0=0.
68     zf=0.5*(fraca(ig,k)+fraca(ig,k+1))
69     zf=0.
70     zf2=1./(1.-zf)
71     ! la première fois on multiplie le coefficient de freinage
72     ! par le module du vent dans la couche en dessous.
73     dua=ua(ig,k-1)-u(ig,k-1)
74     dva=va(ig,k-1)-v(ig,k-1)
75     do iter=1,5
76     ! On choisit une relaxation lineaire.
77     gamma(ig,k)=gamma0
78     ! On choisit une relaxation quadratique.
79     gamma(ig,k)=gamma0*sqrt(dua**2+dva**2)
80     ua(ig,k)=(fm(ig,k)*ua(ig,k-1) &
81     +(zf2*entr(ig,k)+gamma(ig,k))*u(ig,k)) &
82     /(fm(ig,k+1)+detr(ig,k)+entr(ig,k)*zf*zf2 &
83     +gamma(ig,k))
84     va(ig,k)=(fm(ig,k)*va(ig,k-1) &
85     +(zf2*entr(ig,k)+gamma(ig,k))*v(ig,k)) &
86     /(fm(ig,k+1)+detr(ig,k)+entr(ig,k)*zf*zf2 &
87     +gamma(ig,k))
88     ! print*,k,ua(ig,k),va(ig,k),u(ig,k),v(ig,k),dua,dva
89     dua=ua(ig,k)-u(ig,k)
90     dva=va(ig,k)-v(ig,k)
91     ue(ig,k)=(u(ig,k)-zf*ua(ig,k))*zf2
92     ve(ig,k)=(v(ig,k)-zf*va(ig,k))*zf2
93     enddo
94     else
95     ua(ig,k)=u(ig,k)
96     va(ig,k)=v(ig,k)
97     ue(ig,k)=u(ig,k)
98     ve(ig,k)=v(ig,k)
99     gamma(ig,k)=0.
100     endif
101     enddo
102     enddo
103    
104     do k=2,nlay
105     do ig=1,ngrid
106     wud(ig,k)=fm(ig,k)*ue(ig,k)
107     wvd(ig,k)=fm(ig,k)*ve(ig,k)
108     enddo
109     enddo
110     do ig=1,ngrid
111     wud(ig,1)=0.
112     wud(ig,nlay+1)=0.
113     wvd(ig,1)=0.
114     wvd(ig,nlay+1)=0.
115     enddo
116    
117     do k=1,nlay
118     do ig=1,ngrid
119     du(ig,k)=((detr(ig,k)+gamma(ig,k))*ua(ig,k) &
120     -(entr(ig,k)+gamma(ig,k))*ue(ig,k) &
121     -wud(ig,k)+wud(ig,k+1)) &
122     /masse(ig,k)
123     dv(ig,k)=((detr(ig,k)+gamma(ig,k))*va(ig,k) &
124     -(entr(ig,k)+gamma(ig,k))*ve(ig,k) &
125     -wvd(ig,k)+wvd(ig,k+1)) &
126     /masse(ig,k)
127     enddo
128     enddo
129    
130     return
131     end

  ViewVC Help
Powered by ViewVC 1.1.21