/[lmdze]/trunk/dyn3d/exner_hyb.f
ViewVC logotype

Diff of /trunk/dyn3d/exner_hyb.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 3 by guez, Wed Feb 27 13:16:39 2008 UTC revision 10 by guez, Fri Apr 18 14:45:53 2008 UTC
# Line 8  contains Line 8  contains
8    
9      ! From dyn3d/exner_hyb.F, v 1.1.1.1 2004/05/19 12:53:07      ! From dyn3d/exner_hyb.F, v 1.1.1.1 2004/05/19 12:53:07
10    
11      ! Auteurs :  P. Le Van, F. Hourdin.      ! Authors : P. Le Van, F. Hourdin
12    
13      ! Calcule la fonction d'Exner :      ! Calcule la fonction d'Exner :
14      ! pk = Cp * p ** kappa      ! pk = Cp * p ** kappa
15      ! aux milieux des couches.      ! aux milieux des "llm" couches.
16      ! Pk(l) est calculé aux milieux des couches "l", entre les pressions      ! "Pk(l)" est calculé au milieu de la couche "l", entre les pressions
17      ! "p(l)" et "p(l+1)", définies aux interfaces des "llm" couches.      ! "p(l)" et "p(l+1)", définies aux interfaces des couches.
18    
19      ! Au sommet de l'atmosphère :      ! Au sommet de l'atmosphère :
20      ! p(llm+1) = 0.      ! p(llm+1) = 0.
21      ! et ps et pks sont la pression et la fonction d'Exner au sol.      ! "ps" et "pks" sont la pression et la fonction d'Exner au sol.
22    
23      ! À partir des relations :      ! À partir des relations :
24    
25      !   -------- z      !   -------- z
26      !(1) p*dz(pk) = kappa *pk*dz(p)      !(1) p*dz(pk) = kappa * pk * dz(p)
27    
28      !(2) pk(l) = alpha(l)+ beta(l)*pk(l-1)      !(2) pk(l) = alpha(l)+ beta(l) * pk(l-1)
29    
30      ! (voir note de F. Hourdin), on determine successivement, du haut      ! (voir note de F. Hourdin), on détermine successivement, du haut
31      ! vers le bas des couches, les coefficients :      ! vers le bas des couches, les coefficients :
32      ! alpha(llm), beta(llm)..., alpha(l), beta(l)..., alpha(2), beta(2)      ! alpha(llm), beta(llm)..., alpha(l), beta(l)..., alpha(2), beta(2)
33      ! puis "pk(ij, 1)".      ! puis "pk(:, 1)".
34      ! Ensuite, on calcule, du bas vers le haut des couches, "pk(ij, l)"      ! Ensuite, on calcule, du bas vers le haut des couches, "pk(:, l)"
35      ! donné par la relation (2), pour l = 2 à l = llm.      ! donné par la relation (2), pour l = 2 à l = llm.
36    
37      use dimens_m, only: iim, jjm, llm      use dimens_m, only: iim, jjm, llm
# Line 49  contains Line 49  contains
49      ! Variables locales      ! Variables locales
50    
51      real alpha((iim + 1) * (jjm + 1), llm), beta((iim + 1) * (jjm + 1), llm)      real alpha((iim + 1) * (jjm + 1), llm), beta((iim + 1) * (jjm + 1), llm)
52      INTEGER l, ij      INTEGER l
53      REAL unpl2k, dellta      REAL unpl2k, dellta((iim + 1) * (jjm + 1))
54    
55      REAL ppn(iim), pps(iim)      REAL ppn(iim), pps(iim)
     REAL xpn, xps  
     REAL SSUM  
56    
57      !-------------------------------------      !-------------------------------------
58    
59      pks(:) = cpp * (ps(:) / preff)**kappa      pks = cpp * (ps / preff)**kappa
60      ppn(:) = aire_2d(:iim, 1) * pks(:iim)      ppn = aire_2d(:iim, 1) * pks(:iim)
61      pps(:) = aire_2d(:iim, jjm + 1) &      pps = aire_2d(:iim, jjm + 1) &
62           * pks(1 + (iim + 1) * jjm: iim + (iim + 1) * jjm)           * pks(1 + (iim + 1) * jjm: iim + (iim + 1) * jjm)
63      xpn = SSUM(iim, ppn, 1) /apoln      pks(:iim + 1) = SUM(ppn) /apoln
64      xps = SSUM(iim, pps, 1) /apols      pks(1+(iim + 1) * jjm:) = SUM(pps) /apols
     pks(:iim + 1) = xpn  
     pks(1+(iim + 1) * jjm:) = xps  
65    
66      unpl2k = 1. + 2 * kappa      unpl2k = 1. + 2 * kappa
67    
68      ! Calcul des coeff. alpha et beta  pour la couche l = llm      ! Calcul des coefficients alpha et beta pour la couche l = llm :
69      alpha(:, llm) = 0.      alpha(:, llm) = 0.
70      beta (:, llm) = 1./ unpl2k      beta(:, llm) = 1./ unpl2k
71    
72      ! Calcul des coeff. alpha et beta  pour l = llm-1  à l = 2      ! Calcul des coefficients alpha et beta pour l = llm-1 à l = 2 :
73      DO l = llm -1 , 2 , -1      DO l = llm - 1, 2, -1
74         DO ij = 1, (iim + 1) * (jjm + 1)         dellta = p(:, l) * unpl2k + p(:, l+1) * (beta(:, l+1) - unpl2k)
75            dellta = p(ij, l)* unpl2k + p(ij, l+1)* ( beta(ij, l+1)-unpl2k )         alpha(:, l) = - p(:, l+1) / dellta * alpha(:, l+1)
76            alpha(ij, l)  = - p(ij, l+1) / dellta * alpha(ij, l+1)         beta(:, l) = p(:, l) / dellta  
           beta (ij, l)  =   p(ij, l  ) / dellta    
        ENDDO  
77      ENDDO      ENDDO
78    
79      ! Calcul de pk pour la couche 1, près du sol :      ! Calcul de pk pour la couche 1, près du sol :
80      pk(:, 1) = (p(:, 1) * pks(:) - 0.5 * alpha(:, 2) * p(:, 2))  &      pk(:, 1) = (p(:, 1) * pks - 0.5 * alpha(:, 2) * p(:, 2))  &
81           / (p(:, 1) * (1. + kappa) + 0.5 * (beta(:, 2) - unpl2k) * p(:, 2))           / (p(:, 1) * (1. + kappa) + 0.5 * (beta(:, 2) - unpl2k) * p(:, 2))
82    
83      ! Calcul de pk(ij, l) , pour l = 2 à l = llm      ! Calcul de pk(:, l) pour l = 2 à l = llm :
84      DO l = 2, llm      DO l = 2, llm
85         DO   ij   = 1, (iim + 1) * (jjm + 1)         pk(:, l) = alpha(:, l) + beta(:, l) * pk(:, l-1)
           pk(ij, l) = alpha(ij, l) + beta(ij, l) * pk(ij, l-1)  
        ENDDO  
86      ENDDO      ENDDO
87    
88      if (present(pkf)) then      if (present(pkf)) then
89         pkf(:, :) = pk(:, :)         pkf = pk
90         CALL filtreg(pkf, jjm + 1, llm, 2, 1, .TRUE., 1)         CALL filtreg(pkf, jjm + 1, llm, 2, 1, .TRUE., 1)
91      end if      end if
92    
Legend:
Removed from v.3  
changed lines
  Added in v.10
  ViewVC Help
Powered by ViewVC 1.1.21