| 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 | if(prt_level.GE.10) 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 |
|---|
