Changeset 1715


Ignore:
Timestamp:
Jun 8, 2017, 1:03:06 AM (8 years ago)
Author:
jvatant
Message:

Cleaning of the radiative code + coherence between optci and optcv :

+ Get rid of the actually dummy and confusing extra level L_LEVELS+1 for aerosols
+ Harmonisation between optci and optcv
+ Get rid of tauref (called nowhere) and csqtly of ini_radcommon subroutine

JVO

Location:
trunk/LMDZ.GENERIC/libf/phystd
Files:
5 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.GENERIC/libf/phystd/callcorrk.F90

    r1714 r1715  
    123123      REAL*8 taucumi(L_LEVELS,L_NSPECTI,L_NGAUSS)
    124124
    125       REAL*8 tauaero(L_LEVELS+1,naerkind)
     125      REAL*8 tauaero(L_LEVELS,naerkind)
    126126      REAL*8 nfluxtopv,nfluxtopi,nfluxtop,fluxtopvdn
    127127      REAL*8 nfluxoutv_nu(L_NSPECTV)                 ! Outgoing band-resolved VI flux at TOA (W/m2).
    … …  
    189189
    190190        ! test on allocated necessary because of CLFvarying (two calls to callcorrk in physiq)
    191         if(.not.allocated(QXVAER)) allocate(QXVAER(L_LEVELS+1,L_NSPECTV,naerkind))
    192         if(.not.allocated(QSVAER)) allocate(QSVAER(L_LEVELS+1,L_NSPECTV,naerkind))
    193         if(.not.allocated(GVAER)) allocate(GVAER(L_LEVELS+1,L_NSPECTV,naerkind))
    194         if(.not.allocated(QXIAER)) allocate(QXIAER(L_LEVELS+1,L_NSPECTI,naerkind))
    195         if(.not.allocated(QSIAER)) allocate(QSIAER(L_LEVELS+1,L_NSPECTI,naerkind))
    196         if(.not.allocated(GIAER)) allocate(GIAER(L_LEVELS+1,L_NSPECTI,naerkind))
     191        if(.not.allocated(QXVAER)) allocate(QXVAER(L_LEVELS,L_NSPECTV,naerkind))
     192        if(.not.allocated(QSVAER)) allocate(QSVAER(L_LEVELS,L_NSPECTV,naerkind))
     193        if(.not.allocated(GVAER)) allocate(GVAER(L_LEVELS,L_NSPECTV,naerkind))
     194        if(.not.allocated(QXIAER)) allocate(QXIAER(L_LEVELS,L_NSPECTI,naerkind))
     195        if(.not.allocated(QSIAER)) allocate(QSIAER(L_LEVELS,L_NSPECTI,naerkind))
     196        if(.not.allocated(GIAER)) allocate(GIAER(L_LEVELS,L_NSPECTI,naerkind))
    197197
    198198         !!! ALLOCATED instances are necessary because of CLFvarying (strategy to call callcorrk twice in physiq...)
    … …  
    323323     end do !iaer=1,naerkind.
    324324
     325
    325326      ! How much light do we get ?
    326327      do nw=1,L_NSPECTV
    … …  
    439440            ! Test / Correct for freaky s. s. albedo values.
    440441            do iaer=1,naerkind
    441                do k=1,L_LEVELS+1
     442               do k=1,L_LEVELS
    442443
    443444                  do nw=1,L_NSPECTV
    … …  
    482483            ! boundary conditions
    483484            tauaero(1,iaer)          = tauaero(2,iaer)
    484             tauaero(L_LEVELS+1,iaer) = tauaero(L_LEVELS,iaer)
    485485            !tauaero(1,iaer)          = 0.
    486             !tauaero(L_LEVELS+1,iaer) = 0.
    487486           
    488487         end do ! naerkind
    … …  
    574573
    575574      if(kastprof)then
    576      
    577          if(.not.global1d)then ! garde-fou/safeguard added by MT (to be removed in the future)
    578             write(*,*) 'You have to fix mu0, '
    579             write(*,*) 'the cosinus of the solar angle'
    580             stop
    581          endif
    582          
     575
    583576         ! Initial values equivalent to mugaz.
    584577         DO l=1,nlayer

  • trunk/LMDZ.GENERIC/libf/phystd/inifis_mod.F90

    r1682 r1715  
    1111  use init_print_control_mod, only: init_print_control
    1212  use radinc_h, only: ini_radinc_h, naerkind
    13   use radcommon_h, only: ini_radcommon_h
    1413  use datafile_mod, only: datadir
    1514  use comdiurn_h, only: sinlat, coslat, sinlon, coslon
    … …  
    762761  call ini_radinc_h(nlayer)
    763762 
    764   ! allocate "radcommon_h" arrays
    765   call ini_radcommon_h()
    766 
    767763  ! allocate "comsoil_h" arrays
    768764  call ini_comsoil_h(ngrid)

  • trunk/LMDZ.GENERIC/libf/phystd/optci.F90

    r1397 r1715  
    3232
    3333  real*8 DTAUI(L_NLAYRAD,L_NSPECTI,L_NGAUSS)
    34   real*8 DTAUKI(L_LEVELS+1,L_NSPECTI,L_NGAUSS)
     34  real*8 DTAUKI(L_LEVELS,L_NSPECTI,L_NGAUSS)
    3535  real*8 TAUI(L_NLEVRAD,L_NSPECTI,L_NGAUSS)
    3636  real*8 TAUCUMI(L_LEVELS,L_NSPECTI,L_NGAUSS)
    … …  
    4242
    4343  ! for aerosols
    44   real*8  QXIAER(L_LEVELS+1,L_NSPECTI,NAERKIND)
    45   real*8  QSIAER(L_LEVELS+1,L_NSPECTI,NAERKIND)
    46   real*8  GIAER(L_LEVELS+1,L_NSPECTI,NAERKIND)
    47   real*8  TAUAERO(L_LEVELS+1,NAERKIND)
    48   real*8  TAUAEROLK(L_LEVELS+1,L_NSPECTI,NAERKIND)
     44  real*8  QXIAER(L_LEVELS,L_NSPECTI,NAERKIND)
     45  real*8  QSIAER(L_LEVELS,L_NSPECTI,NAERKIND)
     46  real*8  GIAER(L_LEVELS,L_NSPECTI,NAERKIND)
     47  real*8  TAUAERO(L_LEVELS,NAERKIND)
     48  real*8  TAUAEROLK(L_LEVELS,L_NSPECTI,NAERKIND)
    4949  real*8  TAEROS(L_LEVELS,L_NSPECTI,NAERKIND)
    5050
    … …  
    135135  end do                    ! levels
    136136
    137 
     137  ! Spectral dependance of aerosol absorption
    138138  do iaer=1,naerkind
    139139     DO NW=1,L_NSPECTI
    … …  
    147147
    148148     do K=2,L_LEVELS
    149 
    150 ! continuum absorption
    151         DCONT = 0.0d0
     149     
     150        DAERO=SUM(TAEROS(K,NW,1:naerkind)) ! aerosol absorption
     151
     152        DCONT = 0.0d0 ! continuum absorption
    152153
    153154        if(continuum.and.(.not.graybody))then
    … …  
    225226        endif
    226227
    227 ! aerosol absorption
    228         DAERO=SUM(TAEROS(K,NW,1:naerkind))
    229 
    230228        do ng=1,L_NGAUSS-1
    231229
    … …  
    286284  end do
    287285
    288   DTAUKI(L_LEVELS+1,1:L_NSPECTI,1:L_NGAUSS)=0.d0
    289286  !=======================================================================
    290287  !     Now the full treatment for the layers, where besides the opacity
    … …  
    293290  do iaer=1,naerkind
    294291    DO NW=1,L_NSPECTI
    295      DO K=2,L_LEVELS+1
    296            TAUAEROLK(K,NW,IAER) = TAUAERO(K,IAER)*QSIAER(K,NW,IAER)
     292     DO K=2,L_LEVELS
     293           TAUAEROLK(K,NW,IAER) = TAUAERO(K,IAER)*QSIAER(K,NW,IAER) ! effect of scattering albedo
    297294     ENDDO
    298295    ENDDO
    … …  
    300297 
    301298  DO NW=1,L_NSPECTI
    302      DO L=1,L_NLAYRAD
     299     DO L=1,L_NLAYRAD-1
    303300        K              = 2*L+1
    304301        btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind)) + SUM(TAUAEROLK(K+1,NW,1:naerkind))
    305302     END DO ! L vertical loop
     303     
     304     ! Last level
     305     L           = L_NLAYRAD
     306     K           = 2*L+1   
     307     btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind))
     308     
    306309  END DO                    ! NW spectral loop
    307310 
    … …  
    309312  DO NW=1,L_NSPECTI
    310313     NG = L_NGAUSS
    311      DO L=1,L_NLAYRAD
     314     DO L=1,L_NLAYRAD-1
    312315
    313316        K              = 2*L+1
    … …  
    334337
    335338     END DO ! L vertical loop
     339     
     340     ! Last level
     341     
     342     L              = L_NLAYRAD
     343     K              = 2*L+1
     344     DTAUI(L,nw,ng) = DTAUKI(K,NW,NG) ! + 1.e-50
     345
     346     atemp = 0.
     347     if(DTAUI(L,NW,NG) .GT. 1.0D-9) then
     348        do iaer=1,naerkind
     349           atemp = atemp + GIAER(K,NW,IAER)   * TAUAEROLK(K,NW,IAER)
     350        end do
     351        WBARI(L,nw,ng) = btemp(L,nw)  / DTAUI(L,NW,NG)
     352     else
     353        WBARI(L,nw,ng) = 0.0D0
     354        DTAUI(L,NW,NG) = 1.0D-9
     355     endif
     356
     357     if(btemp(L,nw) .GT. 0.0d0) then
     358        cosbi(L,NW,NG) = atemp/btemp(L,nw)
     359     else
     360        cosbi(L,NW,NG) = 0.0D0
     361     end if
     362     
    336363
    337364     ! Now the other Gauss points, if needed.
    … …  
    340367        IF(TAUGSURF(NW,NG) .gt. TLIMIT) THEN
    341368
    342            DO L=1,L_NLAYRAD
     369           DO L=1,L_NLAYRAD-1
    343370              K              = 2*L+1
    344371              DTAUI(L,nw,ng) = DTAUKI(K,NW,NG)+DTAUKI(K+1,NW,NG)! + 1.e-50
    … …  
    355382              cosbi(L,NW,NG) = cosbi(L,NW,L_NGAUSS)
    356383           END DO ! L vertical loop
     384           
     385           ! Last level
     386           L              = L_NLAYRAD
     387           K              = 2*L+1
     388           DTAUI(L,nw,ng) = DTAUKI(K,NW,NG)! + 1.e-50
     389
     390           if(DTAUI(L,NW,NG) .GT. 1.0D-9) then
     391
     392              WBARI(L,nw,ng) = btemp(L,nw)  / DTAUI(L,NW,NG)
     393
     394           else
     395              WBARI(L,nw,ng) = 0.0D0
     396              DTAUI(L,NW,NG) = 1.0D-9
     397           endif
     398
     399           cosbi(L,NW,NG) = cosbi(L,NW,L_NGAUSS)
     400           
    357401        END IF
    358402

  • trunk/LMDZ.GENERIC/libf/phystd/optcv.F90

    r1397 r1715  
    2424  !     
    2525  !     THIS SUBROUTINE SETS THE OPTICAL CONSTANTS IN THE VISUAL 
    26   !     IT CALCUALTES FOR EACH LAYER, FOR EACH SPECRAL INTERVAL IN THE VISUAL
     26  !     IT CALCULATES FOR EACH LAYER, FOR EACH SPECTRAL INTERVAL IN THE VISUAL
    2727  !     LAYER: WBAR, DTAU, COSBAR
    2828  !     LEVEL: TAU
    … …  
    3939
    4040  real*8 DTAUV(L_NLAYRAD,L_NSPECTV,L_NGAUSS)
    41   real*8 DTAUKV(L_LEVELS+1,L_NSPECTV,L_NGAUSS)
     41  real*8 DTAUKV(L_LEVELS,L_NSPECTV,L_NGAUSS)
    4242  real*8 TAUV(L_NLEVRAD,L_NSPECTV,L_NGAUSS)
    4343  real*8 TAUCUMV(L_LEVELS,L_NSPECTV,L_NGAUSS)
    … …  
    4848
    4949  ! for aerosols
    50   real*8  QXVAER(L_LEVELS+1,L_NSPECTV,NAERKIND)
    51   real*8  QSVAER(L_LEVELS+1,L_NSPECTV,NAERKIND)
    52   real*8  GVAER(L_LEVELS+1,L_NSPECTV,NAERKIND)
    53   real*8  TAUAERO(L_LEVELS+1,NAERKIND)
    54   real*8  TAUAEROLK(L_LEVELS+1,L_NSPECTV,NAERKIND)
     50  real*8  QXVAER(L_LEVELS,L_NSPECTV,NAERKIND)
     51  real*8  QSVAER(L_LEVELS,L_NSPECTV,NAERKIND)
     52  real*8  GVAER(L_LEVELS,L_NSPECTV,NAERKIND)
     53  real*8  TAUAERO(L_LEVELS,NAERKIND)
     54  real*8  TAUAEROLK(L_LEVELS,L_NSPECTV,NAERKIND)
    5555  real*8  TAEROS(L_LEVELS,L_NSPECTV,NAERKIND)
    5656
    … …  
    6060  real*8  TAURAY(L_NSPECTV)
    6161  real*8  TRAY(L_LEVELS,L_NSPECTV)
    62   real*8  TRAYAER
    6362  real*8  DPR(L_LEVELS), U(L_LEVELS)
    6463  real*8  LCOEF(4), LKCOEF(L_LEVELS,4)
    … …  
    6665  real*8 taugsurf(L_NSPECTV,L_NGAUSS-1)
    6766  real*8 DCONT,DAERO
     67  real*8 DRAYAER
    6868  double precision wn_cont, p_cont, p_air, T_cont, dtemp, dtempc
    6969  double precision p_cross
    … …  
    123123  end do                    ! levels
    124124
    125 
     125  ! Spectral dependance of aerosol absorption
    126126  do iaer=1,naerkind
    127127     do NW=1,L_NSPECTV
    … …  
    131131     end do
    132132  end do
     133 
     134  ! Rayleigh scattering
    133135  do NW=1,L_NSPECTV
    134136     do K=2,L_LEVELS
    … …  
    142144     do NW=1,L_NSPECTV
    143145
    144         TRAYAER = TRAY(K,NW)
    145         !     TRAYAER is Tau RAYleigh scattering, plus AERosol opacity
     146        DRAYAER = TRAY(K,NW)
     147        !     DRAYAER is Tau RAYleigh scattering, plus AERosol opacity
    146148        do iaer=1,naerkind
    147            TRAYAER = TRAYAER + TAEROS(K,NW,IAER)
     149           DRAYAER = DRAYAER + TAEROS(K,NW,IAER)
    148150        end do
    149151
    … …  
    257259           TAUGSURF(NW,NG) = TAUGSURF(NW,NG) + TAUGAS + DCONT
    258260           DTAUKV(K,nw,ng) = TAUGAS &
    259                              + TRAYAER & ! TRAYAER includes all scattering contributions
     261                             + DRAYAER & ! DRAYAER includes all scattering contributions
    260262                             + DCONT ! For parameterized continuum aborption
    261263
    … …  
    266268
    267269        NG              = L_NGAUSS
    268         DTAUKV(K,nw,ng) = TRAY(K,NW) + DCONT ! For parameterized continuum absorption
    269 
    270         do iaer=1,naerkind
    271            DTAUKV(K,nw,ng) = DTAUKV(K,nw,ng) +  TAEROS(K,NW,IAER)
    272         end do ! a bug was here!
     270        DTAUKV(K,nw,ng) = DRAYAER + DCONT ! Scattering + parameterized continuum absorption
     271
     272!        do iaer=1,naerkind
     273!           DTAUKV(K,nw,ng) = DTAUKV(K,nw,ng) +  TAEROS(K,NW,IAER)
     274!        end do ! a bug was here!
    273275
    274276     end do
    … …  
    282284  do iaer=1,naerkind
    283285    DO NW=1,L_NSPECTV
    284       DO K=2,L_LEVELS   ! AS: shouldn't this be L_LEVELS+1 ? (see optci)
    285            TAUAEROLK(K,NW,IAER) = TAUAERO(K,IAER) * QSVAER(K,NW,IAER)
     286      DO K=2,L_LEVELS
     287           TAUAEROLK(K,NW,IAER) = TAUAERO(K,IAER) * QSVAER(K,NW,IAER) ! effect of scattering albedo
    286288      ENDDO
    287289    ENDDO
    … …  
    293295        atemp(L,NW) = SUM(GVAER(K,NW,1:naerkind) * TAUAEROLK(K,NW,1:naerkind))+SUM(GVAER(K+1,NW,1:naerkind) * TAUAEROLK(K+1,NW,1:naerkind))
    294296        btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind)) + SUM(TAUAEROLK(K+1,NW,1:naerkind))
    295         ctemp(L,NW) = btemp(L,NW) + 0.9999*(TRAY(K,NW) + TRAY(K+1,NW))
     297        ctemp(L,NW) = btemp(L,NW) + 0.9999*(TRAY(K,NW) + TRAY(K+1,NW))  ! JVO 2017 : does this 0.999 is really meaningful ?
    296298        btemp(L,NW) = btemp(L,NW) + TRAY(K,NW) + TRAY(K+1,NW)
    297299        COSBV(L,NW,1:L_NGAUSS) = atemp(L,NW)/btemp(L,NW)
    298300     END DO ! L vertical loop
    299301     
    300      !last level
    301      L              = L_NLAYRAD
    302      K              = 2*L+1
    303      atemp(L,NW)    = SUM(GVAER(K,NW,1:naerkind) * TAUAEROLK(K,NW,1:naerkind))
     302     ! Last level
     303     L           = L_NLAYRAD
     304     K           = 2*L+1
     305     atemp(L,NW) = SUM(GVAER(K,NW,1:naerkind) * TAUAEROLK(K,NW,1:naerkind))
    304306     btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind))
    305      ctemp(L,NW) = btemp(L,NW) + 0.9999*TRAY(K,NW)
     307     ctemp(L,NW) = btemp(L,NW) + 0.9999*TRAY(K,NW) ! JVO 2017 : does this 0.999 is really meaningful ?
    306308     btemp(L,NW) = btemp(L,NW) + TRAY(K,NW)
    307309     COSBV(L,NW,1:L_NGAUSS) = atemp(L,NW)/btemp(L,NW)
    … …  
    320322      END DO ! L vertical loop
    321323
    322         !     No vertical averaging on bottom layer
     324        ! Last level
    323325
    324326        L              = L_NLAYRAD

  • trunk/LMDZ.GENERIC/libf/phystd/radcommon_h.F90

    r1529 r1715  
    126126
    127127      real*8,save :: PTOP
    128       real*8,save,allocatable :: TAUREF(:)
    129128
    130129      real*8,parameter :: UBARI = 0.5D0
    … …  
    132131      real*8,save :: gweight(L_NGAUSS)
    133132!$OMP THREADPRIVATE(QREFvis,QREFir,omegaREFvis,omegaREFir,&     ! gweight read by master in sugas_corrk
    134                 !$OMP tstellar,planckir,PTOP,TAUREF)
     133                !$OMP tstellar,planckir,PTOP)
    135134
    136135!     If the gas optical depth (top to the surface) is less than
    … …  
    156155!$OMP THREADPRIVATE(glat,eclipse)
    157156
    158 contains
    159 
    160       subroutine ini_radcommon_h
    161       use radinc_h, only: L_LEVELS
    162       implicit none
    163      
    164         allocate(TAUREF(L_LEVELS+1))
    165      
    166       end subroutine ini_radcommon_h
    167 
    168157end module radcommon_h
Note: See TracChangeset for help on using the changeset viewer.