lwu.F90

Go to the documentation of this file.

SUBROUTINE lwu &
  &( kidia, kfdia, klon, klev &
  &, paer , pcco2, pdp , ppmb, pqof , ptave, pview, pwv &
  &, pabcu &
  &)

!**** *LWU* - LONGWAVE EFFECTIVE ABSORBER AMOUNTS

!     PURPOSE.
!     --------
!           COMPUTES ABSORBER AMOUNTS INCLUDING PRESSURE AND
!           TEMPERATURE EFFECTS

!**   INTERFACE.
!     ----------

!        EXPLICIT ARGUMENTS :
!        --------------------
!     ==== INPUTS ===
! PAER   : (KLON,6,KLEV)     ; OPTICAL THICKNESS OF THE AEROSOLS
! PCCO2  :                   ; CONCENTRATION IN CO2 (PA/PA)
! PDP    : (KLON,KLEV)       ; LAYER PRESSURE THICKNESS (PA)
! PPMB   : (KLON,KLEV+1)     ; HALF LEVEL PRESSURE
! PQOF   : (KLON,KLEV)       ; CONCENTRATION IN OZONE (PA/PA)
! PTAVE  : (KLON,KLEV)       ; TEMPERATURE
! PWV    : (KLON,KLEV)       ; SPECIFIC HUMIDITY PA/PA
! PVIEW  : (KLON)            ; COSECANT OF VIEWING ANGLE
!     ==== OUTPUTS ===
! PABCU  :(KLON,NUA,3*KLEV+1); EFFECTIVE ABSORBER AMOUNTS

!        IMPLICIT ARGUMENTS :   NONE
!        --------------------

!     METHOD.
!     -------

!          1. COMPUTES THE PRESSURE AND TEMPERATURE WEIGHTED AMOUNTS OF
!     ABSORBERS.

!     EXTERNALS.
!     ----------

!          NONE

!     REFERENCE.
!     ----------

!        SEE RADIATION'S PART OF THE MODEL'S DOCUMENTATION AND
!        ECMWF RESEARCH DEPARTMENT DOCUMENTATION OF THE IFS

!     AUTHOR.
!     -------
!        JEAN-JACQUES MORCRETTE  *ECMWF*

!     MODIFICATIONS.
!     --------------
!        ORIGINAL : 89-07-14
!        JJ Morcrette 97-04-18 Revised Continuum + Clean-up

!-----------------------------------------------------------------------

#include "tsmbkind.h"

USE yomcst   , ONLY : rg
USE yoesw    , ONLY : raer
USE yoelw    , ONLY : nsil     ,nua      ,ng1      ,ng1p1    ,&
            &alwt     ,blwt     ,ro3t     ,rt1      ,tref     ,&
            &rvgco2   ,rvgh2o   ,rvgo3
USE yoerdi   , ONLY : rch4     ,rn2o     ,rcfc11   ,rcfc12
USE yoerdu   , ONLY : r10e     ,repsco   ,repscq


IMPLICIT NONE


!     DUMMY INTEGER SCALARS
integer_m :: kfdia
integer_m :: kidia
integer_m :: klev
integer_m :: klon

!     DUMMY REAL SCALARS
real_b :: pcco2



!-----------------------------------------------------------------------

!*       0.1   ARGUMENTS
!              ---------

real_b :: paer(klon,6,klev), pdp(klon,klev)&
  &,  ppmb(klon,klev+1), pqof(klon,klev)&
  &,  ptave(klon,klev) , pview(klon),  pwv(klon,klev)

real_b :: pabcu(klon,nua,3*klev+1)

!-----------------------------------------------------------------------

!*       0.2   LOCAL ARRAYS
!              ------------
real_b :: zably(klon,7,3*klev+1)  , zdpm(klon,3*klev)&
  &,  zduc(klon, 3*klev+1)    , zfact(klon)&
  &,  zupm(klon,3*klev)
real_b :: zphio(klon),zpsc2(klon) , zpsc3(klon), zpsh1(klon)&
  &,  zpsh2(klon),zpsh3(klon) , zpsh4(klon), zpsh5(klon)&
  &,  zpsh6(klon),zpsio(klon) , ztcon(klon)&
  &,  zphm6(klon),zpsm6(klon) , zphn6(klon), zpsn6(klon)
real_b :: zssig(klon,3*klev+1)    , ztavi(klon)&
  &,  zuaer(klon,nsil)        , zxoz(klon) , zxwv(klon)

