21 |
real, intent(in):: g |
real, intent(in):: g |
22 |
|
|
23 |
REAL zlev(ngrid, klev+1) |
REAL zlev(ngrid, klev+1) |
24 |
! altitude à chaque niveau (interface inférieure de la couche de |
! altitude \`a chaque niveau (interface inf\'erieure de la couche de |
25 |
! même indice) |
! m\^eme indice) |
26 |
|
|
27 |
REAL zlay(ngrid, klev) ! altitude au centre de chaque couche |
REAL zlay(ngrid, klev) ! altitude au centre de chaque couche |
28 |
|
|
29 |
REAL u(ngrid, klev), v(ngrid, klev) |
REAL u(ngrid, klev), v(ngrid, klev) |
30 |
! vitesse au centre de chaque couche (en entrée : la valeur au |
! vitesse au centre de chaque couche (en entr\'ee : la valeur au |
31 |
! début du pas de temps) |
! d\'ebut du pas de temps) |
32 |
|
|
33 |
REAL, intent(in): teta(ngrid, klev) |
REAL, intent(in):: teta(ngrid, klev) |
34 |
! température potentielle au centre de chaque couche (en entrée : |
! temp\'erature potentielle au centre de chaque couche (en entr\'ee : |
35 |
! la valeur au début du pas de temps) |
! la valeur au d\'ebut du pas de temps) |
36 |
|
|
37 |
REAL, intent(in):: cd(:) ! (ngrid) cdrag, valeur au début du pas de temps |
REAL, intent(in):: cd(:) ! (ngrid) cdrag, valeur au d\'ebut du pas de temps |
38 |
|
|
39 |
REAL, intent(inout):: q2(ngrid, klev+1) |
REAL, intent(inout):: q2(ngrid, klev+1) |
40 |
! $q^2$ au bas de chaque couche |
! $q^2$ au bas de chaque couche |
41 |
! En entrée : la valeur au début du pas de temps ; en sortie : la |
! En entr\'ee : la valeur au d\'ebut du pas de temps ; en sortie : la |
42 |
! valeur à la fin du pas de temps. |
! valeur \`a la fin du pas de temps. |
43 |
|
|
44 |
REAL km(ngrid, klev+1) |
REAL km(ngrid, klev+1) |
45 |
! diffusivité turbulente de quantité de mouvement (au bas de |
! diffusivit\'e turbulente de quantit\'e de mouvement (au bas de |
46 |
! chaque couche) (en sortie : la valeur à la fin du pas de temps) |
! chaque couche) (en sortie : la valeur \`a la fin du pas de temps) |
47 |
|
|
48 |
REAL kn(ngrid, klev+1) |
REAL kn(ngrid, klev+1) |
49 |
! diffusivité turbulente des scalaires (au bas de chaque couche) |
! diffusivit\'e turbulente des scalaires (au bas de chaque couche) |
50 |
! (en sortie : la valeur à la fin du pas de temps) |
! (en sortie : la valeur \`a la fin du pas de temps) |
51 |
|
|
52 |
REAL kq(ngrid, klev+1) |
REAL kq(ngrid, klev+1) |
53 |
real ustar(ngrid) |
real ustar(ngrid) |
54 |
|
|
55 |
integer iflag_pbl |
integer, intent(in):: iflag_pbl |
56 |
! iflag_pbl doit valoir entre 6 et 9 |
! iflag_pbl doit valoir entre 6 et 9 |
57 |
! l = 6, on prend systématiquement une longueur d'équilibre |
! l = 6, on prend syst\'ematiquement une longueur d'\'equilibre |
58 |
! iflag_pbl = 6 : MY 2.0 |
! iflag_pbl = 6 : MY 2.0 |
59 |
! iflag_pbl = 7 : MY 2.0.Fournier |
! iflag_pbl = 7 : MY 2.0.Fournier |
60 |
! iflag_pbl = 8 : MY 2.5 |
! iflag_pbl = 8 : MY 2.5 |
69 |
REAL kmpre(ngrid, klev+1), tmp2 |
REAL kmpre(ngrid, klev+1), tmp2 |
70 |
REAL mpre(ngrid, klev+1) |
REAL mpre(ngrid, klev+1) |
71 |
real delta(ngrid, klev+1) |
real delta(ngrid, klev+1) |
72 |
real aa(ngrid, klev+1), aa0, aa1 |
real aa(ngrid, klev+1), aa1 |
73 |
integer, PARAMETER:: nlev = klev+1 |
integer, PARAMETER:: nlev = klev+1 |
74 |
logical:: first = .true. |
logical:: first = .true. |
75 |
integer:: ipas = 0 |
integer:: ipas = 0 |
80 |
real dtetadz(ngrid, klev+1) |
real dtetadz(ngrid, klev+1) |
81 |
real m2cstat, mcstat, kmcstat |
real m2cstat, mcstat, kmcstat |
82 |
real l(ngrid, klev+1) |
real l(ngrid, klev+1) |
83 |
real, save:: l0(ngrid) |
real l0(ngrid) |
84 |
real sq(ngrid), sqz(ngrid), zz(ngrid, klev+1) |
real sq(ngrid), sqz(ngrid), zz(ngrid, klev+1) |
85 |
integer iter |
integer iter |
86 |
real:: ric = 0.195, rifc = 0.191, b1 = 16.6, kap = 0.4 |
real:: ric = 0.195, rifc = 0.191, b1 = 16.6, kap = 0.4 |
|
real rino(ngrid, klev+1), smyam(ngrid, klev), styam(ngrid, klev) |
|
|
real lyam(ngrid, klev), knyam(ngrid, klev) |
|
87 |
|
|
88 |
!----------------------------------------------------------------------- |
!----------------------------------------------------------------------- |
89 |
|
|
130 |
else |
else |
131 |
rif(ig, k) = rifc |
rif(ig, k) = rifc |
132 |
endif |
endif |
133 |
if(rif(ig, k).lt.0.16) then |
if (rif(ig, k).lt.0.16) then |
134 |
alpha(ig, k) = falpha(rif(ig, k)) |
alpha(ig, k) = falpha(rif(ig, k)) |
135 |
sm(ig, k) = fsm(rif(ig, k)) |
sm(ig, k) = fsm(rif(ig, k)) |
136 |
else |
else |
141 |
enddo |
enddo |
142 |
enddo |
enddo |
143 |
|
|
144 |
! Au premier appel, on détermine l et q2 de façon itérative. |
! Au premier appel, on d\'etermine l et q2 de façon it\'erative. |
145 |
! Itération pour déterminer la longueur de mélange |
! It\'eration pour d\'eterminer la longueur de m\'elange |
146 |
|
|
147 |
if (first .or. iflag_pbl == 6) then |
if (first .or. iflag_pbl == 6) then |
148 |
do ig = 1, ngrid |
do ig = 1, ngrid |
197 |
do k = 2, klev |
do k = 2, klev |
198 |
do ig = 1, ngrid |
do ig = 1, ngrid |
199 |
l(ig, k) = fl(zlev(ig, k), l0(ig), q2(ig, k), n2(ig, k)) |
l(ig, k) = fl(zlev(ig, k), l0(ig), q2(ig, k), n2(ig, k)) |
200 |
if(first) then |
if (first) then |
201 |
q2(ig, k) = l(ig, k)**2 * zz(ig, k) |
q2(ig, k) = l(ig, k)**2 * zz(ig, k) |
202 |
endif |
endif |
203 |
enddo |
enddo |
263 |
delta(ig, k) = 1.e-20 |
delta(ig, k) = 1.e-20 |
264 |
endif |
endif |
265 |
km(ig, k) = l(ig, k)*sqrt(q2(ig, k))*sm(ig, k) |
km(ig, k) = l(ig, k)*sqrt(q2(ig, k))*sm(ig, k) |
|
aa0 = (m2(ig, k)-alpha(ig, k)*n2(ig, k)-delta(ig, k)/b1) |
|
266 |
aa1 = (m2(ig, k)*(1.-rif(ig, k))-delta(ig, k)/b1) |
aa1 = (m2(ig, k)*(1.-rif(ig, k))-delta(ig, k)/b1) |
267 |
aa(ig, k) = aa1*dt/(delta(ig, k)*l(ig, k)) |
aa(ig, k) = aa1*dt/(delta(ig, k)*l(ig, k)) |
268 |
qpre = sqrt(q2(ig, k)) |
qpre = sqrt(q2(ig, k)) |
285 |
enddo |
enddo |
286 |
endif |
endif |
287 |
|
|
288 |
! Calcul des coefficients de mélange |
! Calcul des coefficients de m\'elange |
289 |
do k = 2, klev |
do k = 2, klev |
290 |
do ig = 1, ngrid |
do ig = 1, ngrid |
291 |
zq = sqrt(q2(ig, k)) |
zq = sqrt(q2(ig, k)) |
326 |
enddo |
enddo |
327 |
enddo |
enddo |
328 |
|
|
|
! Diagnostique pour stokage |
|
|
|
|
|
rino = rif |
|
|
smyam(:, 1:klev) = sm(:, 1:klev) |
|
|
styam = sm(:, 1:klev)*alpha(:, 1:klev) |
|
|
lyam(1:ngrid, 1:klev) = l(:, 1:klev) |
|
|
|
|
329 |
first = .false. |
first = .false. |
330 |
|
|
331 |
end SUBROUTINE yamada4 |
end SUBROUTINE yamada4 |