Restore spectral scaling to RRTMGP

bld/build-namelist

@@ -755,11 +755,9 @@ if ($rad_pkg =~ /rrtmg/ or $chem =~ /waccm/) {
    # use solar data file as the default for rrtmg and waccm_ma
    add_default($nl, 'solar_irrad_data_file');
 
-   # This option only used by camrt and rrtmg radiation schemes.
-   # The solar spectral scaling is done internal to RRTMGP code.
-   if ($rad_pkg ne 'rrtmgp') {
-       add_default($nl, 'solar_htng_spctrl_scl', 'val'=>'.true.');
-   }
+   # The solar spectral scaling is done based on the distribution from
+   # the solar_irrad_data_file.
+   add_default($nl, 'solar_htng_spctrl_scl', 'val'=>'.true.');
 
 }
 elsif (!$simple_phys) {

doc/ChangeLog

@@ -1,5 +1,88 @@
 ===============================================================
 
+Tag name:
+Originator(s): brianpm, eaton
+Date:
+One-line Summary: Restore spectral scaling to RRTMGP
+Github PR URL:
+
+Purpose of changes (include the issue number and title text for each relevant GitHub issue):
+
+resolve issue #1193 - Restore spectral scaling to RRTMGP
+
+Describe any changes made to build system: none
+
+Describe any changes made to the namelist: none
+
+List any changes to the defaults for the boundary datasets: none
+
+Describe any substantial timing or memory changes: not evaluated
+
+Code reviewed by:
+
+List all files eliminated: none
+
+List all files added and what they do:
+
+src/physics/rrtmgp/rad_solar_var.F90
+. compute scale factors for solar irradiance based on input dataset
+
+List all existing files that have been modified, and describe the changes:
+
+bld/build-namelist
+. change default setting of solar_htng_spctrl_scl to true for rrtmgp
+
+src/physics/rrtmgp/radconstants.F90
+. add module data band2gpt_sw and set using kdist_sw%get_band_lims_gpoint()
+
+src/physics/rrtmgp/radiation.F90
+. radiation_init
+  - add call to rad_solar_var_init
+. radiation_tend
+  - replace code that scales the solar source based on internal RRTMGP
+    spectral distribution by a scaling based on distribution from the 
+    solar_irrad_data_file.
+
+If there were any failures reported from running test_driver.sh on any test
+platform, and checkin with these failures has been OK'd by the gatekeeper,
+then copy the lines from the td.*.status files for the failed tests to the
+appropriate machine below.  All failed tests must be justified.
+
+derecho/intel/aux_cam:
+
+derecho/nvhpc/aux_cam:
+
+izumi/nag/aux_cam:
+
+izumi/gnu/aux_cam:
+
+CAM tag used for the baseline comparison tests if different than previous
+tag:
+
+Summarize any changes to answers, i.e.,
+- what code configurations:
+- what platforms/compilers:
+- nature of change (roundoff; larger than roundoff but same climate; new
+  climate):
+
+If bitwise differences were observed, how did you show they were no worse
+than roundoff?
+
+If this tag changes climate describe the run(s) done to evaluate the new
+climate in enough detail that it(they) could be reproduced, i.e.,
+- source tag (all code used must be in the repository):
+- platform/compilers:
+- configure commandline:
+- build-namelist command (or complete namelist):
+- MSS location of output:
+
+MSS location of control simulations used to validate new climate:
+
+URL for AMWG diagnostics output used to validate new climate:
+
+===============================================================
+===============================================================
+
 Tag name: cam6_4_046
 Originator(s): sjsprecious
 Date: 06 November 2024

src/physics/rrtmgp/rad_solar_var.F90

@@ -0,0 +1,149 @@
+!-------------------------------------------------------------------------------
+! This module uses the solar irradiance data 
+! to provide a spectral scaling factor
+! to approximate the spectral distribution of irradiance
+! when the radiation scheme might use a different solar source function
+!-------------------------------------------------------------------------------
+module rad_solar_var
+
+  use shr_kind_mod ,     only : r8 => shr_kind_r8
+  use radconstants,      only : nswbands, get_sw_spectral_boundaries, band2gpt_sw
+  use solar_irrad_data,  only : sol_irrad, we, nbins, has_spectrum, sol_tsi
+  use solar_irrad_data,  only : do_spctrl_scaling
+  use cam_abortutils,    only : endrun
+  use error_messages,    only : alloc_err
+
+  implicit none
+  save
+
+  private
+  public :: rad_solar_var_init
+  public :: get_variability
+
+  real(r8), allocatable :: irrad(:)           ! solar irradiance at model timestep in each band
+
+  real(r8), allocatable :: radbinmax(:)
+  real(r8), allocatable :: radbinmin(:)
+
+!-------------------------------------------------------------------------------
+contains
+!-------------------------------------------------------------------------------
+
+  subroutine rad_solar_var_init( )
+
+    integer :: ierr
+    integer :: radmax_loc
+
+    if ( do_spctrl_scaling ) then
+
+       if ( .not.has_spectrum ) then
+          call endrun('rad_solar_var_init: solar input file must have irradiance spectrum')
+       endif
+
+       allocate (radbinmax(nswbands),stat=ierr)
+       if (ierr /= 0) then
+          call endrun('rad_solar_var_init: Error allocating space for radbinmax')
+       end if
+
+       allocate (radbinmin(nswbands),stat=ierr)
+       if (ierr /= 0) then
+          call endrun('rad_solar_var_init: Error allocating space for radbinmin')
+       end if
+
+       allocate (irrad(nswbands), stat=ierr)
+       if (ierr /= 0) then
+          call endrun('rad_solar_var_init: Error allocating space for irrad')
+       end if
+
+       call get_sw_spectral_boundaries(radbinmin, radbinmax, 'nm')
+
+       ! Make sure that the far-IR is included, even if radiation grid does not
+       ! extend that far down. 10^5 nm corresponds to a wavenumber of
+       ! 100 cm^-1.
+       radmax_loc = maxloc(radbinmax,1)
+       radbinmax(radmax_loc) = max(100000._r8,radbinmax(radmax_loc))
+
+    endif
+
+  end subroutine rad_solar_var_init
+
+!-------------------------------------------------------------------------------
+!-------------------------------------------------------------------------------
+
+  subroutine get_variability(toa_flux, sfac) 
+
+     ! Arguments 
+     real(r8), intent(in)  :: toa_flux(:,:) ! TOA flux to be scaled (columns,gpts)
+     real(r8), intent(out) :: sfac(:,:)     ! scaling factors (columns,gpts)
+
+     ! Local variables 
+     integer :: i, j, istat, gpt_start, gpt_end, ncols
+     real(r8), allocatable :: scale(:)
+     character(len=*), parameter :: sub = 'get_variability'
+
+     if (do_spctrl_scaling) then 
+
+        ! Determine target irradiance for each band
+        call integrate_spectrum(nbins, nswbands, we, radbinmin, radbinmax, sol_irrad, irrad)
+
+        ncols = size(toa_flux, 1)
+        allocate(scale(ncols), stat=istat)
+        call alloc_err(istat, sub, 'scale', ncols)
+
+        do i = 1, nswbands 
+           gpt_start = band2gpt_sw(1,i) 
+           gpt_end   = band2gpt_sw(2,i) 
+           scale = spread(irrad(i), 1, ncols) / sum(toa_flux(:, gpt_start:gpt_end), dim=2)
+           do j = gpt_start, gpt_end
+              sfac(:,j) = scale
+           end do
+        end do
+
+     else 
+        sfac(:,:) = sol_tsi / spread(sum(toa_flux, 2), 2, size(toa_flux, 2))
+     end if
+  end subroutine get_variability
+
+
+!-------------------------------------------------------------------------------
+! private method.........
+!-------------------------------------------------------------------------------
+
+  subroutine integrate_spectrum( nsrc, ntrg, src_x, min_trg, max_trg, src, trg )
+
+    use mo_util, only : rebin
+
+    implicit none
+
+    !---------------------------------------------------------------
+    !	... dummy arguments
+    !---------------------------------------------------------------
+    integer,  intent(in)  :: nsrc                  ! dimension source array
+    integer,  intent(in)  :: ntrg                  ! dimension target array
+    real(r8), intent(in)  :: src_x(nsrc+1)         ! source coordinates
+    real(r8), intent(in)  :: max_trg(ntrg)         ! target coordinates
+    real(r8), intent(in)  :: min_trg(ntrg)         ! target coordinates
+    real(r8), intent(in)  :: src(nsrc)             ! source array
+    real(r8), intent(out) :: trg(ntrg)             ! target array
+
+    !---------------------------------------------------------------
+    !	... local variables
+    !---------------------------------------------------------------
+    real(r8) :: trg_x(2), targ(1)         ! target coordinates
+    integer  :: i
+
+    do i = 1, ntrg
+
+       trg_x(1) = min_trg(i)
+       trg_x(2) = max_trg(i)
+
+       call rebin( nsrc, 1, src_x, trg_x, src(1:nsrc), targ(:) )
+       ! W/m2/nm --> W/m2
+       trg( i ) = targ(1)*(trg_x(2)-trg_x(1))
+
+    enddo
+
+
+  end subroutine integrate_spectrum
+
+end module rad_solar_var

src/physics/rrtmgp/radconstants.F90

@@ -26,6 +26,9 @@ module radconstants
 
 logical :: wavenumber_boundaries_set = .false.
 
+! First and last g-point for each band.
+integer, public, protected :: band2gpt_sw(2,nswbands)
+
 integer, public, protected :: nswgpts  ! number of SW g-points
 integer, public, protected :: nlwgpts  ! number of LW g-points
 
@@ -104,6 +107,9 @@ subroutine set_wavenumber_bands(kdist_sw, kdist_lw)
    wavenumber_low_shortwave = values(1,:)
    wavenumber_high_shortwave = values(2,:)
 
+   ! First and last g-point for each SW band:
+   band2gpt_sw = kdist_sw%get_band_lims_gpoint()
+
    ! Indices into specific bands
    idx_sw_diag     = get_band_index_by_value('sw', 500.0_r8, 'nm')
    idx_nir_diag    = get_band_index_by_value('sw', 1000.0_r8, 'nm')

src/physics/rrtmgp/radiation.F90

@@ -16,7 +16,6 @@ module radiation
                                pbuf_set_field, pbuf_get_field, pbuf_old_tim_idx
 use camsrfexch,          only: cam_out_t, cam_in_t
 use physconst,           only: cappa, cpair, gravit
-use solar_irrad_data,    only: sol_tsi
 
 use time_manager,        only: get_nstep, is_first_step, is_first_restart_step, &
                                get_curr_calday, get_step_size
@@ -27,6 +26,7 @@ module radiation
 
 use radconstants,        only: nradgas, gasnamelength, gaslist, nswbands, nlwbands, &
                                nswgpts, set_wavenumber_bands
+use rad_solar_var,       only: rad_solar_var_init, get_variability
 
 use cloud_rad_props,     only: cloud_rad_props_init
 
@@ -495,6 +495,7 @@ subroutine radiation_init(pbuf2d)
    ! Set the sw/lw band boundaries in radconstants.  Also sets
    ! indicies of specific bands for diagnostic output and COSP input.
    call set_wavenumber_bands(kdist_sw, kdist_lw)
+   call rad_solar_var_init()
 
    ! The spectral band boundaries need to be set before this init is called.
    call rrtmgp_inputs_init(ktopcam, ktoprad)
@@ -937,8 +938,8 @@ subroutine radiation_tend( &
 
    ! TOA solar flux on RRTMGP g-points
    real(r8), allocatable :: toa_flux(:,:)
-   ! TSI from RRTMGP data (from sum over g-point representation)
-   real(r8) :: tsi_ref
+   ! Scale factors based on spectral distribution from input irradiance dataset
+   real(r8), allocatable :: sfac(:,:)
 
    ! Planck sources for LW.
    type(ty_source_func_lw) :: sources_lw
@@ -1097,6 +1098,7 @@ subroutine radiation_tend( &
 
       allocate( &
          t_sfc(ncol), emis_sfc(nlwbands,ncol), toa_flux(nday,nswgpts),     &
+         sfac(nday,nswgpts),                                               &
          t_rad(ncol,nlay), pmid_rad(ncol,nlay), pint_rad(ncol,nlay+1),     &
          t_day(nday,nlay), pmid_day(nday,nlay), pint_day(nday,nlay+1),     &
          coszrs_day(nday), alb_dir(nswbands,nday), alb_dif(nswbands,nday), &
@@ -1174,8 +1176,8 @@ subroutine radiation_tend( &
                   call stop_on_err(errmsg, sub, 'kdist_sw%gas_optics')
 
                   ! Scale the solar source
-                  tsi_ref = sum(toa_flux(1,:))
-                  toa_flux = toa_flux * sol_tsi * eccf / tsi_ref
+                  call get_variability(toa_flux, sfac)
+                  toa_flux = toa_flux * sfac * eccf
 
                end if
 
@@ -1303,7 +1305,7 @@ subroutine radiation_tend( &
       end if  ! if (dolw)
 
       deallocate( &
-         t_sfc, emis_sfc, toa_flux, t_rad, pmid_rad, pint_rad,  &
+         t_sfc, emis_sfc, toa_flux, sfac, t_rad, pmid_rad, pint_rad,  &
          t_day, pmid_day, pint_day, coszrs_day, alb_dir, alb_dif)
 
       !================!