!     LOCAL INTEGER SCALARS
integer_m :: iae1, iae2, iae3, ic, icp1, ig1, ij, ijpn,&
             &ikip1, ikj, ikjp, ikjpn, ikjr, ikl, ja, jae, &
             &jk, jki, jkk, jl

!     LOCAL REAL SCALARS
real_b :: zalup, zcac8, zcah1, zcah2, zcah3, zcah4,&
          &zcah5, zcah6, zcbc8, zcbh1, zcbh2, zcbh3, &
          &zcbh4, zcbh5, zcbh6, zdiff, zdpmg, zdpmp0, &
          &zfppw, ztx, ztx2, zu6, zup, zupmco2, zupmg, &
          &zupmh2o, zupmo3, zzably


!-----------------------------------------------------------------------

!*         1.    INITIALIZATION
!                --------------

!-----------------------------------------------------------------------


!*         2.    PRESSURE OVER GAUSS SUB-LEVELS
!                ------------------------------

DO jl = kidia,kfdia
  zssig(jl, 1 ) = ppmb(jl,1) * 100._jprb
ENDDO

DO jk = 1 , klev
  ikj=(jk-1)*ng1p1+1
  ikjr = ikj
  ikjp = ikj + ng1p1
  DO jl = kidia,kfdia
    zssig(jl,ikjp)=ppmb(jl,jk+1)* 100._jprb
  ENDDO
  DO ig1=1,ng1
    ikj=ikj+1
    DO jl = kidia,kfdia
      zssig(jl,ikj)= (zssig(jl,ikjr) + zssig(jl,ikjp)) * _half_ &
       &+ rt1(ig1) * (zssig(jl,ikjp) - zssig(jl,ikjr)) * _half_
    ENDDO
  ENDDO
ENDDO

!-----------------------------------------------------------------------


!*         4.    PRESSURE THICKNESS AND MEAN PRESSURE OF SUB-LAYERS
!                --------------------------------------------------

DO jki=1,3*klev
  ikip1=jki+1
  DO jl = kidia,kfdia
    zupm(jl,jki)=(zssig(jl,jki)+zssig(jl,ikip1))*_half_
    zdpm(jl,jki)=(zssig(jl,jki)-zssig(jl,ikip1))/(10._jprb*rg)
  ENDDO
ENDDO

DO jk = 1 , klev
  ikl = klev+1 - jk
  DO jl = kidia,kfdia
    zxwv(jl) = max(pwv(jl,ikl) , repscq )
    zxoz(jl) = max(pqof(jl,ikl) / pdp(jl,ikl) , repsco )
  ENDDO
  ikj=(jk-1)*ng1p1+1
  ikjpn=ikj+ng1
  DO jkk=ikj,ikjpn
    DO jl = kidia,kfdia
      zdpmg = zdpm(jl,jkk)
      zdpmp0 = zdpmg / 101325._jprb
      zupmg = zupm(jl,jkk) * zdpmp0
      zupmco2 = ( zupm(jl,jkk) + rvgco2 ) * zdpmp0
      zupmh2o = ( zupm(jl,jkk) + rvgh2o ) * zdpmp0
      zupmo3  = ( zupm(jl,jkk) + rvgo3  ) * zdpmp0
      zduc(jl,jkk) = zdpmg
      zably(jl,6,jkk) = zxoz(jl) * zdpmg
      zably(jl,7,jkk) = zxoz(jl) * zupmo3
      zu6 = zxwv(jl) * zupmg
      zfppw = 1.6078_jprb * zxwv(jl) / (_one_+0.608_jprb*zxwv(jl))
      zably(jl,1,jkk)  = zxwv(jl) * zupmh2o
      zably(jl,5,jkk) = zu6 * zfppw
      zably(jl,4,jkk) = zu6 * (_one_-zfppw)
      zably(jl,3,jkk)  = pcco2 * zupmco2
      zably(jl,2,jkk)  = pcco2 * zdpmg
    ENDDO
  ENDDO
ENDDO

!-----------------------------------------------------------------------


!*         5.    CUMULATIVE ABSORBER AMOUNTS FROM TOP OF ATMOSPHERE
!                --------------------------------------------------

DO ja = 1, nua
  DO jl = kidia,kfdia
    pabcu(jl,ja,3*klev+1) = _zero_
  ENDDO
ENDDO

DO jk = 1 , klev
  ij=(jk-1)*ng1p1+1
  ijpn=ij+ng1
  ikl=klev+1-jk


