/[lmdze]/trunk/libf/phylmd/Thermcell/thermcell.f90
ViewVC logotype

Diff of /trunk/libf/phylmd/Thermcell/thermcell.f90

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

revision 10 by guez, Fri Apr 18 14:45:53 2008 UTC revision 32 by guez, Tue Apr 6 17:52:58 2010 UTC
# Line 925  c             wqd(ig,k)=fm(ig,k)*0.5*(q( Line 925  c             wqd(ig,k)=fm(ig,k)*0.5*(q(
925    
926        return        return
927        end        end
       subroutine dvthermcell(ngrid,nlay,ptimestep,fm,entr,masse  
      .    ,fraca,larga  
      .    ,u,v,du,dv,ua,va)  
       use dimens_m  
       use dimphy  
       implicit none  
   
 c=======================================================================  
 c  
 c   Calcul du transport verticale dans la couche limite en presence  
 c   de "thermiques" explicitement representes  
 c   calcul du dq/dt une fois qu'on connait les ascendances  
 c  
 c=======================================================================  
   
   
       integer ngrid,nlay  
   
       real ptimestep  
       real masse(ngrid,nlay),fm(ngrid,nlay+1)  
       real fraca(ngrid,nlay+1)  
       real larga(ngrid)  
       real entr(ngrid,nlay)  
       real u(ngrid,nlay)  
       real ua(ngrid,nlay)  
       real du(ngrid,nlay)  
       real v(ngrid,nlay)  
       real va(ngrid,nlay)  
       real dv(ngrid,nlay)  
   
       real qa(klon,klev),detr(klon,klev)  
       real wvd(klon,klev+1),wud(klon,klev+1)  
       real gamma0,gamma(klon,klev+1)  
       real dua,dva  
       integer iter  
   
       integer ig,k  
   
 c   calcul du detrainement  
928    
       do k=1,nlay  
          do ig=1,ngrid  
             detr(ig,k)=fm(ig,k)-fm(ig,k+1)+entr(ig,k)  
          enddo  
       enddo  
   
 c   calcul de la valeur dans les ascendances  
       do ig=1,ngrid  
          ua(ig,1)=u(ig,1)  
          va(ig,1)=v(ig,1)  
       enddo  
   
       do k=2,nlay  
          do ig=1,ngrid  
             if ((fm(ig,k+1)+detr(ig,k))*ptimestep.gt.  
      s         1.e-5*masse(ig,k)) then  
 c   On itère sur la valeur du coeff de freinage.  
 c              gamma0=rho(ig,k)*(zlev(ig,k+1)-zlev(ig,k))  
                gamma0=masse(ig,k)  
      s         *sqrt( 0.5*(fraca(ig,k+1)+fraca(ig,k)) )  
      s         *0.5/larga(ig)  
 c              gamma0=0.  
 c   la première fois on multiplie le coefficient de freinage  
 c   par le module du vent dans la couche en dessous.  
                dua=ua(ig,k-1)-u(ig,k-1)  
                dva=va(ig,k-1)-v(ig,k-1)  
                do iter=1,5  
                   gamma(ig,k)=gamma0*sqrt(dua**2+dva**2)  
                   ua(ig,k)=(fm(ig,k)*ua(ig,k-1)  
      s               +(entr(ig,k)+gamma(ig,k))*u(ig,k))  
      s               /(fm(ig,k+1)+detr(ig,k)+gamma(ig,k))  
                   va(ig,k)=(fm(ig,k)*va(ig,k-1)  
      s               +(entr(ig,k)+gamma(ig,k))*v(ig,k))  
      s               /(fm(ig,k+1)+detr(ig,k)+gamma(ig,k))  
 c                 print*,k,ua(ig,k),va(ig,k),u(ig,k),v(ig,k),dua,dva  
                   dua=ua(ig,k)-u(ig,k)  
                   dva=va(ig,k)-v(ig,k)  
                enddo  
             else  
                ua(ig,k)=u(ig,k)  
                va(ig,k)=v(ig,k)  
                gamma(ig,k)=0.  
             endif  
          enddo  
       enddo  
   
       do k=2,nlay  
          do ig=1,ngrid  
             wud(ig,k)=fm(ig,k)*u(ig,k)  
             wvd(ig,k)=fm(ig,k)*v(ig,k)  
          enddo  
       enddo  
       do ig=1,ngrid  
          wud(ig,1)=0.  
          wud(ig,nlay+1)=0.  
          wvd(ig,1)=0.  
          wvd(ig,nlay+1)=0.  
       enddo  
   
       do k=1,nlay  
          do ig=1,ngrid  
             du(ig,k)=((detr(ig,k)+gamma(ig,k))*ua(ig,k)  
      s               -(entr(ig,k)+gamma(ig,k))*u(ig,k)  
      s               -wud(ig,k)+wud(ig,k+1))  
      s               /masse(ig,k)  
             dv(ig,k)=((detr(ig,k)+gamma(ig,k))*va(ig,k)  
      s               -(entr(ig,k)+gamma(ig,k))*v(ig,k)  
      s               -wvd(ig,k)+wvd(ig,k+1))  
      s               /masse(ig,k)  
          enddo  
       enddo  
   
       return  
       end  
929        subroutine dqthermcell2(ngrid,nlay,ptimestep,fm,entr,masse,frac        subroutine dqthermcell2(ngrid,nlay,ptimestep,fm,entr,masse,frac
930       .    ,q,dq,qa)       .    ,q,dq,qa)
931        use dimens_m        use dimens_m

Legend:
Removed from v.10  
changed lines
  Added in v.32

  ViewVC Help
Powered by ViewVC 1.1.21