source: LMDZ4/trunk/libf/phylmd/thermcell_closure.F90 @ 1046

Last change on this file since 1046 was 1026, checked in by lmdzadmin, 16 years ago

Modifs thermiques
FH/IM

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 2.4 KB
Line 
1      SUBROUTINE thermcell_closure(ngrid,nlay,r_aspect,ptimestep,rho,  &
2     &   zlev,lalim,alim_star,alim_star_tot,zmax_sec,wmax_sec,zmax,wmax,f,lev_out)
3
4!-------------------------------------------------------------------------
5!thermcell_closure: fermeture, determination de f
6!-------------------------------------------------------------------------
7      IMPLICIT NONE
8
9#include "iniprint.h"
10#include "thermcell.h"
11      INTEGER ngrid,nlay
12      INTEGER ig,k       
13      REAL r_aspect,ptimestep
14      integer lev_out                           ! niveau pour les print
15
16      INTEGER lalim(ngrid)
17      REAL alim_star(ngrid,nlay)
18      REAL alim_star_tot(ngrid)
19      REAL rho(ngrid,nlay)
20      REAL zlev(ngrid,nlay)
21      REAL zmax(ngrid),zmax_sec(ngrid)
22      REAL wmax(ngrid),wmax_sec(ngrid)
23      real zdenom
24
25      REAL alim_star2(ngrid)
26
27      REAL f(ngrid)
28
29      do ig=1,ngrid
30         alim_star2(ig)=0.
31      enddo
32      do ig=1,ngrid
33         if (alim_star(ig,1).LT.1.e-10) then
34            f(ig)=0.
35         else   
36             do k=1,lalim(ig)
37                alim_star2(ig)=alim_star2(ig)+alim_star(ig,k)**2  &
38     &                    /(rho(ig,k)*(zlev(ig,k+1)-zlev(ig,k)))
39             enddo
40             zdenom=max(500.,zmax(ig))*r_aspect*alim_star2(ig)
41             if (zdenom<1.e-14) then
42                print*,'ig=',ig
43                print*,'alim_star2',alim_star2(ig)
44                print*,'zmax',zmax(ig)
45                print*,'r_aspect',r_aspect
46                print*,'zdenom',zdenom
47                print*,'alim_star',alim_star(ig,:)
48                print*,'zmax_sec',zmax_sec(ig)
49                print*,'wmax_sec',wmax_sec(ig)
50                stop
51             endif
52             if ((zmax_sec(ig).gt.1.e-10).and.(iflag_thermals_ed.eq.0)) then
53             f(ig)=wmax_sec(ig)*alim_star_tot(ig)/(max(500.,zmax_sec(ig))*r_aspect  &
54     &             *alim_star2(ig))
55!            f(ig)=f(ig)+(f0(ig)-f(ig))*exp((-ptimestep/  &
56!    &                     zmax_sec(ig))*wmax_sec(ig))
57             print*,'closure dry',f(ig),wmax_sec(ig),alim_star_tot(ig),zmax_sec(ig)
58             else
59             f(ig)=wmax(ig)*alim_star_tot(ig)/zdenom
60!            f(ig)=f(ig)+(f0(ig)-f(ig))*exp((-ptimestep/  &
61!     &                     zmax(ig))*wmax(ig))
62             print*,'closure moist',f(ig),wmax(ig),alim_star_tot(ig),zmax(ig)
63             endif
64         endif
65!         f0(ig)=f(ig)
66      enddo
67      if (prt_level.ge.1) print*,'apres fermeture'
68
69!
70      return
71      end
Note: See TracBrowser for help on using the repository browser.