LCOV - code coverage report
Current view: top level - physics/cam - chem_surfvals.F90 (source / functions) Hit Total Coverage
Test: coverage.info Lines: 127 273 46.5 %
Date: 2024-12-17 22:39:59 Functions: 5 8 62.5 %

          Line data    Source code
       1             : 
       2             : module chem_surfvals
       3             : 
       4             : !-----------------------------------------------------------------------------------
       5             : ! Purpose: Provides greenhouse gas (ghg) values at the Earth's surface.
       6             : !          These values may be time dependent.
       7             : !
       8             : ! Author: Brian Eaton (assembled module from existing scattered code pieces)
       9             : !-----------------------------------------------------------------------------------
      10             : 
      11             :    use shr_kind_mod,   only: r8=>shr_kind_r8
      12             :    use spmd_utils,     only: masterproc
      13             :    use time_manager,   only: get_curr_date, get_start_date, is_end_curr_day, &
      14             :                              timemgr_datediff, get_curr_calday
      15             :    use cam_abortutils, only: endrun
      16             :    use netcdf
      17             :    use error_messages, only: handle_ncerr  
      18             :    use cam_logfile,    only: iulog
      19             :    use m_types,        only: time_ramp
      20             :    use constituents,   only: pcnst
      21             : 
      22             : !-----------------------------------------------------------------------
      23             : !- module boilerplate --------------------------------------------------
      24             : !-----------------------------------------------------------------------
      25             :    implicit none
      26             :    private                   ! Make default access private
      27             :    save
      28             : 
      29             :    public ::&
      30             :       chem_surfvals_readnl,  &! read namelist input
      31             :       chem_surfvals_init,    &! initialize options that depend on namelist input
      32             :       chem_surfvals_set,     &! set ghg surface values when scenario_ghg is 'RAMPED' or 'CHEM_LBC_FILE'
      33             :       chem_surfvals_get,     &! return surface values for: CO2VMR, CO2MMR, CH4VMR
      34             :                               ! N2OVMR, F11VMR, and F12VMR
      35             :       chem_surfvals_co2_rad   ! return co2 for radiation
      36             : 
      37             :    public :: flbc_list
      38             : 
      39             : ! Private module data
      40             : 
      41             :    ! Default values for namelist variables -- now set by build-namelist
      42             :    real(r8) :: o2mmr = .23143_r8               ! o2 mass mixing ratio
      43             :    real(r8) :: co2vmr_rad = -1.0_r8            ! co2 vmr override for radiation
      44             :    real(r8) :: co2vmr = -1.0_r8                ! co2   volume mixing ratio 
      45             :    real(r8) :: n2ovmr = -1.0_r8                ! n2o   volume mixing ratio 
      46             :    real(r8) :: ch4vmr = -1.0_r8                ! ch4   volume mixing ratio 
      47             :    real(r8) :: f11vmr = -1.0_r8                ! cfc11 volume mixing ratio 
      48             :    real(r8) :: f12vmr = -1.0_r8                ! cfc12 volume mixing ratio 
      49             :    character(len=16) :: scenario_ghg = 'FIXED' ! 'FIXED','RAMPED', 'RAMP_CO2_ONLY', 'CHEM_LBC_FILE'
      50             :    integer  :: rampYear_ghg = 0                ! ramped gases fixed at this year (if > 0)
      51             :    character(len=256) :: bndtvghg = 'NONE'     ! filename for ramped data
      52             :    integer  :: ramp_co2_start_ymd = 0          ! start date for co2 ramping (yyyymmdd)
      53             :    real(r8) :: ramp_co2_annual_rate = 1.0_r8      ! % amount of co2 ramping per yr; default is 1% 
      54             :    real(r8) :: ramp_co2_cap = -9999.0_r8          ! co2 ramp cap if rate>0, floor otherwise 
      55             :                                                ! as multiple or fraction of inital value
      56             :                                                ! ex. 4.0 => cap at 4x initial co2 setting 
      57             :    integer  :: ghg_yearStart_model = 0         ! model start year
      58             :    integer  :: ghg_yearStart_data  = 0         ! data  start year   
      59             : 
      60             :    logical  :: ghg_use_calendar                ! true => data year = model year
      61             :    logical :: doRamp_ghg    ! true => turn on ramping for ghg
      62             :    logical :: ramp_just_co2 ! true => ramping to be done just for co2 and not other ghg's
      63             :    integer :: fixYear_ghg   ! year at which Ramped gases are fixed
      64             :    integer :: co2_start     ! date at which co2 begins ramping
      65             :    real(r8) :: co2_daily_factor    ! daily multiplier to achieve annual rate of co2 ramp
      66             :    real(r8) :: co2_limit    ! value of co2vmr where ramping ends
      67             :    real(r8) :: co2_base     ! initial co2 volume mixing ratio, before any ramping
      68             :    integer :: ntim = -1               ! number of yearly data values
      69             :    integer,  allocatable :: yrdata(:) ! yearly data values
      70             :    real(r8), allocatable :: co2(:)    ! co2 mixing ratios in ppmv 
      71             :    real(r8), allocatable :: ch4(:)    ! ppbv
      72             :    real(r8), allocatable :: n2o(:)    ! ppbv
      73             :    real(r8), allocatable :: f11(:)    ! pptv
      74             :    real(r8), allocatable :: f12(:)    ! pptv
      75             :    real(r8), allocatable :: adj(:)    ! unitless adjustment factor for f11 & f12
      76             :    
      77             :    ! fixed lower boundary 
      78             :    
      79             :    character(len=256) :: flbc_file = 'NONE'
      80             :    character(len=16)  :: flbc_list(pcnst) = ''
      81             :    type(time_ramp)    :: flbc_timing     != time_ramp( "CYCLICAL",  19970101, 0 )
      82             : 
      83             : !=========================================================================================
      84             : contains
      85             : !=========================================================================================
      86             : 
      87        1536 : subroutine chem_surfvals_readnl(nlfile)
      88             : 
      89             :    ! Read chem_surfvals_nl namelist group.
      90             : 
      91             :    use namelist_utils,  only: find_group_name
      92             :    use units,           only: getunit, freeunit
      93             :    use spmd_utils,      only: mpicom, mstrid=>masterprocid, mpi_logical, mpi_integer, &
      94             :                               mpi_real8, mpi_character
      95             : 
      96             :    character(len=*), intent(in) :: nlfile  ! filepath for file containing namelist input
      97             : 
      98             :    ! Local variables
      99             :    integer :: unitn, ierr, i
     100             :    character(len=*), parameter :: sub = 'chem_surfvals_readnl'
     101             :    
     102             :    character(len=8)   :: flbc_type = 'CYCLICAL'     ! 'CYCLICAL' | 'SERIAL' | 'FIXED'
     103             :    integer            :: flbc_cycle_yr = 0
     104             :    integer            :: flbc_fixed_ymd = 0
     105             :    integer            :: flbc_fixed_tod = 0
     106             : 
     107             :    namelist /chem_surfvals_nl/ co2vmr, n2ovmr, ch4vmr, f11vmr, f12vmr, &
     108             :                                co2vmr_rad, scenario_ghg, rampyear_ghg, bndtvghg, &
     109             :                                ramp_co2_start_ymd, ramp_co2_annual_rate, ramp_co2_cap, &
     110             :                                ghg_yearStart_model, ghg_yearStart_data
     111             : 
     112             :    namelist /chem_surfvals_nl/ flbc_type, flbc_cycle_yr, flbc_fixed_ymd, flbc_fixed_tod, flbc_list, flbc_file
     113             : 
     114             :    !-----------------------------------------------------------------------------
     115             : 
     116        1536 :    if (masterproc) then
     117           2 :       unitn = getunit()
     118           2 :       open( unitn, file=trim(nlfile), status='old' )
     119           2 :       call find_group_name(unitn, 'chem_surfvals_nl', status=ierr)
     120           2 :       if (ierr == 0) then
     121           2 :          read(unitn, chem_surfvals_nl, iostat=ierr)
     122           2 :          if (ierr /= 0) then
     123           0 :             call endrun(sub // ':: ERROR reading namelist')
     124             :          end if
     125             :       end if
     126           2 :       close(unitn)
     127           2 :       call freeunit(unitn)
     128             :    end if
     129             : 
     130        1536 :    call mpi_bcast(co2vmr, 1, mpi_real8, mstrid, mpicom, ierr)
     131        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: co2vmr")
     132        1536 :    call mpi_bcast(n2ovmr, 1, mpi_real8, mstrid, mpicom, ierr)
     133        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: n2ovmr")
     134        1536 :    call mpi_bcast(ch4vmr, 1, mpi_real8, mstrid, mpicom, ierr)
     135        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: ch4vmr")
     136        1536 :    call mpi_bcast(f11vmr, 1, mpi_real8, mstrid, mpicom, ierr)
     137        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: f11vmr")
     138        1536 :    call mpi_bcast(f12vmr, 1, mpi_real8, mstrid, mpicom, ierr)
     139        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: f12vmr")
     140        1536 :    call mpi_bcast(co2vmr_rad, 1, mpi_real8, mstrid, mpicom, ierr)
     141        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: co2vmr_rad")
     142        1536 :    call mpi_bcast(scenario_ghg, len(scenario_ghg), mpi_character, mstrid, mpicom, ierr)
     143        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: scenario_ghg")
     144        1536 :    call mpi_bcast(rampyear_ghg, 1, mpi_integer, mstrid, mpicom, ierr)
     145        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: rampyear_ghg")
     146        1536 :    call mpi_bcast(bndtvghg, len(bndtvghg), mpi_character, mstrid, mpicom, ierr)
     147        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: bndtvghg")
     148        1536 :    call mpi_bcast(ramp_co2_start_ymd, 1, mpi_integer, mstrid, mpicom, ierr)
     149        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: ramp_co2_start_ymd")
     150        1536 :    call mpi_bcast(ramp_co2_annual_rate, 1, mpi_real8, mstrid, mpicom, ierr)
     151        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: ramp_co2_annual_rate")
     152        1536 :    call mpi_bcast(ramp_co2_cap, 1, mpi_real8, mstrid, mpicom, ierr)
     153        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: ramp_co2_cap")
     154        1536 :    call mpi_bcast(ghg_yearstart_model, 1, mpi_integer, mstrid, mpicom, ierr)
     155        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: ghg_yearstart_model")
     156        1536 :    call mpi_bcast(ghg_yearstart_data, 1, mpi_integer, mstrid, mpicom, ierr)
     157        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: ghg_yearstart_data")
     158        1536 :    call mpi_bcast(flbc_type, len(flbc_type), mpi_character, mstrid, mpicom, ierr)
     159        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: flbc_type")
     160        1536 :    call mpi_bcast(flbc_cycle_yr, 1, mpi_integer, mstrid, mpicom, ierr)
     161        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: flbc_cycle_yr")
     162        1536 :    call mpi_bcast(flbc_fixed_ymd, 1, mpi_integer, mstrid, mpicom, ierr)
     163        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: flbc_fixed_ymd")
     164        1536 :    call mpi_bcast(flbc_fixed_tod, 1, mpi_integer, mstrid, mpicom, ierr)
     165        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: flbc_fixed_tod")
     166        1536 :    call mpi_bcast(flbc_list, len(flbc_list(1))*pcnst, mpi_character, mstrid, mpicom, ierr)
     167        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: flbc_list")
     168        1536 :    call mpi_bcast(flbc_file, len(flbc_file), mpi_character, mstrid, mpicom, ierr)
     169        1536 :    if (ierr /= 0) call endrun(sub//" mpi_bcast: flbc_file")
     170             : 
     171        1536 :    flbc_timing%type      = flbc_type
     172        1536 :    flbc_timing%cycle_yr  = flbc_cycle_yr
     173        1536 :    flbc_timing%fixed_ymd = flbc_fixed_ymd
     174        1536 :    flbc_timing%fixed_tod = flbc_fixed_tod
     175             : 
     176        1536 :    if ( (bndtvghg.ne.'NONE') .and. (flbc_file.ne.'NONE') ) then
     177           0 :       call endrun(sub//': Cannot specify both bndtvghg and flbc_file ')
     178             :    endif
     179             : 
     180        1536 :    if (masterproc) then
     181           2 :       write(iulog,*) ' '
     182           2 :       write(iulog,*) sub//': Settings for control of GHG surface values '
     183           2 :       write(iulog,*) '  scenario_ghg = '//trim(scenario_ghg)
     184             :       
     185           2 :       if (scenario_ghg == 'FIXED' .or. scenario_ghg == 'RAMP_CO2_ONLY') then
     186             : 
     187           0 :          if (scenario_ghg == 'RAMP_CO2_ONLY') then
     188           0 :             write(iulog,*) '  CO2 will be ramped as follows:'
     189           0 :             write(iulog,*) '    Initial co2vmr       = ', co2vmr
     190           0 :             write(iulog,*) '    ramp_co2_annual_rate = ', ramp_co2_annual_rate
     191           0 :             write(iulog,*) '    ramp_co2_cap         = ', ramp_co2_cap
     192           0 :             if (ramp_co2_start_ymd == 0) then
     193           0 :                write(iulog,*) '    ramp starts at initial run time '
     194             :             else
     195           0 :                write(iulog,*) '    ramp_co2_start_ymd   = ', ramp_co2_start_ymd
     196             :             end if
     197             :  
     198             :          else
     199           0 :             write(iulog,*) '  CO2 will be fixed:'
     200           0 :             write(iulog,*) '    co2vmr = ', co2vmr
     201             :          end if
     202             : 
     203           0 :          write(iulog,*) '  Other GHG values will be fixed as follows:'
     204           0 :          write(iulog,*) '    n2ovmr = ', n2ovmr
     205           0 :          write(iulog,*) '    ch4vmr = ', ch4vmr
     206           0 :          write(iulog,*) '    f11vmr = ', f11vmr
     207           0 :          write(iulog,*) '    f12vmr = ', f12vmr
     208             : 
     209           2 :       else if (scenario_ghg == 'RAMPED') then
     210           0 :          write(iulog,*) '  GHG values time interpolated from global annual averages:'
     211           0 :          write(iulog,*) '    bndtvghg = ', trim(bndtvghg)
     212           0 :          if (rampYear_ghg > 0) then
     213           0 :             write(iulog,*) '  Data from bndtvghg FIXED at year ', rampYear_ghg
     214             :          end if
     215             : 
     216           2 :       else if (scenario_ghg == 'CHEM_LBC_FILE') then
     217           2 :          write(iulog,*) '  GHG values time interpolated from LBC file:'
     218           2 :          write(iulog,*) '    flbc_file = ', trim(flbc_file)
     219           2 :          write(iulog,*) '    flbc_type = ', trim(flbc_type)
     220           2 :          if (flbc_type == 'CYCLICAL') then
     221           0 :             write(iulog,*) '    flbc_cycle_yr = ', flbc_cycle_yr
     222           2 :          else if (flbc_type == 'FIXED') then
     223           0 :             write(iulog,*) '    flbc_fixed_ymd = ', flbc_fixed_ymd
     224           0 :             write(iulog,*) '    flbc_fixed_tod = ', flbc_fixed_tod
     225             :          end if
     226           2 :          write(iulog,*) '  Species from LBC file:'
     227          14 :          do i = 1, pcnst
     228          14 :             if (flbc_list(i) == ' ') exit
     229          14 :             write(iulog,*) '    '//trim(flbc_list(i))
     230             :          end do
     231             : 
     232             :       else
     233             :          call endrun (sub//': scenario_ghg must be set to either FIXED, RAMPED, RAMP_CO2_ONLY, &
     234           0 :                    & or CHEM_LBC_FILE')
     235             : 
     236             :       end if
     237             : 
     238           2 :       if (co2vmr_rad > 0._r8) then
     239           0 :          write(iulog,*) '  *** The CO2 prescribed by scenario_ghg has been OVERRIDDEN by ***'
     240           0 :          write(iulog,*) '      co2vmr_rad = ', co2vmr_rad
     241             :       end if
     242             : 
     243           2 :       write(iulog,*) ' '
     244             :    end if
     245             : 
     246        1536 : end subroutine chem_surfvals_readnl
     247             : 
     248             : !================================================================================================
     249             : 
     250        1536 : subroutine chem_surfvals_init()
     251             : 
     252             : !----------------------------------------------------------------------- 
     253             : ! 
     254             : ! Purpose: 
     255             : ! Initialize the ramp options that are controlled by namelist input.
     256             : ! Set surface values at initial time.
     257             : ! N.B. This routine must be called after the time manager has been initialized
     258             : !      since chem_surfvals_set calls time manager methods.
     259             : ! 
     260             : ! Author: B. Eaton - merged code from parse_namelist and rampnl_ghg.
     261             : ! 
     262             : !-----------------------------------------------------------------------
     263             : 
     264             :    use infnan,       only: posinf, assignment(=)
     265             :    use mo_flbc,      only: flbc_inti
     266             :    use phys_control, only: use_simple_phys
     267             : 
     268             :    !---------------------------Local variables-----------------------------
     269             :    integer :: yr, mon, day, ncsec
     270             :    character(len=*), parameter :: sub = 'chem_surfvals_init'
     271             :    !-----------------------------------------------------------------------
     272             : 
     273        1536 :    if (use_simple_phys) return
     274             : 
     275        1536 :    if (scenario_ghg == 'FIXED') then
     276           0 :       doRamp_ghg = .false.
     277           0 :       ramp_just_co2 = .false.
     278             : 
     279        1536 :    else if (scenario_ghg == 'RAMPED') then
     280           0 :       doRamp_ghg = .true.
     281           0 :       ramp_just_co2 = .false.
     282           0 :       call ghg_ramp_read
     283             : 
     284           0 :       fixYear_ghg = rampYear_ghg
     285             : 
     286           0 :       call chem_surfvals_set()
     287             : 
     288        1536 :    else if (scenario_ghg == 'RAMP_CO2_ONLY') then
     289             : 
     290           0 :       if(ramp_co2_start_ymd == 0) then
     291             :          ! by default start the ramp at the initial run time
     292           0 :          call get_start_date(yr, mon, day, ncsec)
     293           0 :          ramp_co2_start_ymd = yr*10000 + mon*100 + day
     294             :       end if
     295           0 :       co2_start = ramp_co2_start_ymd
     296             : 
     297           0 :       if(ramp_co2_annual_rate <= -100.0_r8) then
     298           0 :          write(iulog,*) 'RAMP_CO2:  invalid ramp_co2_annual_rate= ',ramp_co2_annual_rate
     299           0 :          call endrun (sub//': RAMP_CO2_ANNUAL_RATE must be greater than -100.0')
     300             :       end if
     301             : 
     302           0 :       doRamp_ghg = .true.
     303           0 :       ramp_just_co2 = .true.
     304           0 :       co2_base = co2vmr        ! save initial setting 
     305             : 
     306           0 :       co2_daily_factor = (ramp_co2_annual_rate*0.01_r8+1.0_r8)**(1.0_r8/365.0_r8)
     307             : 
     308           0 :       if (ramp_co2_cap > 0.0_r8) then  
     309           0 :          co2_limit = ramp_co2_cap * co2_base
     310             :       else                                  ! if no cap/floor specified, provide default
     311           0 :          if (ramp_co2_annual_rate < 0.0_r8) then
     312           0 :             co2_limit = 0.0_r8
     313             :          else
     314           0 :             co2_limit = posinf
     315             :          end if
     316             :       end if
     317             : 
     318           0 :       if ((ramp_co2_annual_rate<0.0_r8 .and. co2_limit>co2_base) .or. &
     319             :           (ramp_co2_annual_rate>0.0_r8 .and. co2_limit<co2_base)) then
     320           0 :          write(iulog,*) sub//': ramp_co2_cap is unreachable'
     321           0 :          write(iulog,*) sub//': ramp_co2_annual_rate= ',ramp_co2_annual_rate,' ramp_co2_cap= ',ramp_co2_cap
     322           0 :          call endrun(sub//':  ramp_co2_annual_rate and ramp_co2_cap incompatible')
     323             :       end if
     324             : 
     325           0 :       call chem_surfvals_set()
     326             : 
     327        1536 :    else if (scenario_ghg == 'CHEM_LBC_FILE') then
     328             :       ! set by lower boundary conditions file
     329        1536 :       call flbc_inti( flbc_file, flbc_list, flbc_timing, co2vmr, ch4vmr, n2ovmr, f11vmr, f12vmr )
     330        1536 :       call chem_surfvals_set()
     331             : 
     332             :    endif
     333             : 
     334        1536 :    if (masterproc) then
     335           2 :       write(iulog,*) ' '
     336           2 :       write(iulog,*) 'chem_surfvals_init: Initial ghg surface values:'
     337           2 :       write(iulog,*) '  co2 volume mixing ratio = ', chem_surfvals_co2_rad(vmr_in=.true.)
     338           2 :       write(iulog,*) '  ch4 volume mixing ratio = ', ch4vmr
     339           2 :       write(iulog,*) '  n2o volume mixing ratio = ', n2ovmr
     340           2 :       write(iulog,*) '  f11 volume mixing ratio = ', f11vmr
     341           2 :       write(iulog,*) '  f12 volume mixing ratio = ', f12vmr
     342           2 :       write(iulog,*) ' '
     343             :    end if
     344             : 
     345        1536 : end subroutine chem_surfvals_init
     346             : 
     347             : !=========================================================================================
     348             : 
     349           0 : subroutine ghg_ramp_read()
     350             : 
     351             : !----------------------------------------------------------------------- 
     352             : ! 
     353             : ! Purpose: 
     354             : ! Read ramped greenhouse gas surface data.  
     355             : ! 
     356             : ! Author: T. Henderson
     357             : ! 
     358             : !-----------------------------------------------------------------------
     359             : 
     360        1536 :    use ioFileMod, only: getfil
     361             : #if ( defined SPMD )
     362             :    use mpishorthand, only: mpicom, mpiint, mpir8
     363             : #endif
     364             :    character(len=*), parameter :: subname = 'ghg_ramp_read'
     365             : 
     366             : !---------------------------Local variables-----------------------------
     367             :    integer :: ncid
     368             :    integer :: co2_id
     369             :    integer :: ch4_id
     370             :    integer :: n2o_id
     371             :    integer :: f11_id
     372             :    integer :: f12_id
     373             :    integer :: adj_id
     374             :    integer :: date_id
     375             :    integer :: time_id
     376             :    integer :: ierror
     377             :    character(len=256) :: locfn          ! netcdf local filename to open
     378             : 
     379           0 :    if (masterproc) then
     380           0 :      call getfil (bndtvghg, locfn, 0)
     381           0 :      call handle_ncerr( nf90_open (trim(locfn), NF90_NOWRITE, ncid),subname,__LINE__)
     382             : 
     383           0 :      write(iulog,*)'GHG_RAMP_READ:  reading ramped greenhouse gas surface data from file ',trim(locfn)
     384             : 
     385           0 :      call handle_ncerr( nf90_inq_varid( ncid, 'date', date_id ),subname,__LINE__)
     386           0 :      call handle_ncerr( nf90_inq_varid( ncid, 'CO2', co2_id ),subname,__LINE__)
     387           0 :      call handle_ncerr( nf90_inq_varid( ncid, 'CH4', ch4_id ),subname,__LINE__)
     388           0 :      call handle_ncerr( nf90_inq_varid( ncid, 'N2O', n2o_id ),subname,__LINE__)
     389           0 :      call handle_ncerr( nf90_inq_varid( ncid, 'f11', f11_id ),subname,__LINE__)
     390           0 :      call handle_ncerr( nf90_inq_varid( ncid, 'f12', f12_id ),subname,__LINE__)
     391           0 :      call handle_ncerr( nf90_inq_varid( ncid, 'adj', adj_id ),subname,__LINE__)
     392           0 :      call handle_ncerr( nf90_inq_dimid( ncid, 'time', time_id ),subname,__LINE__)
     393           0 :      call handle_ncerr( nf90_inquire_dimension( ncid, time_id, len=ntim ),subname,__LINE__)
     394             : 
     395             :    endif
     396             : #if (defined SPMD )
     397           0 :    call mpibcast (ntim, 1, mpiint, 0, mpicom)
     398             : #endif
     399             :    ! these arrays are never deallocated
     400             :    allocate ( yrdata(ntim), co2(ntim), ch4(ntim), n2o(ntim),    &
     401           0 :                  f11(ntim), f12(ntim), adj(ntim), stat=ierror )
     402           0 :    if (ierror /= 0) then
     403           0 :      write(iulog,*)'GHG_RAMP_READ:  ERROR, allocate() failed!'
     404           0 :      call endrun
     405             :    endif
     406           0 :    if (masterproc) then
     407           0 :      call handle_ncerr( nf90_get_var (ncid, date_id, yrdata ),subname,__LINE__)
     408           0 :      yrdata = yrdata / 10000
     409           0 :      call handle_ncerr( nf90_get_var (ncid, co2_id, co2 ),subname,__LINE__)
     410           0 :      call handle_ncerr( nf90_get_var (ncid, ch4_id, ch4 ),subname,__LINE__)
     411           0 :      call handle_ncerr( nf90_get_var (ncid, n2o_id, n2o ),subname,__LINE__)
     412           0 :      call handle_ncerr( nf90_get_var (ncid, f11_id, f11 ),subname,__LINE__)
     413           0 :      call handle_ncerr( nf90_get_var (ncid, f12_id, f12 ),subname,__LINE__)
     414           0 :      call handle_ncerr( nf90_get_var (ncid, adj_id, adj ),subname,__LINE__)
     415           0 :      call handle_ncerr( nf90_close (ncid),subname,__LINE__)
     416           0 :      write(iulog,*)'GHG_RAMP_READ:  successfully read ramped greenhouse gas surface data from years ',&
     417           0 :         yrdata(1),' through ',yrdata(ntim)
     418             :    endif
     419             : #if (defined SPMD )
     420           0 :    call mpibcast (co2, ntim, mpir8, 0, mpicom)
     421           0 :    call mpibcast (ch4, ntim, mpir8, 0, mpicom)
     422           0 :    call mpibcast (n2o, ntim, mpir8, 0, mpicom)
     423           0 :    call mpibcast (f11, ntim, mpir8, 0, mpicom)
     424           0 :    call mpibcast (f12, ntim, mpir8, 0, mpicom)
     425           0 :    call mpibcast (adj, ntim, mpir8, 0, mpicom)
     426           0 :    call mpibcast (yrdata, ntim, mpiint, 0, mpicom)
     427             : #endif
     428             : 
     429           0 :    return
     430             : 
     431             : end subroutine ghg_ramp_read
     432             : 
     433             : !=========================================================================================
     434             : 
     435    15420432 : function chem_surfvals_get(name)
     436             :   use physconst,    only: mwdry, mwco2
     437             : 
     438             :   character(len=*), intent(in) :: name
     439             : 
     440             :   real(r8) :: rmwco2 
     441             :   real(r8) :: chem_surfvals_get
     442             : 
     443    15420432 :   rmwco2 = mwco2/mwdry    ! ratio of molecular weights of co2 to dry air
     444     1495368 :   select case (name)
     445             :   case ('CO2VMR')
     446     1495368 :      chem_surfvals_get = co2vmr
     447             :   case ('CO2MMR')
     448           0 :      chem_surfvals_get = rmwco2 * co2vmr
     449             :   case ('N2OVMR')
     450     3107424 :      chem_surfvals_get = n2ovmr
     451             :   case ('CH4VMR')
     452     3107424 :      chem_surfvals_get = ch4vmr
     453             :   case ('F11VMR')
     454     3107424 :      chem_surfvals_get = f11vmr
     455             :   case ('F12VMR')
     456     3107424 :      chem_surfvals_get = f12vmr
     457             :   case ('O2MMR')
     458     1495368 :      chem_surfvals_get = o2mmr
     459             :   case default
     460    15420432 :      call endrun('chem_surfvals_get does not know name')
     461             :   end select
     462             : 
     463    15420432 : end function chem_surfvals_get
     464             : 
     465             : 
     466             : !=========================================================================================
     467             : 
     468     1618250 : function chem_surfvals_co2_rad(vmr_in)
     469             :  
     470             :    ! Return the value of CO2 (as mmr) that is radiatively active.
     471             : 
     472             :    ! This method is used by ghg_data to set the prescribed value of CO2 in
     473             :    ! the physics buffer.  If the user has set the co2vmr_rad namelist
     474             :    ! variable then that value will override either the value set by the
     475             :    ! co2vmr namelist variable, or the values time interpolated from a
     476             :    ! dataset.
     477             :    
     478             :    ! This method is also used by cam_history to write the radiatively active
     479             :    ! CO2 to the history file.  The optional argument allows returning the
     480             :    ! value as vmr.
     481             : 
     482             :    use physconst,    only: mwdry, mwco2
     483             : 
     484             :    ! Arguments
     485             :    logical, intent(in), optional :: vmr_in  ! return CO2 as vmr
     486             : 
     487             :    ! Return value
     488             :    real(r8) :: chem_surfvals_co2_rad
     489             : 
     490             :    ! Local variables
     491             :    real(r8) :: convert_vmr      ! convert vmr to desired output
     492             :    !-----------------------------------------------------------------------
     493             : 
     494             :    ! by default convert vmr to mmr
     495     1618250 :    convert_vmr = mwco2/mwdry    ! ratio of molecular weights of co2 to dry air
     496     1618250 :    if (present(vmr_in)) then
     497             :       ! if request return vmr
     498      122882 :       if (vmr_in) convert_vmr = 1.0_r8
     499             :    end if
     500             : 
     501     1618250 :    if (co2vmr_rad > 0._r8) then
     502           0 :       chem_surfvals_co2_rad = convert_vmr * co2vmr_rad
     503             :    else                           
     504     1618250 :       chem_surfvals_co2_rad = convert_vmr * co2vmr     
     505             :    end if
     506             : 
     507     1618250 : end function chem_surfvals_co2_rad
     508             : 
     509             : !=========================================================================================
     510             : 
     511      372480 : subroutine chem_surfvals_set()
     512             : 
     513             :    use ppgrid,         only: begchunk, endchunk
     514             :    use mo_flbc,        only: flbc_gmean_vmr, flbc_chk
     515             :    use scamMod,        only: single_column, scmiop_flbc_inti, use_camiop
     516             : 
     517             : !---------------------------Local variables-----------------------------
     518             : 
     519             :    integer  :: yr, mon, day, ncsec ! components of a date
     520             :    integer  :: ncdate              ! current date in integer format [yyyymmdd]
     521             :    
     522      372480 :    if ( doRamp_ghg ) then
     523           0 :       if(ramp_just_co2) then
     524           0 :          call chem_surfvals_set_co2()
     525             :       else
     526           0 :          call chem_surfvals_set_all()
     527             :       end if
     528      372480 :    elseif (scenario_ghg == 'CHEM_LBC_FILE') then
     529             :       ! set mixing ratios from cam-chem/waccm lbc file 
     530      372480 :       call flbc_chk()
     531      372480 :       if (single_column .and. use_camiop) then
     532           0 :          call scmiop_flbc_inti( co2vmr, ch4vmr, n2ovmr, f11vmr, f12vmr )
     533             :       else
     534             :          ! set by lower boundary conditions file
     535      372480 :          call flbc_gmean_vmr(co2vmr,ch4vmr,n2ovmr,f11vmr,f12vmr)
     536             :       endif
     537             :    endif
     538             : 
     539      372480 :    if (masterproc .and. is_end_curr_day()) then
     540          12 :       call get_curr_date(yr, mon, day, ncsec)
     541          12 :       ncdate = yr*10000 + mon*100 + day
     542          12 :       write(iulog,*) 'chem_surfvals_set: ncdate= ',ncdate,' co2vmr=',co2vmr
     543             : 
     544          12 :       if (.not. ramp_just_co2 .and. mon==1 .and. day==1) then
     545           2 :          write(iulog,*) 'chem_surfvals_set: ch4vmr=', ch4vmr, ' n2ovmr=', n2ovmr, &
     546           4 :                         ' f11vmr=', f11vmr, ' f12vmr=', f12vmr
     547             :       end if
     548             : 
     549             :    end if
     550             : 
     551      372480 :    return
     552      372480 : end subroutine chem_surfvals_set
     553             : 
     554             : !=========================================================================================
     555             : 
     556           0 : subroutine chem_surfvals_set_all()
     557             : !----------------------------------------------------------------------- 
     558             : ! 
     559             : ! Purpose: 
     560             : ! Computes greenhouse gas volume mixing ratios via interpolation of
     561             : ! yearly input data.
     562             : ! 
     563             : ! Author: B. Eaton - updated ramp_ghg for use in chem_surfvals module
     564             : ! 
     565             : !-----------------------------------------------------------------------
     566      372480 :    use interpolate_data, only: get_timeinterp_factors
     567             : 
     568             : !---------------------------Local variables-----------------------------
     569             : 
     570             :    integer yrmodel           ! model year
     571             :    integer nyrm              ! year index
     572             :    integer nyrp              ! year index
     573             :    integer :: yr, mon, day   ! components of a date
     574             :    integer :: ncdate         ! current date in integer format [yyyymmdd]
     575             :    integer :: ncsec          ! current time of day [seconds]
     576             : 
     577             :    real(r8) :: calday            ! current calendar day
     578             :    real(r8) doymodel             ! model day of year
     579             :    real(r8) doydatam             ! day of year for input data yrdata(nyrm)
     580             :    real(r8) doydatap             ! day or year for input data yrdata(nyrp)
     581             :    real(r8) deltat               ! delta time
     582             :    real(r8) fact1, fact2         ! time interpolation factors
     583             :    real(r8) cfcscl               ! cfc scale factor for f11
     584             : 
     585             :    integer yearRan_model         ! model ran year
     586             : !
     587             : ! ---------------------------------------------------------------------
     588             : !
     589           0 :    calday = get_curr_calday()
     590           0 :    call get_curr_date(yr, mon, day, ncsec)
     591           0 :    ncdate = yr*10000 + mon*100 + day
     592             : !
     593             : ! determine ghg_use_calendar      
     594             : !
     595           0 :    if ( ghg_yearStart_model > 0 .and. ghg_yearStart_data > 0 ) then
     596           0 :       ghg_use_calendar = .false.
     597             :    else
     598           0 :       ghg_use_calendar = .true.
     599             :    end if
     600             : !
     601             : ! determine index into input data
     602             : !
     603           0 :    if ( fixYear_ghg > 0) then
     604           0 :       yrmodel  = fixYear_ghg
     605           0 :       nyrm = fixYear_ghg - yrdata(1) + 1
     606             :    else
     607           0 :       if ( ghg_use_calendar) then
     608           0 :          yrmodel  = yr          
     609           0 :          nyrm = yr - yrdata(1) + 1
     610             :       else 
     611           0 :          yearRan_model = yr - ghg_yearStart_model
     612           0 :          if ( yearRan_model < 0 ) then
     613           0 :             call endrun('chem_surfvals_set_all: incorrect ghg_yearStart_model')
     614             :          endif
     615           0 :          yrmodel  = yearRan_model + ghg_yearStart_data
     616             :  
     617           0 :          nyrm = ghg_yearStart_data + yearRan_model - yrdata(1) + 1
     618             :       end if
     619             :    end if
     620             : 
     621           0 :    nyrp       = nyrm + 1
     622             : !
     623             : ! if current date is before yrdata(1), quit
     624             : !
     625           0 :    if (nyrm < 1) then
     626           0 :       write(iulog,*)'chem_surfvals_set_all: data time index is out of bounds'
     627           0 :       write(iulog,*)'nyrm = ',nyrm,' nyrp= ',nyrp, ' ncdate= ', ncdate
     628           0 :       call endrun
     629             :    endif
     630             : !
     631             : ! if current date later than yrdata(ntim), call endrun.
     632             : ! if want to use ntim values - uncomment the following lines
     633             : ! below and comment the call to endrun and previous write
     634             : !
     635           0 :    if (nyrp > ntim) then
     636           0 :       call endrun ('chem_surfvals_set_all: error - current date is past the end of valid data')
     637             : !         write(iulog,*)'chem_surfvals_set_all: using ghg data for ',yrdata(ntim)
     638             : !         co2vmr = co2(ntim)*1.e-06
     639             : !         ch4vmr = ch4(ntim)*1.e-09
     640             : !         n2ovmr = n2o(ntim)*1.e-09
     641             : !         f11vmr = f11(ntim)*1.e-12*(1.+cfcscl)
     642             : !         f12vmr = f12(ntim)*1.e-12
     643             : !         co2mmr = rmwco2 * co2vmr
     644             : !         return
     645             :    endif
     646             : !
     647             : ! determine time interpolation factors, check sanity
     648             : ! of interpolation factors to within 32-bit roundoff
     649             : ! assume that day of year is 1 for all input data
     650             : !
     651           0 :    doymodel = yrmodel*365._r8    + calday
     652           0 :    doydatam = yrdata(nyrm)*365._r8 + 1._r8
     653           0 :    doydatap = yrdata(nyrp)*365._r8 + 1._r8
     654             : 
     655             :    call get_timeinterp_factors(.false.,2,doydatam,doydatap, doymodel, &
     656           0 :         fact1, fact2,'chem_surfvals')
     657             : 
     658             : !
     659             : ! do time interpolation:
     660             : !   co2     in ppmv
     661             : !   n2o,ch4 in ppbv
     662             : !   f11,f12 in pptv
     663             : !
     664           0 :    co2vmr = (co2(nyrm)*fact1 + co2(nyrp)*fact2)*1.e-06_r8
     665           0 :    ch4vmr = (ch4(nyrm)*fact1 + ch4(nyrp)*fact2)*1.e-09_r8
     666           0 :    n2ovmr = (n2o(nyrm)*fact1 + n2o(nyrp)*fact2)*1.e-09_r8
     667             : 
     668           0 :    cfcscl = (adj(nyrm)*fact1 + adj(nyrp)*fact2)
     669           0 :    f11vmr = (f11(nyrm)*fact1 + f11(nyrp)*fact2)*1.e-12_r8*(1._r8+cfcscl)
     670           0 :    f12vmr = (f12(nyrm)*fact1 + f12(nyrp)*fact2)*1.e-12_r8
     671             : 
     672           0 :    return
     673             : end subroutine chem_surfvals_set_all
     674             : 
     675             : !=========================================================================================
     676             : 
     677           0 : subroutine chem_surfvals_set_co2()
     678             : !----------------------------------------------------------------------- 
     679             : ! 
     680             : ! Purpose: 
     681             : ! Computes co2 greenhouse gas volume mixing ratio via ramping info 
     682             : ! provided in namelist var's
     683             : ! 
     684             : ! Author: B. Eaton - updated ramp_ghg for use in chem_surfvals module
     685             : ! 
     686             : !-----------------------------------------------------------------------
     687             :    use shr_kind_mod, only: r8 => shr_kind_r8
     688             : 
     689             : !---------------------------Local variables-----------------------------
     690             : 
     691             :    real(r8) :: daydiff             ! number of days of co2 ramping
     692             :    integer  :: yr, mon, day, ncsec ! components of a date
     693             :    integer  :: ncdate              ! current date in integer format [yyyymmdd]
     694             : !-----------------------------------------------------------------------
     695             : 
     696           0 :    call get_curr_date(yr, mon, day, ncsec)
     697           0 :    ncdate = yr*10000 + mon*100 + day
     698             : 
     699           0 :    call timemgr_datediff(co2_start, 0, ncdate, ncsec, daydiff)
     700             : 
     701           0 :    if (daydiff > 0.0_r8) then
     702             : 
     703           0 :       co2vmr = co2_base*(co2_daily_factor)**daydiff
     704             : 
     705           0 :       if(co2_daily_factor < 1.0_r8) then
     706           0 :          co2vmr = max(co2vmr,co2_limit)
     707             :       else
     708           0 :          co2vmr = min(co2vmr,co2_limit)
     709             :       end if
     710             :    end if
     711             : 
     712           0 :    return
     713             : end subroutine chem_surfvals_set_co2
     714             : 
     715             : 
     716             : !=========================================================================================
     717             : 
     718             : end module chem_surfvals

Generated by: LCOV version 1.14