Ignore:
Timestamp:
Jul 24, 2023, 3:09:24 AM (15 months ago)
Author:
emillour
Message:

Generic GCM:

Chemistry: correction in photolysis online with level instead of layer

YJ

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.GENERIC/libf/aeronostd/photolysis_online.F

    r2553 r3009  
    22
    33      subroutine photolysis_online(nlayer, alt, press, temp,
    4      $           zmmean, rm,
     4     $           mmean, rm,
    55     $           tau, sza, dist_sol, v_phot, e_phot, ig, ngrid,
    66     $           nreact)
     
    3131      integer, intent(in) :: ngrid  ! number of atmospheric columns
    3232
    33       real,    intent(in), dimension(nlayer) :: press, temp, zmmean ! pressure (hpa)/temperature (k)/mean molecular mass (g.mol-1)
     33      real,    intent(in), dimension(nlayer) :: press, temp, mmean ! pressure (hpa)/temperature (k)/mean molecular mass (g.mol-1)
    3434      real,    intent(in), dimension(nlayer) :: alt                 ! altitude (km)
    3535      real,    intent(in), dimension(nlayer,nesp) :: rm             ! mixing ratios
     
    5252!     atmosphere
    5353
    54       real, dimension(nw)             :: albedo_chim                ! Surface albedo calculated on chemistry wavelenght grid
    55       real, dimension(nlayer)         :: colinc                     ! air column increment (molecule.cm-2)
    56       real, dimension(nlayer)         :: airlev                     ! air density at each specified altitude (molec/cc)
    57       real, dimension(nlayer)         :: edir, edn, eup             ! normalised irradiances
    58       real, dimension(nlayer)         :: fdir, fdn, fup             ! normalised actinic fluxes
    59       real, dimension(nlayer)         :: saflux                     ! total actinic flux
    60       real, dimension(nlayer,nw)      :: dtrl                       ! rayleigh optical depth
    61       real, dimension(nlayer,nw)      :: dtaer                      ! aerosol optical depth
    62       real, dimension(nlayer,nw)      :: omaer                      ! aerosol single scattering albedo
    63       real, dimension(nlayer,nw)      :: gaer                       ! aerosol asymmetry parameter
    64       real, dimension(nlayer,nw)      :: dtcld                      ! cloud optical depth
    65       real, dimension(nlayer,nw)      :: omcld                      ! cloud single scattering albedo
    66       real, dimension(nlayer,nw)      :: gcld                       ! cloud asymmetry parameter
    67       real, dimension(nlayer,nw)      :: dagas                      ! total gas optical depth
    68       real, dimension(nlayer,nw,nabs) :: dtgas                      ! optical depth for each gas
    69 
    70       integer, dimension(0:nlayer)        :: nid
    71       real,    dimension(0:nlayer,nlayer) :: dsdh
     54      real, dimension(nw)               :: albedo_chim              ! Surface albedo calculated on chemistry wavelenght grid
     55      real, dimension(nlayer+1)         :: colinc                   ! air column increment (molecule.cm-2)
     56      real, dimension(nlayer+1)         :: airlev                   ! air density at each specified altitude (molec/cc)
     57      real, dimension(nlayer+1)         :: edir, edn, eup           ! normalised irradiances
     58      real, dimension(nlayer+1)         :: fdir, fdn, fup           ! normalised actinic fluxes
     59      real, dimension(nlayer+1)         :: saflux                   ! total actinic flux
     60      real, dimension(nlayer+1,nw)      :: dtrl                     ! rayleigh optical depth
     61      real, dimension(nlayer+1,nw)      :: dtaer                    ! aerosol optical depth
     62      real, dimension(nlayer+1,nw)      :: omaer                    ! aerosol single scattering albedo
     63      real, dimension(nlayer+1,nw)      :: gaer                     ! aerosol asymmetry parameter
     64      real, dimension(nlayer+1,nw)      :: dtcld                    ! cloud optical depth
     65      real, dimension(nlayer+1,nw)      :: omcld                    ! cloud single scattering albedo
     66      real, dimension(nlayer+1,nw)      :: gcld                     ! cloud asymmetry parameter
     67      real, dimension(nlayer+1,nw)      :: dagas                    ! total gas optical depth
     68      real, dimension(nlayer+1,nw,nabs) :: dtgas                    ! optical depth for each gas
     69      real, dimension(nlayer+1)         :: zpress                   ! pressure (hpa)/altitude (km)/temperature (k)/mean molecular mass (g.mol-1)
     70      real, dimension(nlayer+1)         :: zalt                     ! pressure (hpa)/altitude (km)/temperature (k)/mean molecular mass (g.mol-1)
     71      real, dimension(nlayer+1)         :: ztemp                    ! pressure (hpa)/altitude (km)/temperature (k)/mean molecular mass (g.mol-1)
     72      real, dimension(nlayer+1)         :: zmmean                   ! pressure (hpa)/altitude (km)/temperature (k)/mean molecular mass (g.mol-1)
     73
     74      integer, dimension(0:nlayer+1)          :: nid
     75      real,    dimension(0:nlayer+1,nlayer+1) :: dsdh
    7276
    7377      integer :: i, ilay, iw, ialt, iphot, ispe, ij, igas, ireact
     78      integer :: ilev, nlev
    7479      real    :: deltaj
    7580
     81!==== define working vertical grid for the uv radiative code
     82
     83      nlev = nlayer + 1
     84
     85      do ilev = 1,nlev-1
     86         zpress(ilev) = press(ilev)
     87         zalt(ilev)   = alt(ilev)
     88         ztemp(ilev)  = temp(ilev)
     89         zmmean(ilev) = mmean(ilev)
     90      end do
     91
     92      zpress(nlev) = 0.                  ! top of the atmosphere
     93      zalt(nlev)   = zalt(nlev-1) + (zalt(nlev-1) - zalt(nlev-2))
     94      ztemp(nlev)  = ztemp(nlev-1)
     95      zmmean(nlev) = zmmean(nlev-1)
     96
    7697!==== air column increments and rayleigh optical depth
    7798
    78       call setair(nlayer, nw, wl, wc, press, temp, zmmean, colinc, dtrl,
     99      call setair(nlev, nw, wl, wc, zpress, ztemp, zmmean, colinc, dtrl,
    79100     $            airlev)
    80101
     
    89110      igas   = 0
    90111      ireact = 0
     112      dtgas(:,:,:) = 0.
    91113
    92114      do while(iphot<nb_phot_hv_max)
     
    166188!==== set aerosol properties and optical depth
    167189
    168       call setaer(nlayer,alt,tau,nw,dtaer,omaer,gaer)
     190      call setaer(nlev,zalt,tau,nw,dtaer,omaer,gaer)
    169191
    170192!==== set cloud properties and optical depth
    171193
    172       call setcld(nlayer,nw,dtcld,omcld,gcld)
     194      call setcld(nlev,nw,dtcld,omcld,gcld)
    173195
    174196!==== slant path lengths in spherical geometry
    175197
    176       call sphers(nlayer,alt,sza,dsdh,nid)
     198      call sphers(nlev,zalt,sza,dsdh,nid)
    177199
    178200!==== solar flux at mars
     
    198220!     dir = direct beam, dn = down-welling diffuse, up = up-welling diffuse
    199221
    200          call rtlink(nlayer, nw, iw, albedo_chim(iw),
     222         call rtlink(nlev, nw, iw, albedo_chim(iw),
    201223     $               sza, dsdh, nid, dtrl,
    202224     $               dagas, dtcld, omcld, gcld, dtaer, omaer, gaer,
     
    276298         colinc(ilev) = avocado*0.1*dp/(zmmean(ilev)*g)
    277299      end do
    278       colinc(nlev) = avocado*0.1*press(nlev)*100./(zmmean(nlev)*g)
     300      colinc(nlev) = 0.0
    279301
    280302      do iw = 1, nw - 1
     
    308330!==============================================================================
    309331
    310       subroutine setaer(nlayer,alt,tau,nw,dtaer,omaer,gaer)
     332      subroutine setaer(nlev,zalt,tau,nw,dtaer,omaer,gaer)
    311333
    312334!-----------------------------------------------------------------------------*
     
    320342!     input
    321343
    322       integer :: nlayer                              ! number of vertical layers
    323       integer :: nw                                  ! number of wavelength grid points
    324       real, dimension(nlayer) :: alt                 ! altitude (km)
    325       real :: tau                                    ! integrated aerosol optical depth at the surface
     344      integer :: nlev                              ! number of vertical layers
     345      integer :: nw                                ! number of wavelength grid points
     346      real, dimension(nlev) :: zalt                ! altitude (km)
     347      real :: tau                                  ! integrated aerosol optical depth at the surface
    326348
    327349!     output
    328350
    329       real, dimension(nlayer,nw) :: dtaer            ! aerosol optical depth
    330       real, dimension(nlayer,nw) :: omaer            ! aerosol single scattering albedo
    331       real, dimension(nlayer,nw) :: gaer             ! aerosol asymmetry parameter
     351      real, dimension(nlev,nw) :: dtaer            ! aerosol optical depth
     352      real, dimension(nlev,nw) :: omaer            ! aerosol single scattering albedo
     353      real, dimension(nlev,nw) :: gaer             ! aerosol asymmetry parameter
    332354
    333355!     local
    334356
    335       integer :: ilay, iw
    336       real, dimension(nlayer) :: aer                                 ! dust extinction
     357      integer :: ilev, iw
     358      real, dimension(nlev) :: aer                 ! dust extinction
    337359      real :: omega, g, scaleh, gamma
    338360      real :: dz, tautot, q0
     
    350372
    351373      tautot = 0.
    352       do ilay = 1, nlayer-1
    353          dz = alt(ilay+1) - alt(ilay)
    354          tautot = tautot + exp(gamma*(1. - exp(alt(ilay)/scaleh)))*dz
     374      do ilev = 1, nlev-1
     375         dz = zalt(ilev+1) - zalt(ilev)
     376         tautot = tautot + exp(gamma*(1. - exp(zalt(ilev)/scaleh)))*dz
    355377      end do
    356378     
    357379      q0 = tau/tautot
    358       do ilay = 1, nlayer-1
    359          dz = alt(ilay+1) - alt(ilay)
    360          dtaer(ilay,:) = q0*exp(gamma*(1. - exp(alt(ilay)/scaleh)))*dz
    361          omaer(ilay,:) = omega
    362          gaer(ilay,:)  = g
     380      do ilev = 1, nlev-1
     381         dz = zalt(ilev+1) - zalt(ilev)
     382         dtaer(ilev,:) = q0*exp(gamma*(1. - exp(zalt(ilev)/scaleh)))*dz
     383         omaer(ilev,:) = omega
     384         gaer(ilev,:)  = g
    363385      end do
    364386
     
    367389!==============================================================================
    368390
    369       subroutine setcld(nlayer,nw,dtcld,omcld,gcld)
     391      subroutine setcld(nlev,nw,dtcld,omcld,gcld)
    370392
    371393!-----------------------------------------------------------------------------*
     
    379401!     input
    380402
    381       integer :: nlayer                              ! number of vertical layers
    382       integer :: nw                                  ! number of wavelength grid points
     403      integer :: nlev                              ! number of vertical layers
     404      integer :: nw                                ! number of wavelength grid points
    383405
    384406!     output
    385407
    386       real, dimension(nlayer,nw) :: dtcld            ! cloud optical depth
    387       real, dimension(nlayer,nw) :: omcld            ! cloud single scattering albedo
    388       real, dimension(nlayer,nw) :: gcld             ! cloud asymmetry parameter
     408      real, dimension(nlev,nw) :: dtcld            ! cloud optical depth
     409      real, dimension(nlev,nw) :: omcld            ! cloud single scattering albedo
     410      real, dimension(nlev,nw) :: gcld             ! cloud asymmetry parameter
    389411
    390412!     local
    391413
    392       integer :: ilay, iw
     414      integer :: ilev, iw
    393415
    394416!     dtcld : optical depth
     
    396418!     gcld  : asymmetry factor
    397419
    398       do ilay = 1, nlayer - 1
     420      do ilev = 1, nlev - 1
    399421         do iw = 1, nw - 1
    400             dtcld(ilay,iw) = 0.       ! no clouds for the moment
    401             omcld(ilay,iw) = 0.99
    402             gcld(ilay,iw)  = 0.85
     422            dtcld(ilev,iw) = 0.       ! no clouds for the moment
     423            omcld(ilev,iw) = 0.99
     424            gcld(ilev,iw)  = 0.85
    403425         end do
    404426      end do
     
    11501172
    11511173      use chimiedata_h,  only: albedo_snow_chim, albedo_co2_ice_chim
    1152       use slab_ice_h,    only: h_alb_ice, alb_ice_min, alb_ice_max
     1174      use slab_ice_h,    only: alb_ocean, h_alb_ice,
     1175     &                         alb_ice_min, alb_ice_max
    11531176      use tracer_h,      only: igcm_h2o_ice, igcm_co2_ice
    1154       use callkeys_mod,  only: ok_slab_ocean, co2cond, alb_ocean
     1177      use callkeys_mod,  only: ok_slab_ocean, co2cond
    11551178      use phys_state_var_mod, only: albedo_bareground,
    11561179     &                              rnat, qsurf, sea_ice,
Note: See TracChangeset for help on using the changeset viewer.