!*         5.1  CUMULATIVE AEROSOL AMOUNTS FROM TOP OF ATMOSPHERE
!               --------------------------------------------------
! --            NB: 'PAER' AEROSOLS ARE ENTERED FROM TOP TO BOTTOM

  iae1=3*klev+1-ij
  iae2=3*klev+1-(ij+1)
  iae3=3*klev+1-ijpn
  DO jae=1,6
    DO jl = kidia,kfdia
      zuaer(jl,jae) =&
       &(raer(jae,1)*paer(jl,1,jk)+raer(jae,2)*paer(jl,2,jk)&
       &+raer(jae,3)*paer(jl,3,jk)+raer(jae,4)*paer(jl,4,jk)&
       &+raer(jae,5)*paer(jl,5,jk)+raer(jae,6)*paer(jl,5,jk))&
       &/(zduc(jl,iae1)+zduc(jl,iae2)+zduc(jl,iae3))
    ENDDO
  ENDDO



!*         5.2  INTRODUCES TEMPERATURE EFFECTS ON ABSORBER AMOUNTS
!               --------------------------------------------------

  DO jl = kidia,kfdia
    ztavi(jl)=ptave(jl,ikl)
    zfact(jl)=_one_-ztavi(jl)/296._jprb
    ztcon(jl)=exp(6.08_jprb*(296._jprb/ztavi(jl)-_one_))
!     ZTCON(JL)=EXP(6.08*ZFACT(JL))
    ztx=ztavi(jl)-tref
    ztx2=ztx*ztx
    zzably = zably(jl,1,iae1)+zably(jl,1,iae2)+zably(jl,1,iae3)
    zup=min( max( _half_*r10e*log( zzably ) + 5._jprb, _zero_), 6.0_jprb)
    zcah1=alwt(1,1)+zup*(alwt(1,2)+zup*(alwt(1,3)))
    zcbh1=blwt(1,1)+zup*(blwt(1,2)+zup*(blwt(1,3)))
    zpsh1(jl)=exp( zcah1 * ztx + zcbh1 * ztx2 )
    zcah2=alwt(2,1)+zup*(alwt(2,2)+zup*(alwt(2,3)))
    zcbh2=blwt(2,1)+zup*(blwt(2,2)+zup*(blwt(2,3)))
    zpsh2(jl)=exp( zcah2 * ztx + zcbh2 * ztx2 )
    zcah3=alwt(3,1)+zup*(alwt(3,2)+zup*(alwt(3,3)))
    zcbh3=blwt(3,1)+zup*(blwt(3,2)+zup*(blwt(3,3)))
    zpsh3(jl)=exp( zcah3 * ztx + zcbh3 * ztx2 )
    zcah4=alwt(4,1)+zup*(alwt(4,2)+zup*(alwt(4,3)))
    zcbh4=blwt(4,1)+zup*(blwt(4,2)+zup*(blwt(4,3)))
    zpsh4(jl)=exp( zcah4 * ztx + zcbh4 * ztx2 )
    zcah5=alwt(5,1)+zup*(alwt(5,2)+zup*(alwt(5,3)))
    zcbh5=blwt(5,1)+zup*(blwt(5,2)+zup*(blwt(5,3)))
    zpsh5(jl)=exp( zcah5 * ztx + zcbh5 * ztx2 )
    zcah6=alwt(6,1)+zup*(alwt(6,2)+zup*(alwt(6,3)))
    zcbh6=blwt(6,1)+zup*(blwt(6,2)+zup*(blwt(6,3)))
    zpsh6(jl)=exp( zcah6 * ztx + zcbh6 * ztx2 )
    zphm6(jl)=exp(-5.81e-4_jprb * ztx - 1.13e-6_jprb * ztx2 )
    zpsm6(jl)=exp(-5.57e-4_jprb * ztx - 3.30e-6_jprb * ztx2 )
    zphn6(jl)=exp(-3.46e-5_jprb * ztx + 2.05e-7_jprb * ztx2 )
    zpsn6(jl)=exp( 3.70e-3_jprb * ztx - 2.30e-6_jprb * ztx2 )
  ENDDO

  DO jl = kidia,kfdia
    ztavi(jl)=ptave(jl,ikl)
    ztx=ztavi(jl)-tref
    ztx2=ztx*ztx
    zzably = zably(jl,3,iae1)+zably(jl,3,iae2)+zably(jl,3,iae3)
    zalup = r10e * log( zzably )
    zup   = max( _zero_ , 5.0_jprb + _half_ * zalup )
    zpsc2(jl) = (ztavi(jl)/tref) ** zup
    zcac8=alwt(8,1)+zup*(alwt(8,2)+zup*(alwt(8,3)))
    zcbc8=blwt(8,1)+zup*(blwt(8,2)+zup*(blwt(8,3)))
    zpsc3(jl)=exp( zcac8 * ztx + zcbc8 * ztx2 )
    zphio(jl) = exp( ro3t(1) * ztx + ro3t(2) * ztx2)
    zpsio(jl) = exp( _two_* (ro3t(3)*ztx+ro3t(4)*ztx2))
  ENDDO

  DO jkk=ij,ijpn
    ic=3*klev+1-jkk
    icp1=ic+1
    DO jl = kidia,kfdia
      zdiff = pview(jl)
