1 |
|
|
module exner_hyb_m |
2 |
|
|
|
3 |
|
|
IMPLICIT NONE |
4 |
|
|
|
5 |
|
|
contains |
6 |
|
|
|
7 |
|
2306 |
SUBROUTINE exner_hyb ( ngrid, ps, p, pks, pk, pkf ) |
8 |
|
|
|
9 |
|
|
! Auteurs : P.Le Van , Fr. Hourdin . |
10 |
|
|
! .......... |
11 |
|
|
! |
12 |
|
|
! .... ngrid, ps,p sont des argum.d'entree au sous-prog ... |
13 |
|
|
! .... pks,pk,pkf sont des argum.de sortie au sous-prog ... |
14 |
|
|
! |
15 |
|
|
! ************************************************************************ |
16 |
|
|
! Calcule la fonction d'Exner pk = Cp * (p/preff) ** kappa , aux milieux des |
17 |
|
|
! couches . Pk(l) sera calcule aux milieux des couches l ,entre les |
18 |
|
|
! pressions p(l) et p(l+1) ,definis aux interfaces des llm couches . |
19 |
|
|
! ************************************************************************ |
20 |
|
|
! .. N.B : Au sommet de l'atmosphere, p(llm+1) = 0. , et ps et pks sont |
21 |
|
|
! la pression et la fonction d'Exner au sol . |
22 |
|
|
! |
23 |
|
|
! -------- z |
24 |
|
|
! A partir des relations ( 1 ) p*dz(pk) = kappa *pk*dz(p) et |
25 |
|
|
! ( 2 ) pk(l) = alpha(l)+ beta(l)*pk(l-1) |
26 |
|
|
! ( voir note de Fr.Hourdin ) , |
27 |
|
|
! |
28 |
|
|
! on determine successivement , du haut vers le bas des couches, les |
29 |
|
|
! coef. alpha(llm),beta(llm) .,.,alpha(l),beta(l),,,alpha(2),beta(2), |
30 |
|
|
! puis pk(ij,1). Ensuite ,on calcule,du bas vers le haut des couches, |
31 |
|
|
! pk(ij,l) donne par la relation (2), pour l = 2 a l = llm . |
32 |
|
|
! |
33 |
|
|
! |
34 |
|
|
! |
35 |
|
|
USE comconst_mod, ONLY: jmp1, cpp, kappa, r |
36 |
|
|
USE comvert_mod, ONLY: preff |
37 |
|
|
|
38 |
|
|
IMPLICIT NONE |
39 |
|
|
|
40 |
|
|
include "dimensions.h" |
41 |
|
|
include "paramet.h" |
42 |
|
|
include "comgeom.h" |
43 |
|
|
|
44 |
|
|
INTEGER ngrid |
45 |
|
|
REAL p(ngrid,llmp1),pk(ngrid,llm) |
46 |
|
|
real, optional:: pkf(ngrid,llm) |
47 |
|
4612 |
REAL ps(ngrid),pks(ngrid), alpha(ngrid,llm),beta(ngrid,llm) |
48 |
|
|
|
49 |
|
|
! .... variables locales ... |
50 |
|
|
|
51 |
|
|
INTEGER l, ij |
52 |
|
|
REAL unpl2k,dellta |
53 |
|
|
|
54 |
|
|
logical,save :: firstcall=.true. |
55 |
|
|
character(len=*),parameter :: modname="exner_hyb" |
56 |
|
|
|
57 |
|
|
! Sanity check |
58 |
✓✓ |
2306 |
if (firstcall) then |
59 |
|
|
! sanity checks for Shallow Water case (1 vertical layer) |
60 |
|
|
if (llm.eq.1) then |
61 |
|
|
if (kappa.ne.1) then |
62 |
|
|
call abort_gcm(modname, & |
63 |
|
|
"kappa!=1 , but running in Shallow Water mode!!",42) |
64 |
|
|
endif |
65 |
|
|
if (cpp.ne.r) then |
66 |
|
|
call abort_gcm(modname, & |
67 |
|
|
"cpp!=r , but running in Shallow Water mode!!",42) |
68 |
|
|
endif |
69 |
|
|
endif ! of if (llm.eq.1) |
70 |
|
|
|
71 |
|
1 |
firstcall=.false. |
72 |
|
|
endif ! of if (firstcall) |
73 |
|
|
|
74 |
|
|
! Specific behaviour for Shallow Water (1 vertical layer) case: |
75 |
|
|
if (llm.eq.1) then |
76 |
|
|
|
77 |
|
|
! Compute pks(:),pk(:),pkf(:) |
78 |
|
|
|
79 |
|
|
DO ij = 1, ngrid |
80 |
|
|
pks(ij) = (cpp/preff) * ps(ij) |
81 |
|
|
pk(ij,1) = .5*pks(ij) |
82 |
|
|
ENDDO |
83 |
|
|
|
84 |
|
|
if (present(pkf)) then |
85 |
|
|
pkf = pk |
86 |
|
|
CALL filtreg ( pkf, jmp1, llm, 2, 1, .TRUE., 1 ) |
87 |
|
|
end if |
88 |
|
|
|
89 |
|
|
! our work is done, exit routine |
90 |
|
|
return |
91 |
|
|
endif ! of if (llm.eq.1) |
92 |
|
|
|
93 |
|
|
! General case: |
94 |
|
|
|
95 |
|
2306 |
unpl2k = 1.+ 2.* kappa |
96 |
|
|
|
97 |
|
|
! ------------- |
98 |
|
|
! Calcul de pks |
99 |
|
|
! ------------- |
100 |
|
|
|
101 |
✓✓ |
2513540 |
DO ij = 1, ngrid |
102 |
|
2513540 |
pks(ij) = cpp * ( ps(ij)/preff ) ** kappa |
103 |
|
|
ENDDO |
104 |
|
|
|
105 |
|
|
! .... Calcul des coeff. alpha et beta pour la couche l = llm .. |
106 |
|
|
! |
107 |
✓✓ |
2513540 |
DO ij = 1, ngrid |
108 |
|
2511234 |
alpha(ij,llm) = 0. |
109 |
|
2513540 |
beta (ij,llm) = 1./ unpl2k |
110 |
|
|
ENDDO |
111 |
|
|
! |
112 |
|
|
! ... Calcul des coeff. alpha et beta pour l = llm-1 a l = 2 ... |
113 |
|
|
! |
114 |
✓✓ |
87628 |
DO l = llm -1 , 2 , -1 |
115 |
|
|
! |
116 |
✓✓ |
93003286 |
DO ij = 1, ngrid |
117 |
|
92915658 |
dellta = p(ij,l)* unpl2k + p(ij,l+1)* ( beta(ij,l+1)-unpl2k ) |
118 |
|
92915658 |
alpha(ij,l) = - p(ij,l+1) / dellta * alpha(ij,l+1) |
119 |
|
93000980 |
beta (ij,l) = p(ij,l ) / dellta |
120 |
|
|
ENDDO |
121 |
|
|
ENDDO |
122 |
|
|
|
123 |
|
|
! *********************************************************************** |
124 |
|
|
! ..... Calcul de pk pour la couche 1 , pres du sol .... |
125 |
|
|
! |
126 |
✓✓ |
2513540 |
DO ij = 1, ngrid |
127 |
|
|
pk(ij,1) = ( p(ij,1)*pks(ij) - 0.5*alpha(ij,2)*p(ij,2) ) / & |
128 |
|
2513540 |
( p(ij,1)* (1.+kappa) + 0.5*( beta(ij,2)-unpl2k )* p(ij,2) ) |
129 |
|
|
ENDDO |
130 |
|
|
! |
131 |
|
|
! ..... Calcul de pk(ij,l) , pour l = 2 a l = llm ........ |
132 |
|
|
! |
133 |
✓✓ |
89934 |
DO l = 2, llm |
134 |
✓✓ |
95516826 |
DO ij = 1, ngrid |
135 |
|
95514520 |
pk(ij,l) = alpha(ij,l) + beta(ij,l) * pk(ij,l-1) |
136 |
|
|
ENDDO |
137 |
|
|
ENDDO |
138 |
|
|
|
139 |
✓✗ |
2306 |
if (present(pkf)) then |
140 |
|
|
! calcul de pkf |
141 |
✓✓✓✓
|
98030366 |
pkf = pk |
142 |
|
2306 |
CALL filtreg ( pkf, jmp1, llm, 2, 1, .TRUE., 1 ) |
143 |
|
|
end if |
144 |
|
|
|
145 |
✓✗ |
2306 |
END SUBROUTINE exner_hyb |
146 |
|
|
|
147 |
|
|
end module exner_hyb_m |