source: trunk/LMDZ.GENERIC/libf/phystd/thermcell_closure.F90 @ 2613

Last change on this file since 2613 was 2177, checked in by aboissinot, 5 years ago

Cleanup thermal plums model subroutines
In thermcell_flux, "bidouilles" are modified:

  • now the plumes stop when the updraft fraction is greater than alpha_max
  • e > e_max is no longer permitted
  • b <= incoming mass flux is checked last
File size: 3.5 KB
Line 
1!
2!
3!
4SUBROUTINE thermcell_closure(ngrid,nlay,ptimestep,rho,zlev,                   &
5                             lmax,alim_star,zmin,zmax,wmax,f)
6     
7     
8!===============================================================================
9!  Purpose: fermeture, determination de f
10!
11! Modification 7 septembre 2009
12! 1. On enleve alim_star_tot des arguments pour le recalculer et etre ainis
13! coherent avec l'integrale au numerateur.
14! 2. On ne garde qu'une version des couples wmax,zmax et wmax_sec,zmax_sec
15! l'idee etant que le choix se fasse a l'appel de thermcell_closure
16! 3. Vectorisation en mettant les boucles en l a l'exterieur avec des if
17!===============================================================================
18     
19      USE thermcell_mod
20     
21      IMPLICIT NONE
22     
23     
24!===============================================================================
25! Declaration
26!===============================================================================
27     
28!     Inputs:
29!     -------
30     
31      INTEGER, INTENT(in) :: ngrid
32      INTEGER, INTENT(in) :: nlay
33      INTEGER, INTENT(in) :: lmax(ngrid)
34     
35      REAL, INTENT(in) :: ptimestep
36      REAL, INTENT(in) :: rho(ngrid,nlay)
37      REAL, INTENT(in) :: zlev(ngrid,nlay)
38      REAL, INTENT(in) :: alim_star(ngrid,nlay)
39      REAL, INTENT(in) :: zmin(ngrid)
40      REAL, INTENT(in) :: zmax(ngrid)
41      REAL, INTENT(in) :: wmax(ngrid)
42     
43!     Outputs:
44!     --------
45     
46      REAL, INTENT(out) :: f(ngrid)
47     
48!     Local:
49!     ------
50     
51      INTEGER ig, l
52      INTEGER llmax
53     
54      REAL alim_star_tot(ngrid)
55      REAL alim_star2(ngrid)
56      REAL plume_height(ngrid)
57     
58!===============================================================================
59! Initialization
60!===============================================================================
61     
62      alim_star2(:) = 0.
63      alim_star_tot(:) = 0.
64     
65      f(:) = 0.
66     
67      llmax = 1
68     
69      DO ig=1,ngrid
70         plume_height(ig) = zmax(ig) - zmin(ig)
71      ENDDO
72     
73!===============================================================================
74! Closure
75!===============================================================================
76     
77!-------------------------------------------------------------------------------
78! Indice vertical max atteint par les thermiques sur le domaine
79!-------------------------------------------------------------------------------
80     
81      DO ig=1,ngrid
82         IF (lmax(ig) > llmax) THEN
83            llmax = lmax(ig)
84         ENDIF
85      ENDDO
86     
87!-------------------------------------------------------------------------------
88! Calcul des integrales sur la verticale de alim_star et de alim_star^2/(rho dz)
89!-------------------------------------------------------------------------------
90     
91      DO l=1,llmax-1
92         DO ig=1,ngrid
93               alim_star2(ig) = alim_star2(ig) + alim_star(ig,l)**2           &
94               &              / (rho(ig,l) * (zlev(ig,l+1) - zlev(ig,l)))        ! => integration is ok because alim_star = a* dz
95               alim_star_tot(ig) = alim_star_tot(ig) + alim_star(ig,l)
96         ENDDO
97      ENDDO
98     
99      DO ig=1,ngrid
100         IF ((alim_star2(ig) > 0.).and.(plume_height(ig) > 0.)) THEN
101            f(ig) = wmax(ig) * alim_star_tot(ig)                              &  ! => normalization is ok
102            &     / (plume_height(ig) * r_aspect_thermals * alim_star2(ig))
103         ELSE
104            f(ig) = 0.
105         ENDIF
106      ENDDO
107     
108     
109RETURN
110END
Note: See TracBrowser for help on using the repository browser.