LCOV - code coverage report
Current view: top level - physics/waccm - nlte_lw.F90 (source / functions) Hit Total Coverage
Test: coverage.info Lines: 95 145 65.5 %
Date: 2025-03-14 01:23:43 Functions: 3 5 60.0 %

          Line data    Source code
       1             : module nlte_lw
       2             : 
       3             : !
       4             : ! interface for calculation of non-LTE heating rates
       5             : !
       6             :   use shr_kind_mod,       only: r8 => shr_kind_r8
       7             :   use spmd_utils,         only: masterproc
       8             :   use ppgrid,             only: pcols, pver
       9             :   use pmgrid,             only: plev
      10             :   use rad_constituents,   only: rad_cnst_get_gas, rad_cnst_get_info
      11             : 
      12             :   use nlte_fomichev,      only: nlte_fomichev_init, nlte_fomichev_calc, nocooling, o3pcooling
      13             :   use nlte_aliarms,       only: nlte_aliarms_init, nlte_aliarms_calc
      14             : 
      15             :   use waccm_forcing,      only: waccm_forcing_init, waccm_forcing_adv,  get_cnst
      16             :   use cam_logfile,        only: iulog
      17             : 
      18             :   implicit none
      19             :   private
      20             :   save
      21             : 
      22             : ! Public interfaces
      23             :   public &
      24             :        nlte_register,      &
      25             :        nlte_init,          &
      26             :        nlte_timestep_init, &
      27             :        nlte_tend
      28             : 
      29             : ! Private module data
      30             : 
      31             : ! namelist variables
      32             :   logical :: nlte_use_mo              ! Determines which constituents are used from NLTE calculations
      33             : !  = .true. uses MOZART constituents
      34             : !  = .false. uses constituents from bnd dataset cftgcm
      35             : 
      36             :   logical :: nlte_use_aliarms = .false.
      37             :   integer :: nlte_aliarms_every_X = 0
      38             : 
      39             :   logical :: use_data_o3
      40             :   logical :: use_waccm_forcing = .false.
      41             : 
      42             :   real(r8) :: o3_mw = -huge(1.0_r8)        ! O3 molecular weight
      43             : 
      44             : ! indexes of required constituents in model constituent array
      45             :   integer :: ico2 = -1                      ! CO2 index
      46             :   integer :: io1 = -1                       ! O index
      47             :   integer :: io2 = -1                       ! O2 index
      48             :   integer :: io3 = -1                       ! O3 index
      49             :   integer :: ih = -1                        ! H index
      50             :   integer :: ino = -1                       ! NO index
      51             :   integer :: qrlaliarms_idx = -1
      52             : 
      53             : ! merge limits for data ozone
      54             :   integer :: nbot_mlt = huge(1)              ! bottom of pure tgcm range
      55             :   integer :: ntop_cam = huge(1)              ! bottom of merge range
      56             :   real(r8):: wt_o3_mrg(pver) = -huge(1.0_r8) ! merge weights for cam o3
      57             : 
      58             : !================================================================================================
      59             : contains
      60             : !================================================================================================
      61             : 
      62           0 :   subroutine nlte_register()
      63             :   use physics_buffer,   only: pbuf_add_field, dtype_r8
      64             : 
      65           0 :   call pbuf_add_field('qrlaliarms',  'global', dtype_r8, (/pcols,pver/),qrlaliarms_idx)
      66             : 
      67           0 :   end subroutine nlte_register
      68             : 
      69             : !================================================================================================
      70             : 
      71        1536 :   subroutine nlte_init (pref_mid, max_pressure_lw, nlte_use_mo_in, nlte_limit_co2, nlte_use_aliarms_in, nlte_aliarms_every_X_in)
      72             : !
      73             : ! Initialize the nlte parameterizations and tgcm forcing data, if required
      74             : !------------------------------------------------------------------------
      75           0 :     use constituents, only: cnst_mw, cnst_get_ind
      76             :     use physconst,    only: mwco2
      77             :     use cam_history,  only: add_default, addfld
      78             :     use mo_waccm_hrates,  only: has_hrates
      79             :     use phys_control, only: phys_getopts
      80             : 
      81             :     real(r8),         intent(in) :: pref_mid(plev)
      82             :     real(r8),         intent(in) :: max_pressure_lw
      83             :     logical,          intent(in) :: nlte_use_mo_in
      84             :     logical,          intent(in) :: nlte_limit_co2
      85             :     logical,          intent(in) :: nlte_use_aliarms_in
      86             :     integer,          intent(in) :: nlte_aliarms_every_X_in
      87             : 
      88             : 
      89             :     real(r8) :: o1_mw = -huge(1.0_r8)      ! O molecular weight
      90             :     real(r8) :: o2_mw = -huge(1.0_r8)      ! O2 molecular weight
      91             :     real(r8) :: co2_mw = -huge(1.0_r8)     ! CO2 molecular weight
      92             :     real(r8) :: n2_mw = -huge(1.0_r8)      ! N2 molecular weight
      93             :     real(r8) :: no_mw = -huge(1.0_r8)      ! NO molecular weight
      94             :     real(r8) :: psh(pver)                  ! pressure scale height
      95             :     real(r8) :: pshmn                      ! lower range of merge
      96             :     real(r8) :: pshmx                      ! upper range of merge
      97             :     real(r8) :: pshdd                      ! scale
      98             :     integer  :: k                          ! index
      99             :     logical  :: rad_use_data_o3
     100             :     logical  :: history_waccm
     101             : !----------------------------------------------------------------------------------------
     102             : 
     103        1536 :     call phys_getopts(history_waccm_out=history_waccm)
     104             : 
     105             : ! Set flag to use mozart (or tgcm) consituents and flag to use ALI-ARMS scheme
     106        1536 :     nlte_use_mo          = nlte_use_mo_in
     107        1536 :     nlte_use_aliarms     = nlte_use_aliarms_in
     108        1536 :     nlte_aliarms_every_X = nlte_aliarms_every_X_in
     109             : 
     110             :     ! ask rad_constituents module whether the O3 used in the climate
     111             :     ! calculation is from data
     112        1536 :     call rad_cnst_get_info(0, use_data_o3=rad_use_data_o3)
     113             : 
     114             :     ! Use data ozone if nlte_use_mo=false, or if nlte_use_mo=true and the flag to use data ozone
     115             :     ! for the interactive radiation calculation has been set to .true. in the rad_constituents module
     116        1536 :     use_data_o3 = .false.
     117        1536 :     if ( .not. nlte_use_mo  .or. &
     118           0 :          (nlte_use_mo .and. rad_use_data_o3) ) use_data_o3 = .true.
     119             : 
     120             : ! Define merge weights for data ozone
     121        1536 :     if (use_data_o3) then
     122           0 :        pshmn=7.0_r8
     123           0 :        pshmx=8.5_r8
     124           0 :        pshdd=1.0_r8
     125             : 
     126           0 :        nbot_mlt = 0
     127           0 :        ntop_cam = 0
     128           0 :        do k = 1, plev
     129           0 :           psh(k) = log(1e5_r8/pref_mid(k))
     130           0 :           if (psh(k) >= pshmx) nbot_mlt = k
     131           0 :           if (psh(k) >= pshmn) ntop_cam = k+1
     132             :        end do
     133             : 
     134           0 :        wt_o3_mrg(:) = 0._r8
     135           0 :        do k = nbot_mlt+1, ntop_cam-1
     136           0 :           wt_o3_mrg(k) = 1._r8 - tanh( (psh(k)-pshmn)/pshdd )
     137             :        enddo
     138           0 :        write(iulog,*) 'NLTE data ozone merge range is ', nbot_mlt+1, ntop_cam-1
     139           0 :        write(iulog,*) 'NLTE data ozone merge weights are ', wt_o3_mrg(nbot_mlt+1 : ntop_cam-1)
     140             : 
     141           0 :        call addfld ('O3MRG',(/ 'lev' /), 'A','mol/mol','merged (eUV+CAM) O3 vmr')
     142             : 
     143             :     end if
     144             : 
     145             : ! Get molecular weights and constituent indexes
     146        1536 :     if (nlte_use_mo)  then
     147             : 
     148        1536 :        call cnst_get_ind( 'CO2', ico2 )
     149        1536 :        call cnst_get_ind( 'O',   io1 )
     150        1536 :        call cnst_get_ind( 'O2',  io2 )
     151        1536 :        call cnst_get_ind( 'O3',  io3 )
     152        1536 :        call cnst_get_ind( 'H',   ih  )
     153        1536 :        call cnst_get_ind( 'NO',  ino )
     154             : 
     155        1536 :        co2_mw= cnst_mw(ico2)
     156        1536 :        o1_mw = cnst_mw(io1)
     157        1536 :        o2_mw = cnst_mw(io2)
     158        1536 :        o3_mw = cnst_mw(io3)
     159        1536 :        no_mw = cnst_mw(ino)
     160        1536 :        n2_mw = 28._r8
     161             : 
     162             :     else
     163             : 
     164           0 :        co2_mw = mwco2
     165           0 :        o1_mw  = 16._r8
     166           0 :        o2_mw  = 32._r8
     167           0 :        o3_mw  = 48._r8
     168           0 :        no_mw  = 30._r8
     169           0 :        n2_mw  = 28._r8
     170             : 
     171             :     end if
     172             : 
     173        1536 :     use_waccm_forcing = use_data_o3 .or. (.not.nlte_use_mo) .or. (.not. has_hrates)
     174             : 
     175             : ! Initialize Fomichev parameterization
     176        1536 :     call nlte_fomichev_init (co2_mw, n2_mw, o1_mw, o2_mw, o3_mw, no_mw, nlte_limit_co2)
     177             : 
     178             : ! Initialize ALI-ARMS parameterization
     179        1536 :     if (nlte_use_aliarms) then
     180           0 :        call nlte_aliarms_init (max_pressure_lw,co2_mw,n2_mw,o1_mw,o2_mw)
     181             :     end if
     182             : 
     183             : ! Initialize waccm forcing data
     184        1536 :     if (use_waccm_forcing) then
     185           0 :        call waccm_forcing_init ()
     186             :     endif
     187             : 
     188        1536 :     if (masterproc) then
     189             : 
     190           2 :        if (nlte_use_mo) then
     191           2 :           write(iulog,*) 'NLTE constituents are obtained from the MOZART chemistry module'
     192             :        else
     193           0 :           write(iulog,*) 'NLTE constituents are obtained from boundary dataset'
     194             :        endif
     195             :     end if
     196             : 
     197        3072 :     call addfld ('QRLNLTE',(/ 'lev' /), 'A','K/s','Non-LTE LW heating (includes QNO and QO3P)')
     198        3072 :     call addfld ('QNO',    (/ 'lev' /), 'A','K/s','NO cooling')
     199        3072 :     call addfld ('QCO2',   (/ 'lev' /), 'A','K/s','CO2 cooling')
     200        3072 :     call addfld ('QO3',    (/ 'lev' /), 'A','K/s','O3 cooling')
     201        3072 :     call addfld ('QHC2S',  (/ 'lev' /), 'A','K/s','Cooling to Space')
     202        3072 :     call addfld ('QO3P',   (/ 'lev' /), 'A','K/s','O3P cooling')
     203             : 
     204             : ! add output to default output for primary history tapes
     205        1536 :     if (history_waccm) then
     206        1536 :        call add_default ('QRLNLTE', 1, ' ')
     207        1536 :        call add_default ('QNO ', 1, ' ')
     208        1536 :        call add_default ('QCO2', 1, ' ')
     209        1536 :        call add_default ('QO3',  1, ' ')
     210        1536 :        call add_default ('QHC2S',1, ' ')
     211        1536 :        call add_default ('QO3P ', 1, ' ')
     212             :     end if
     213             : 
     214        1536 :   end subroutine nlte_init
     215             : 
     216             : !=======================================================================
     217             : 
     218       32256 :   subroutine nlte_timestep_init(state, pbuf2d)
     219        1536 :     use physics_types,  only: physics_state
     220             :     use ppgrid,         only: begchunk, endchunk
     221             :     use physics_buffer, only: physics_buffer_desc
     222             : 
     223             : !
     224             : ! Time interpolation of waccm forcing fields to the current time
     225             : !
     226             : !------------------------------------------------------------------------
     227             : 
     228             :     type(physics_state), intent(in):: state(begchunk:endchunk)
     229             :     type(physics_buffer_desc), pointer :: pbuf2d(:,:)
     230             : 
     231             : 
     232             : !---------------------------Local workspace--------------------------------------
     233             : 
     234       16128 :     if (use_waccm_forcing) then
     235           0 :        call waccm_forcing_adv (state, pbuf2d)
     236             :     endif
     237             : 
     238       16128 :     return
     239       16128 :   end subroutine nlte_timestep_init
     240             : 
     241             : !================================================================================================
     242             : !================================================================================================
     243             : 
     244       80640 :   subroutine nlte_tend(state, pbuf, qrlf)
     245             : !
     246             : ! Driver for nlte calculations
     247             : !-------------------------------------------------------------------------
     248       16128 :     use physconst,       only: mwdry
     249             :     use air_composition, only: cpairv
     250             :     use physics_types, only: physics_state
     251             :     use physics_buffer, only : physics_buffer_desc
     252             :     use perf_mod,      only: t_startf, t_stopf
     253             :     use cam_history,   only: outfld
     254             :     use physics_buffer,only: pbuf_get_field
     255             :     use time_manager,  only: get_nstep
     256             : 
     257             : ! Arguments
     258             :     type(physics_state), target, intent(in)  :: state   ! Physics state variables
     259             : 
     260             :     type(physics_buffer_desc), pointer :: pbuf(:)
     261             : 
     262             :     real(r8), intent(out) :: qrlf(pcols,pver)   ! nlte longwave heating rate
     263             : 
     264             : ! Local workspace for waccm
     265             :     integer :: lchnk              ! chunk identifier
     266             :     integer :: ncol               ! no. of columns in chunk
     267             : 
     268             :     real(r8) :: nocool (pcols,pver)  ! NO cooling (K/s)
     269             :     real(r8) :: o3pcool (pcols,pver) ! O3P cooling (K/s)
     270             :     real(r8) :: qout (pcols,pver)    ! temp for outfld
     271             :     real(r8) :: co2cool(pcols,pver), o3cool(pcols,pver), c2scool(pcols,pver) ! (K/s)
     272             : 
     273       80640 :     real(r8), pointer :: qrlaliarms(:,:) ! ALI-ARMS NLTE CO2 cooling rate (K/s)
     274             : 
     275             :     real(r8) :: qrlfomichev(pcols,pver) ! Fomichev cooling rate ! (K/s)
     276             : 
     277       80640 :     real(r8), pointer, dimension(:,:) :: xco2mmr  ! CO2 mmr
     278       80640 :     real(r8), pointer, dimension(:,:) :: xommr    ! O   mmr
     279       80640 :     real(r8), pointer, dimension(:,:) :: xo2mmr   ! O2  mmr
     280       80640 :     real(r8), pointer, dimension(:,:) :: xo3mmr   ! O3  mmr
     281       80640 :     real(r8), pointer, dimension(:,:) :: xhmmr    ! H  mmr
     282       80640 :     real(r8), pointer, dimension(:,:) :: xnommr   ! NO mmr
     283       80640 :     real(r8), pointer, dimension(:,:) :: xn2mmr   ! N2  mmr
     284             : 
     285             :     real(r8), target :: n2mmr (pcols,pver)   ! N2  mmr
     286             :     real(r8), target :: o3mrg(pcols,pver)    ! merged O3
     287       80640 :     real(r8), pointer, dimension(:,:) :: to3mmr  ! O3 mmr   (tgcm)
     288             : 
     289             :     integer :: k
     290             :     integer :: nstep
     291             : 
     292             : !------------------------------------------------------------------------
     293             : 
     294       80640 :     lchnk = state%lchnk
     295       80640 :     ncol  = state%ncol
     296             : 
     297             : ! Get radiatively active ozone
     298       80640 :     call rad_cnst_get_gas(0, 'O3', state, pbuf,  xo3mmr)
     299       80640 :     if (use_data_o3) then
     300           0 :        call get_cnst (lchnk, o3=to3mmr)
     301           0 :        call merge_o3 (ncol, xo3mmr, to3mmr, o3mrg)
     302           0 :        qout(:ncol,:) = o3mrg(:ncol,:)*mwdry/o3_mw
     303           0 :        call outfld ('O3MRG', qout, pcols,lchnk)
     304           0 :        xo3mmr => o3mrg(:,:)
     305             :     end if
     306             : 
     307       80640 :     if (nlte_use_mo) then
     308             : 
     309             : ! Get relevant constituents from the chemistry module
     310       80640 :        xco2mmr => state%q(:,:,ico2)
     311       80640 :        xommr   => state%q(:,:,io1)
     312       80640 :        xo2mmr  => state%q(:,:,io2)
     313       80640 :        xhmmr   => state%q(:,:,ih)
     314       80640 :        xnommr  => state%q(:,:,ino)
     315             : 
     316             :     else
     317             : 
     318           0 :        call get_cnst (lchnk, co2=xco2mmr, o1=xommr, o2=xo2mmr, no=xnommr, h=xhmmr)
     319             : 
     320             :     endif
     321             : 
     322     5725440 :     do k = 1,pver
     323   173940480 :        n2mmr (:ncol,k) = 1._r8 - (xommr(:ncol,k) + xo2mmr(:ncol,k) + xhmmr(:ncol,k))
     324             :     enddo
     325       80640 :     xn2mmr  => n2mmr(:,:)
     326             : 
     327             : ! do non-LTE cooling rate calculations
     328             : 
     329       80640 :     call t_startf('nlte_fomichev_calc')
     330             :     call nlte_fomichev_calc (lchnk,ncol,state%pmid,state%pint,state%t, &
     331       80640 :          xo2mmr,xommr,xo3mmr,xn2mmr,xco2mmr,qrlfomichev,co2cool,o3cool,c2scool)
     332       80640 :     call t_stopf('nlte_fomichev_calc')
     333             : 
     334             : 
     335             :  !  Call the optional ALI-ARMS.  Note that this does not replace the fomichev
     336             :  !  call as the other individual cooling rates from fomichev still need to be calculated
     337             : 
     338       80640 :     if (nlte_use_aliarms) then
     339             : 
     340           0 :        call t_startf('nlte_aliarms_calc')
     341             : 
     342           0 :        call pbuf_get_field(pbuf, qrlaliarms_idx, qrlaliarms )
     343             :        ! Only run ALI-ARMS every nlte_aliarms_every_X timesteps
     344           0 :        nstep = get_nstep()
     345           0 :        if (MOD(nstep, nlte_aliarms_every_X) == 0) then
     346           0 :           call nlte_aliarms_calc (lchnk,ncol,state%zm, state%pmid,state%t,xo2mmr,xommr,xn2mmr,xco2mmr,qrlaliarms)
     347             :        end if
     348             : 
     349             :        ! Apply the ALI-ARMS heating rate to the qrlf summation
     350           0 :        qrlf(:ncol,:) = o3cool(:ncol,:) + qrlaliarms(:ncol,:) * cpairv(:ncol,:,lchnk)
     351             : 
     352           0 :        call t_stopf('nlte_aliarms_calc')
     353             : 
     354             :     else
     355    87091200 :        qrlf(:ncol,:) = qrlfomichev(:ncol,:)
     356             :     end if
     357             : 
     358             : 
     359             : ! do NO cooling
     360       80640 :     call nocooling (ncol, state%t, state%pmid, xnommr,xommr,xo2mmr,xo3mmr,xn2mmr,nocool)
     361             : 
     362             : ! do O3P cooling
     363       80640 :     call o3pcooling (ncol, state%t, xommr, o3pcool)
     364             : 
     365     5725440 :     do k = 1,pver
     366    87010560 :        qrlf(:ncol,k) = qrlf(:ncol,k) + nocool(:ncol,k) + o3pcool(:ncol,k)
     367             :     end do
     368             : 
     369    87010560 :     qout(:ncol,:) = nocool(:ncol,:)/cpairv(:ncol,:,lchnk)
     370       80640 :     call outfld ('QNO'    , qout, pcols, lchnk)
     371    87010560 :     qout(:ncol,:) = o3pcool(:ncol,:)/cpairv(:ncol,:,lchnk)
     372       80640 :     call outfld ('QO3P'    , qout, pcols, lchnk)
     373    87010560 :     qout(:ncol,:) = qrlf(:ncol,:)/cpairv(:ncol,:,lchnk)
     374       80640 :     call outfld ('QRLNLTE', qout, pcols, lchnk)
     375             : 
     376    87010560 :     qout(:ncol,:) = co2cool(:ncol,:)/cpairv(:ncol,:,lchnk)
     377       80640 :     call outfld ('QCO2', qout, pcols, lchnk)
     378    87010560 :     qout(:ncol,:) = o3cool(:ncol,:)/cpairv(:ncol,:,lchnk)
     379       80640 :     call outfld ('QO3', qout, pcols, lchnk)
     380    87010560 :     qout(:ncol,:) = c2scool(:ncol,:)/cpairv(:ncol,:,lchnk)
     381       80640 :     call outfld ('QHC2S', qout, pcols, lchnk)
     382             : 
     383      161280 :   end subroutine nlte_tend
     384             : 
     385             : !======================================================================================
     386             : 
     387           0 :   subroutine merge_o3 (ncol, o3cam, o3mlt, o3mrg)
     388             : !
     389             : ! Merges CAM O3 (usually climatology) with mesosphere/lower thermosphere O3 (usually TIME/GCM)
     390             : !
     391             : !------------------Input arguments----------------------------------------------
     392             : 
     393             :     integer,  intent(in)    :: ncol                  ! number of atmospheric columns
     394             :     real(r8), intent(in)    :: o3mlt(pcols,pver)     ! MLT O3 mmr
     395             :     real(r8), intent(in)    :: o3cam(pcols,pver)     ! CAM O3 mmr
     396             :     real(r8), intent(out)   :: o3mrg(pcols,pver)     ! merged product
     397             : 
     398             : !---------------------------Local Workspace--------------------------------------------
     399             : 
     400             :     integer k                                        ! index
     401             : 
     402             : !-------------------------------------------------------------------------------------
     403             : 
     404             : ! combine ozone profiles of TIME/GCM with CAM
     405             : 
     406             : ! load TIME/GCM above NBOT_MLT
     407           0 :     do k = 1, nbot_mlt
     408           0 :        o3mrg(:ncol,k) = o3mlt(:ncol,k)
     409             :     end do
     410             : 
     411             : ! merge
     412           0 :     do k=nbot_mlt+1,ntop_cam-1
     413           0 :        o3mrg(:ncol,k) = (1._r8 - wt_o3_mrg(k)) * o3cam(:ncol,k) + wt_o3_mrg(k) * o3mlt(:ncol,k)
     414             :     end do
     415             : 
     416             : ! load CAM below NTOP_CAM
     417           0 :     do k=ntop_cam,pver
     418           0 :        o3mrg(:ncol,k) = o3cam(:ncol,k)
     419             :     end do
     420             : 
     421       80640 :   end subroutine merge_o3
     422             : 
     423             : end module nlte_lw

Generated by: LCOV version 1.14