source: LMDZ4/branches/V3_test/libf/dyn3dpar/exner_hyb_p.F

Last change on this file was 709, checked in by Laurent Fairhead, 18 years ago

Nouvelles versions de la dynamique YM
LF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 3.7 KB
Line 
1      SUBROUTINE  exner_hyb_p ( ngrid, ps, p,alpha,beta, pks, pk, pkf )
2c
3c     Auteurs :  P.Le Van  , Fr. Hourdin  .
4c    ..........
5c
6c    ....  ngrid, ps,p             sont des argum.d'entree  au sous-prog ...
7c    .... alpha,beta, pks,pk,pkf   sont des argum.de sortie au sous-prog ...
8c
9c   ************************************************************************
10c    Calcule la fonction d'Exner pk = Cp * p ** kappa , aux milieux des
11c    couches .   Pk(l) sera calcule aux milieux  des couches l ,entre les
12c    pressions p(l) et p(l+1) ,definis aux interfaces des llm couches .
13c   ************************************************************************
14c  .. N.B : Au sommet de l'atmosphere,  p(llm+1) = 0. , et ps et pks sont
15c    la pression et la fonction d'Exner  au  sol  .
16c
17c                                 -------- z                                   
18c    A partir des relations  ( 1 ) p*dz(pk) = kappa *pk*dz(p)      et
19c                            ( 2 ) pk(l) = alpha(l)+ beta(l)*pk(l-1)
20c    ( voir note de Fr.Hourdin )  ,
21c
22c    on determine successivement , du haut vers le bas des couches, les
23c    coef. alpha(llm),beta(llm) .,.,alpha(l),beta(l),,,alpha(2),beta(2),
24c    puis pk(ij,1). Ensuite ,on calcule,du bas vers le haut des couches, 
25c     pk(ij,l)  donne  par la relation (2),  pour l = 2 a l = llm .
26c
27c
28      USE parallel
29      IMPLICIT NONE
30c
31#include "dimensions.h"
32#include "paramet.h"
33#include "comconst.h"
34#include "comgeom.h"
35#include "comvert.h"
36#include "serre.h"
37
38      INTEGER  ngrid
39      REAL p(ngrid,llmp1),pk(ngrid,llm),pkf(ngrid,llm)
40      REAL ps(ngrid),pks(ngrid), alpha(ngrid,llm),beta(ngrid,llm)
41
42c    .... variables locales   ...
43
44      INTEGER l, ij
45      REAL unpl2k,dellta
46
47      REAL ppn(iim),pps(iim)
48      REAL xpn, xps
49      REAL SSUM
50      EXTERNAL SSUM
51      INTEGER ije,ijb,jje,jjb
52c
53     
54      unpl2k    = 1.+ 2.* kappa
55c
56      ijb=ij_begin
57      ije=ij_end
58     
59      DO   ij  = ijb, ije
60        pks(ij) = cpp * ( ps(ij)/preff ) ** kappa
61      ENDDO
62
63      if (pole_nord) then
64        DO  ij   = 1, iim
65          ppn(ij) = aire(   ij   ) * pks(  ij     )
66        ENDDO
67        xpn      = SSUM(iim,ppn,1) /apoln
68 
69        DO ij   = 1, iip1
70          pks(   ij     )  =  xpn
71        ENDDO
72      endif
73     
74      if (pole_sud) then
75        DO  ij   = 1, iim
76          pps(ij) = aire(ij+ip1jm) * pks(ij+ip1jm )
77        ENDDO
78        xps      = SSUM(iim,pps,1) /apols
79 
80        DO ij   = 1, iip1
81          pks( ij+ip1jm )  =  xps
82        ENDDO
83      endif
84
85c
86c
87c    .... Calcul des coeff. alpha et beta  pour la couche l = llm ..
88c
89      DO     ij      = ijb,ije
90       alpha(ij,llm) = 0.
91       beta (ij,llm) = 1./ unpl2k
92      ENDDO
93c
94c     ... Calcul des coeff. alpha et beta  pour l = llm-1  a l = 2 ...
95c
96      DO l = llm -1 , 2 , -1
97c
98        DO ij = ijb, ije
99        dellta = p(ij,l)* unpl2k + p(ij,l+1)* ( beta(ij,l+1)-unpl2k )
100        alpha(ij,l)  = - p(ij,l+1) / dellta * alpha(ij,l+1)
101        beta (ij,l)  =   p(ij,l  ) / dellta   
102        ENDDO
103c
104      ENDDO
105c
106c  ***********************************************************************
107c     .....  Calcul de pk pour la couche 1 , pres du sol  ....
108c
109      DO   ij   = ijb, ije
110       pk(ij,1) = ( p(ij,1)*pks(ij) - 0.5*alpha(ij,2)*p(ij,2) )  /
111     *    (  p(ij,1)* (1.+kappa) + 0.5*( beta(ij,2)-unpl2k )* p(ij,2) )
112      ENDDO
113c
114c    ..... Calcul de pk(ij,l) , pour l = 2 a l = llm  ........
115c
116      DO l = 2, llm
117        DO   ij   = ijb, ije
118         pk(ij,l) = alpha(ij,l) + beta(ij,l) * pk(ij,l-1)
119        ENDDO
120      ENDDO
121c
122c
123c      CALL SCOPY   ( ngrid * llm, pk, 1, pkf, 1 )
124      pkf(ijb:ije,1:llm)=pk(ijb:ije,1:llm)
125     
126      jjb=jj_begin
127      jje=jj_end
128      CALL filtreg_p ( pkf,jjb,jje, jmp1, llm, 2, 1, .TRUE., 1 )
129     
130
131      RETURN
132      END
Note: See TracBrowser for help on using the repository browser.