| 1 | SUBROUTINE thermcell_closure(ngrid,nlay,r_aspect,ptimestep,rho, & |
|---|
| 2 | & zlev,lalim,alim_star,zmax_sec,wmax_sec,zmax,wmax,f,f0,lev_out) |
|---|
| 3 | |
|---|
| 4 | !------------------------------------------------------------------------- |
|---|
| 5 | !thermcell_closure: fermeture, determination de f |
|---|
| 6 | !------------------------------------------------------------------------- |
|---|
| 7 | IMPLICIT NONE |
|---|
| 8 | |
|---|
| 9 | #include "iniprint.h" |
|---|
| 10 | INTEGER ngrid,nlay |
|---|
| 11 | INTEGER ig,k |
|---|
| 12 | REAL r_aspect,ptimestep |
|---|
| 13 | integer lev_out ! niveau pour les print |
|---|
| 14 | |
|---|
| 15 | INTEGER lalim(ngrid) |
|---|
| 16 | REAL alim_star(ngrid,nlay) |
|---|
| 17 | REAL rho(ngrid,nlay) |
|---|
| 18 | REAL zlev(ngrid,nlay) |
|---|
| 19 | REAL zmax(ngrid),zmax_sec(ngrid) |
|---|
| 20 | REAL wmax(ngrid),wmax_sec(ngrid) |
|---|
| 21 | real zdenom |
|---|
| 22 | |
|---|
| 23 | REAL alim_star2(ngrid) |
|---|
| 24 | |
|---|
| 25 | REAL f(ngrid) |
|---|
| 26 | REAL f0(ngrid) |
|---|
| 27 | |
|---|
| 28 | do ig=1,ngrid |
|---|
| 29 | alim_star2(ig)=0. |
|---|
| 30 | enddo |
|---|
| 31 | do ig=1,ngrid |
|---|
| 32 | if (alim_star(ig,1).LT.1.e-10) then |
|---|
| 33 | f(ig)=0. |
|---|
| 34 | else |
|---|
| 35 | do k=1,lalim(ig) |
|---|
| 36 | alim_star2(ig)=alim_star2(ig)+alim_star(ig,k)**2 & |
|---|
| 37 | & /(rho(ig,k)*(zlev(ig,k+1)-zlev(ig,k))) |
|---|
| 38 | enddo |
|---|
| 39 | zdenom=max(500.,zmax(ig))*r_aspect*alim_star2(ig) |
|---|
| 40 | if (zdenom<1.e-14) then |
|---|
| 41 | print*,'ig=',ig |
|---|
| 42 | print*,'alim_star2',alim_star2(ig) |
|---|
| 43 | print*,'zmax',zmax(ig) |
|---|
| 44 | print*,'r_aspect',r_aspect |
|---|
| 45 | print*,'zdenom',zdenom |
|---|
| 46 | print*,'alim_star',alim_star(ig,:) |
|---|
| 47 | print*,'zmax_sec',zmax_sec(ig) |
|---|
| 48 | print*,'wmax_sec',wmax_sec(ig) |
|---|
| 49 | stop |
|---|
| 50 | endif |
|---|
| 51 | if ((zmax_sec(ig).gt.1.e-10).and.(1.eq.1)) then |
|---|
| 52 | f(ig)=wmax_sec(ig)/(max(500.,zmax_sec(ig))*r_aspect & |
|---|
| 53 | & *alim_star2(ig)) |
|---|
| 54 | f(ig)=f(ig)+(f0(ig)-f(ig))*exp((-ptimestep/ & |
|---|
| 55 | & zmax_sec(ig))*wmax_sec(ig)) |
|---|
| 56 | else |
|---|
| 57 | f(ig)=wmax(ig)/zdenom |
|---|
| 58 | f(ig)=f(ig)+(f0(ig)-f(ig))*exp((-ptimestep/ & |
|---|
| 59 | & zmax(ig))*wmax(ig)) |
|---|
| 60 | endif |
|---|
| 61 | endif |
|---|
| 62 | f0(ig)=f(ig) |
|---|
| 63 | enddo |
|---|
| 64 | if (prt_level.ge.1) print*,'apres fermeture' |
|---|
| 65 | |
|---|
| 66 | ! |
|---|
| 67 | return |
|---|
| 68 | end |
|---|
