Changeset 3959 for trunk/LMDZ.PLUTO/libf


Ignore:
Timestamp:
Nov 12, 2025, 6:12:16 PM (4 weeks ago)
Author:
debatzbr
Message:

Pluto PCM: Take clouds into account in radiative transfer.
BBT

Location:
trunk/LMDZ.PLUTO/libf/phypluto
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.PLUTO/libf/phypluto/aeropacity.F90

    r3802 r3959  
    1313       use comcstfi_mod, only: g, pi, mugaz, avocado
    1414       use geometry_mod, only: latitude
    15        use callkeys_mod, only: kastprof, callmufi
     15       use callkeys_mod, only: kastprof, callmufi, callmuclouds
    1616       use mp2m_diagnostics
    1717       implicit none
     
    3232!     Input
    3333!     -----
    34 !     ngrid             Number of horizontal gridpoints
    35 !     nlayer            Number of layers
    36 !     nq                Number of tracers
    37 !     pplev             Pressure (Pa) at each layer boundary
    38 !     pq                Aerosol mixing ratio
    39 !     reffrad(ngrid,nlayer,naerkind)         Aerosol effective radius
    40 !     QREFvis3d(ngrid,nlayer,naerkind) \ 3d extinction coefficients
     34!     ngrid                            Number of horizontal gridpoints
     35!     nlayer                           Number of layers
     36!     nq                               Number of tracers
     37!     pplev                            Pressure (Pa) at each layer boundary
     38!     pq                               Aerosol mixing ratio
     39!     reffrad(ngrid,nlayer,naerkind)   Aerosol effective radius
     40!     QREFvis3d(ngrid,nlayer,naerkind) \ 3D extinction coefficients
    4141!     QREFir3d(ngrid,nlayer,naerkind)  / at reference wavelengths
    4242!
     
    4848!=======================================================================
    4949
    50       INTEGER,INTENT(IN) :: ngrid  ! number of atmospheric columns
    51       INTEGER,INTENT(IN) :: nlayer ! number of atmospheric layers
    52       INTEGER,INTENT(IN) :: nq     ! number of tracers
    53       REAL,INTENT(IN) :: pplay(ngrid,nlayer) ! mid-layer pressure (Pa)
    54       REAL,INTENT(IN) :: pplev(ngrid,nlayer+1) ! inter-layer pressure (Pa)
    55       REAL,INTENT(IN) :: zzlev(ngrid,nlayer)   ! Altitude at the layer boundaries.
    56       REAL,INTENT(IN) :: pq(ngrid,nlayer,nq) ! tracers (.../kg_of_air)
    57       REAL,INTENT(IN) :: pt(ngrid,nlayer) ! mid-layer temperature (K)
    58       REAL,INTENT(OUT) :: dtau_aer(ngrid,nlayer,naerkind) ! aerosol optical depth
    59       REAL,INTENT(IN) :: reffrad(ngrid,nlayer,naerkind)   ! aerosol effective radius
    60       REAL,INTENT(IN) :: nueffrad(ngrid,nlayer,naerkind)  ! aerosol effective variance
    61       REAL,INTENT(IN) :: QREFvis3d(ngrid,nlayer,naerkind) ! extinction coefficient in the visible
    62       REAL,INTENT(IN) :: QREFir3d(ngrid,nlayer,naerkind)
    63       REAL,INTENT(OUT):: tau_col(ngrid) !column integrated visible optical depth
     50      INTEGER,INTENT(IN) :: ngrid                            ! Number of atmospheric columns
     51      INTEGER,INTENT(IN) :: nlayer                           ! Number of atmospheric layers
     52      INTEGER,INTENT(IN) :: nq                               ! Number of tracers
     53      REAL,INTENT(IN)    :: pplay(ngrid,nlayer)              ! Mid-layer pressure (Pa)
     54      REAL,INTENT(IN)    :: pplev(ngrid,nlayer+1)            ! Inter-layer pressure (Pa)
     55      REAL,INTENT(IN)    :: zzlev(ngrid,nlayer)              ! Altitude at the layer boundaries.
     56      REAL,INTENT(IN)    :: pq(ngrid,nlayer,nq)              ! Tracers (.../kg_of_air)
     57      REAL,INTENT(IN)    :: pt(ngrid,nlayer)                 ! Mid-layer temperature (K)
     58      REAL,INTENT(OUT)   :: dtau_aer(ngrid,nlayer,naerkind)  ! Aerosol optical depth
     59      REAL,INTENT(IN)    :: reffrad(ngrid,nlayer,naerkind)   ! Aerosol effective radius
     60      REAL,INTENT(IN)    :: nueffrad(ngrid,nlayer,naerkind)  ! Aerosol effective variance
     61      REAL,INTENT(IN)    :: QREFvis3d(ngrid,nlayer,naerkind) ! Extinction coefficient in the visible
     62      REAL,INTENT(IN)    :: QREFir3d(ngrid,nlayer,naerkind)  ! Extinction coefficient in the infrared
     63      REAL,INTENT(OUT)   :: tau_col(ngrid)                   ! Column integrated visible optical depth
    6464
    6565      real aerosol0, obs_tau_col_aurora, pm
     
    7474      real m0as(ngrid,nlayer)
    7575      real m0af(ngrid,nlayer)
     76      real m0ccn(ngrid,nlayer)
    7677
    7778      INTEGER l,ig,iq,iaer,ia
     
    135136            endwhere
    136137           
    137             ! write(*,*) 'dtau_as  :', MINVAL(dtau_aer(:,:,1)), '-', MAXVAL(dtau_aer(:,:,1))
    138             ! write(*,*) 'dtau_af  :', MINVAL(dtau_aer(:,:,2)), '-', MAXVAL(dtau_aer(:,:,2))
     138            ! Cloud drops
     139            if (callmuclouds) then
     140               m0ccn(:,:) = pq(:,:,micro_indx(5)) * (pplev(:,1:nlayer) - pplev(:,2:nlayer+1)) / g
     141               sig = 0.2
     142               where ((m0ccn(:,:) >= 1e-8).and.(mp2m_rc_cld(:,:) >= 5e-9))
     143                  dtau_aer(:,:,3) = m0ccn(:,:) * QREFvis3d(:,:,3) * pi * mp2m_rc_cld(:,:)**2 * exp(2*sig**2)
     144               elsewhere
     145                  dtau_aer(:,:,3) = 0d0
     146               endwhere
     147            endif ! end callmuclouds
    139148         
    140149         else
  • trunk/LMDZ.PLUTO/libf/phypluto/callcorrk.F90

    r3955 r3959  
    3737                              methane,carbox,cooling,nlte,strobel,&
    3838                              ch4fix,vmrch4_proffix,vmrch4fix,&
    39                               callmufi,triton
     39                              callmufi,callmuclouds,triton
    4040      use optcv_mod, only: optcv
    4141      use optci_mod, only: optci
     
    492492            endwhere
    493493            nueffrad(:,:,2) = exp(sig**2) - 1
     494            ! Cloud drops
     495            if (callmuclouds) then
     496               sig = 0.2
     497               where (mp2m_rc_cld(:,:) >= 5e-9)
     498                  reffrad(:,:,3) = mp2m_rc_cld(:,:) * exp(5.*sig**2 / 2.)
     499               elsewhere
     500                  reffrad(:,:,3) = 0d0
     501               endwhere
     502               nueffrad(:,:,3) = exp(sig**2) - 1
     503            endif ! end callmuclouds
    494504
    495505         else
     
    593603!-----------------------------------------------------------------------
    594604
    595 
    596       ! AF24: for now only consider one aerosol (=haze)
    597605      if (optichaze) then
    598606        do iaer=1,naerkind
     
    676684        end do ! naerkind
    677685
    678         ! Test / Correct for freaky s. s. albedo values.
     686        ! Test / Correct for freaky single scattering albedo values.
    679687        do iaer=1,naerkind
    680688            do k=1,L_LEVELS
  • trunk/LMDZ.PLUTO/libf/phypluto/callkeys_mod.F90

    r3949 r3959  
    217217!! Variables for aerosol absorption
    218218      real,save :: Fabs_aers_VI, Fabs_aerf_VI, Fabs_aers_IR, Fabs_aerf_IR
    219 !$OMP THREADPRIVATE(Fabs_aers_VI, Fabs_aerf_VI, Fabs_aers_IR, Fabs_aerf_IR)
     219!$OMP THREADPRIVATE(Fabs_aers_VI,Fabs_aerf_VI,Fabs_aers_IR,Fabs_aerf_IR)
     220!! Variables for cloud drop absorption
     221      real,save :: Fabs_cldd_VI, Fabs_cldd_IR
     222!$OMP THREADPRIVATE(Fabs_cldd_VI,Fabs_cldd_IR)
    220223
    221224      integer,save :: iddist
  • trunk/LMDZ.PLUTO/libf/phypluto/datafile_mod.F90

    r3949 r3959  
    2222      character(len=300),save :: aersprop_file
    2323      character(len=300),save :: aerfprop_file
     24      character(len=300),save :: clddprop_file
    2425
    2526      ! surfdir stores planetary topography, albedo, etc. (surface.nc files)
  • trunk/LMDZ.PLUTO/libf/phypluto/inifis_mod.F90

    r3951 r3959  
    1414  use radii_mod, only: radfixed
    1515  use datafile_mod, only: datadir,hazeprop_file,hazerad_file,hazemmr_file,hazedens_file, &
    16                           config_mufi, mugasflux_file, aersprop_file, aerfprop_file
     16                          config_mufi, mugasflux_file, aersprop_file, aerfprop_file, clddprop_file
    1717  use comdiurn_h, only: sinlat, coslat, sinlon, coslon
    1818  use comgeomfi_h, only: totarea, totarea_planet
     
    801801     if (is_master) write(*,*) trim(rname)//" aerfprop_file = ",trim(aerfprop_file)
    802802
     803     if (is_master) write(*,*) "Cloud drop optical properties datafile"
     804     clddprop_file="optprop_rannou_r2-200nm_nu003.dat"  ! default file
     805     call getin_p("clddprop_file",clddprop_file)
     806     if (is_master) write(*,*) trim(rname)//" clddprop_file = ",trim(clddprop_file)
     807
    803808     if (is_master) write(*,*) "Use haze production from CH4 photolysis or production rate?"
    804809     call_haze_prod_pCH4=.false. ! default value
     
    879884     call getin_p("Fabs_aerf_IR",Fabs_aerf_IR)
    880885     if (is_master) write(*,*)" Fabs_aerf_IR = ",Fabs_aerf_IR
     886
     887     if (is_master) write(*,*) "Cloud drop absorption correction in VI?"
     888     Fabs_cldd_VI=1. ! default value
     889     call getin_p("Fabs_cldd_VI",Fabs_cldd_VI)
     890     if (is_master) write(*,*)" Fabs_cldd_VI = ",Fabs_cldd_VI
     891
     892     if (is_master) write(*,*) "Cloud drop absorption correction in IR?"
     893     Fabs_cldd_IR=1. ! default value
     894     call getin_p("Fabs_cldd_IR",Fabs_cldd_IR)
     895     if (is_master) write(*,*)" Fabs_cldd_IR = ",Fabs_cldd_IR
    881896
    882897     ! Pluto haze model
     
    14541469      call abort_physic(rname, 'if microphysics is on, naerkind must be > 1!', 1)
    14551470     endif
    1456      ! if ((callmufi).and.(callmuclouds).and..not.(naerkind.gt.2)) then
    1457      !  call abort_physic(rname, 'if microphysical clouds are on, naerkind must be > 2!', 1)
    1458      ! endif
     1471     if ((callmufi).and..not.(callmuclouds).and.(naerkind.gt.2)) then
     1472      call abort_physic(rname, 'Warning: here microphysical clouds are on, naerkind must be = 2!', 1)
     1473     endif
     1474     if ((callmufi).and.(callmuclouds).and..not.(naerkind.gt.2)) then
     1475      call abort_physic(rname, 'if microphysical clouds are on, naerkind must be > 2!', 1)
     1476     endif
    14591477     if (.not.(callmufi.or.haze).and.(optichaze)) then
    14601478      call abort_physic(rname, 'if microphysics and haze are off, optichaze must be deactivated!', 1)
  • trunk/LMDZ.PLUTO/libf/phypluto/optci.F90

    r3889 r3959  
    1515                     igas_CH4, igas_N2
    1616  use comcstfi_mod, only: g, r, mugaz
    17   use callkeys_mod, only: kastprof,continuum,graybody,callmufi,Fabs_aers_IR,Fabs_aerf_IR
     17  use callkeys_mod, only: kastprof,continuum,graybody,callmufi,callmuclouds, &
     18                          Fabs_aers_IR,Fabs_aerf_IR,Fabs_cldd_IR
    1819  use recombin_corrk_mod, only: corrk_recombin, gasi_recomb
    1920  use tpindex_mod, only: tpindex
     
    143144  lkcoef(:,:)   = 0.0
    144145
    145   if(callmufi) then
     146  if (callmufi) then
    146147   Fabs_aer(1) = Fabs_aers_IR
    147148   Fabs_aer(2) = Fabs_aerf_IR
     149   if (callmuclouds) then
     150      Fabs_aer(3) = Fabs_cldd_IR
     151   endif
    148152  else
    149153   Fabs_aer(:) = 1.0
    150   endif
     154  endif ! end callmufi
    151155
    152156  do K=2,L_LEVELS
     
    186190        do K=2,L_LEVELS
    187191           TAEROS(K,NW,IAER) = TAUAERO(K,IAER) * QXIAER(K,NW,IAER)
    188            ! Aerosol absorption correction --> impact on opacity.
     192           ! Aerosol/Cloud absorption correction --> impact on opacity.
    189193           if (callmufi .and. (TAEROS(K,NW,IAER).gt.0.d0)) then
    190194            TAEROS(K,NW,IAER) = TAEROS(K,NW,IAER) * &
     
    343347        ! Aerosol absorption correction --> impact on single scattering albedo.
    344348        if (callmufi) then
     349           ! Spherical aerosols
    345350           if ((TAEROS(K,NW,1).gt.0.d0) .and. TAEROS(K+1,NW,1).gt.0.d0) then
    346351              btemp(L,NW) = TAUAEROLK(K,NW,1) / &
     
    353358              btemp(L,NW) = TAUAEROLK(K,NW,1) + TAUAEROLK(K+1,NW,1)
    354359           endif
     360           ! Fractal aerosols
    355361           if ((TAEROS(K,NW,2).gt.0.d0) .and. TAEROS(K+1,NW,2).gt.0.d0) then
    356362              btemp(L,NW) = btemp(L,NW) + &
     
    364370              btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,2) + TAUAEROLK(K+1,NW,2)
    365371           endif
     372           ! Cloud absorption correction --> impact on single scattering albedo.
     373           if (callmuclouds) then
     374              if ((TAEROS(K,NW,3).gt.0.d0) .and. TAEROS(K+1,NW,3).gt.0.d0) then
     375               btemp(L,NW) = btemp(L,NW) + &
     376                           (TAUAEROLK(K,NW,3) / &
     377                           (QSIAER(K,NW,3)/QXIAER(K,NW,3) + &
     378                              Fabs_aer(3)*(1.-QSIAER(K,NW,3)/QXIAER(K,NW,3)))) + &
     379                           (TAUAEROLK(K+1,NW,3) / &
     380                           (QSIAER(K+1,NW,3)/QXIAER(K+1,NW,3) + &
     381                              Fabs_aer(3)*(1.-QSIAER(K+1,NW,3)/QXIAER(K+1,NW,3))))
     382              else
     383               btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,3) + TAUAEROLK(K+1,NW,3)
     384              endif
     385           endif ! end callmuclouds
    366386        else
    367387           btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind)) + SUM(TAUAEROLK(K+1,NW,1:naerkind))
     
    375395     ! Aerosol absorption correction --> impact on single scattering albedo.
    376396     if (callmufi) then
     397      ! Spherical aerosols
    377398      if (TAEROS(K,NW,1).gt.0.d0) then
    378399         btemp(L,NW) = TAUAEROLK(K,NW,1) / &
     
    382403         btemp(L,NW) = TAUAEROLK(K,NW,1)
    383404      endif
     405      ! Fractal aerosols
    384406      if (TAEROS(K,NW,2).gt.0.d0) then
    385407         btemp(L,NW) = btemp(L,NW) + &
     
    390412         btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,2)
    391413      endif
     414      ! Cloud absorption correction --> impact on single scattering albedo.
     415      if (callmuclouds) then
     416         if (TAEROS(K,NW,3).gt.0.d0) then
     417            btemp(L,NW) = btemp(L,NW) + &
     418                          (TAUAEROLK(K,NW,3) / &
     419                          (QSIAER(K,NW,3)/QXIAER(K,NW,3) + &
     420                           Fabs_aer(3)*(1.-QSIAER(K,NW,3)/QXIAER(K,NW,3))))
     421         else
     422            btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,3)
     423         endif
     424      endif ! end callmuclouds
    392425     else
    393426      btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind))
  • trunk/LMDZ.PLUTO/libf/phypluto/optcv.F90

    r3947 r3959  
    1515                     igas_CH4, igas_N2
    1616  use comcstfi_mod, only: g, r, mugaz
    17   use callkeys_mod, only: kastprof,continuum,graybody,callgasvis,callmufi,Fabs_aers_VI,Fabs_aerf_VI
     17  use callkeys_mod, only: kastprof,continuum,graybody,callgasvis,callmufi,callmuclouds, &
     18                          Fabs_aers_VI,Fabs_aerf_VI,Fabs_cldd_VI
    1819  use recombin_corrk_mod, only: corrk_recombin, gasv_recomb
    1920  use tpindex_mod, only: tpindex
     
    146147  wbarv_aer(:,:,:) = 1.0
    147148
    148   if(callmufi) then
     149  if (callmufi) then
    149150   Fabs_aer(1) = Fabs_aers_VI
    150151   Fabs_aer(2) = Fabs_aerf_VI
     152   if (callmuclouds) then
     153      Fabs_aer(3) = Fabs_cldd_VI
     154   endif
    151155  else
    152156   Fabs_aer(:) = 1.0
     
    185189        do K=5,L_LEVELS
    186190           TAEROS(K,NW,IAER) = TAUAERO(K,IAER) * QXVAER(K,NW,IAER)
    187            ! Aerosol absorption correction --> impact on opacity.
     191           ! Aerosol/Cloud absorption correction --> impact on opacity.
    188192           if (callmufi .and. (TAEROS(K,NW,IAER).gt.0.d0)) then
    189193            TAEROS(K,NW,IAER) = TAEROS(K,NW,IAER) * &
     
    376380            btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,2) + TAUAEROLK(K+1,NW,2)
    377381         endif
     382
     383         ! Cloud absorption correction --> impact on single scattering albedo.
     384         if (callmuclouds) then
     385            if ((TAEROS(K,NW,3).gt.0.d0) .and. TAEROS(K+1,NW,3).gt.0.d0) then
     386               btemp(L,NW) = btemp(L,NW) + &
     387                           (TAUAEROLK(K,NW,3) / &
     388                           (QSVAER(K,NW,3)/QXVAER(K,NW,3) + &
     389                              Fabs_aer(3)*(1.-QSVAER(K,NW,3)/QXVAER(K,NW,3)))) + &
     390                           (TAUAEROLK(K+1,NW,3) / &
     391                           (QSVAER(K+1,NW,3)/QXVAER(K+1,NW,3) + &
     392                              Fabs_aer(3)*(1.-QSVAER(K+1,NW,3)/QXVAER(K+1,NW,3))))
     393            else
     394               btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,3) + TAUAEROLK(K+1,NW,3)
     395            endif
     396         endif ! end callmuclouds
    378397      else
    379398         btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind)) + SUM(TAUAEROLK(K+1,NW,1:naerkind))
     
    409428         btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,2)
    410429      endif
     430      ! Cloud absorption correction --> impact on single scattering albedo.
     431      if (callmuclouds) then
     432         if (TAEROS(K,NW,3).gt.0.d0) then
     433            btemp(L,NW) = btemp(L,NW) + &
     434                          (TAUAEROLK(K,NW,3) / &
     435                          (QSVAER(K,NW,3)/QXVAER(K,NW,3) + &
     436                           Fabs_aer(3)*(1.-QSVAER(K,NW,3)/QXVAER(K,NW,3))))
     437         else
     438            btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,3)
     439         endif
     440      endif ! end callmuclouds
    411441     else
    412442      btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind))
     
    421451    DO NW=1,L_NSPECTV
    422452     DO L=1,L_NLAYRAD-1
    423 
    424453        K              = 2*L+1
    425454        DTAUV(L,nw,ng) = DTAUKV(K,NW,NG) + DTAUKV(K+1,NW,NG)
    426455        WBARV(L,nw,ng) = ctemp(L,NW) / DTAUV(L,nw,ng)
    427 
    428456      END DO ! L vertical loop
    429457
    430         ! Last level
    431 
    432         L              = L_NLAYRAD
    433         K              = 2*L+1
    434         DTAUV(L,nw,ng) = DTAUKV(K,NW,NG)
    435 
    436         WBARV(L,NW,NG) = ctemp(L,NW) / DTAUV(L,NW,NG)
    437 
    438      END DO                 ! NW spectral loop
    439   END DO                    ! NG Gauss loop
     458      ! Last level
     459      !-----------
     460      L = L_NLAYRAD
     461      K = 2*L+1
     462           
     463      DTAUV(L,nw,ng) = DTAUKV(K,NW,NG)
     464      WBARV(L,NW,NG) = ctemp(L,NW) / DTAUV(L,NW,NG)
     465     END DO ! NW spectral loop
     466  END DO ! NG Gauss loop
    440467
    441468  ! Aerosols extinction optical depths
  • trunk/LMDZ.PLUTO/libf/phypluto/suaer_corrk.F90

    r3613 r3959  
    1111      use radcommon_h, only: blamv,blami,lamrefir,lamrefvis
    1212      use datafile_mod, only: datadir, aerdir, &
    13                               hazeprop_file, aersprop_file, aerfprop_file
     13                              hazeprop_file, aersprop_file, aerfprop_file, clddprop_file
    1414
    1515      ! outputs
     
    1818      use radcommon_h, only: qrefvis,qrefir,omegarefir !,omegarefvis
    1919      use aerosol_mod, only: iaero_haze
    20       use callkeys_mod, only: tplanet, callmufi
     20      use callkeys_mod, only: tplanet, callmufi, callmuclouds
    2121      use tracer_h, only: noms
    2222
     
    127127
    128128!--------------------------------------------------------------
    129       ! allocate file_id, as naerkind is a variable
     129      ! allocate file_id, as naerkind is a variable (VIS & IR)
    130130      allocate(file_id(naerkind,2))
    131 
    132131
    133132      if (callmufi) then
     
    137136            write(*,*)'Suaer fractal aerosols optical properties, using: ', &
    138137                       TRIM(aerfprop_file)
     138            if (callmuclouds) then
     139               write(*,*)'Suaer cloud drop optical properties, using: ', &
     140                          TRIM(clddprop_file)
     141            endif
    139142         endif
    140143         ! Visible
    141144         file_id(1,1)=TRIM(aersprop_file)
    142145         file_id(2,1)=TRIM(aerfprop_file)
     146         if (callmuclouds) then
     147            file_id(3,1)=TRIM(clddprop_file)
     148         endif
    143149         ! Infrared
    144150         file_id(1,2)=file_id(1,1)
    145151         file_id(2,2)=file_id(2,1)
     152         if (callmuclouds) then
     153            file_id(3,2)=file_id(3,1)
     154         endif
    146155
    147156         do iaer=1,naerkind
Note: See TracChangeset for help on using the changeset viewer.