!- H2O continuum      
      pabcu(jl,10,ic)=pabcu(jl,10,icp1)+ zably(jl,4,ic)          *zdiff
      pabcu(jl,11,ic)=pabcu(jl,11,icp1)+ zably(jl,5,ic)*ztcon(jl)*zdiff
!- O3      
      pabcu(jl,12,ic)=pabcu(jl,12,icp1)+ zably(jl,6,ic)*zphio(jl)*zdiff
      pabcu(jl,13,ic)=pabcu(jl,13,icp1)+ zably(jl,7,ic)*zpsio(jl)*zdiff
!- CO2
      pabcu(jl,7,ic)=pabcu(jl,7,icp1)+ zably(jl,3,ic)*zpsc2(jl)*zdiff
      pabcu(jl,8,ic)=pabcu(jl,8,icp1)+ zably(jl,3,ic)*zpsc3(jl)*zdiff
      pabcu(jl,9,ic)=pabcu(jl,9,icp1)+ zably(jl,3,ic)*zpsc3(jl)*zdiff
!- H2O
      pabcu(jl,1,ic)=pabcu(jl,1,icp1)+ zably(jl,1,ic)*zpsh1(jl)
      pabcu(jl,2,ic)=pabcu(jl,2,icp1)+ zably(jl,1,ic)*zpsh2(jl)
      pabcu(jl,3,ic)=pabcu(jl,3,icp1)+ zably(jl,1,ic)*zpsh5(jl)*zdiff
      pabcu(jl,4,ic)=pabcu(jl,4,icp1)+ zably(jl,1,ic)*zpsh3(jl)
      pabcu(jl,5,ic)=pabcu(jl,5,icp1)+ zably(jl,1,ic)*zpsh4(jl)
      pabcu(jl,6,ic)=pabcu(jl,6,icp1)+ zably(jl,1,ic)*zpsh6(jl)*zdiff
!- aerosols
      pabcu(jl,14,ic)=pabcu(jl,14,icp1)+ zuaer(jl,1)    *zduc(jl,ic)*zdiff
      pabcu(jl,15,ic)=pabcu(jl,15,icp1)+ zuaer(jl,2)    *zduc(jl,ic)*zdiff
      pabcu(jl,16,ic)=pabcu(jl,16,icp1)+ zuaer(jl,3)    *zduc(jl,ic)*zdiff
      pabcu(jl,17,ic)=pabcu(jl,17,icp1)+ zuaer(jl,4)    *zduc(jl,ic)*zdiff
      pabcu(jl,18,ic)=pabcu(jl,18,icp1)+ zuaer(jl,5)    *zduc(jl,ic)*zdiff
!- CH4
      pabcu(jl,19,ic)=pabcu(jl,19,icp1)&
       &+ zably(jl,2,ic)*rch4/pcco2*zphm6(jl)*zdiff
      pabcu(jl,20,ic)=pabcu(jl,20,icp1)&
       &+ zably(jl,3,ic)*rch4/pcco2*zpsm6(jl)*zdiff
!- N2O
      pabcu(jl,21,ic)=pabcu(jl,21,icp1)&
       &+ zably(jl,2,ic)*rn2o/pcco2*zphn6(jl)*zdiff
      pabcu(jl,22,ic)=pabcu(jl,22,icp1)&
       &+ zably(jl,3,ic)*rn2o/pcco2*zpsn6(jl)*zdiff
!- CFC11
      pabcu(jl,23,ic)=pabcu(jl,23,icp1)&
       &+ zably(jl,2,ic)*rcfc11/pcco2        *zdiff
!- CFC12
      pabcu(jl,24,ic)=pabcu(jl,24,icp1)&
       &+ zably(jl,2,ic)*rcfc12/pcco2        *zdiff
    ENDDO
  ENDDO

ENDDO
!      print *,'END OF LWU'

!-----------------------------------------------------------------------

RETURN
END SUBROUTINE lwu