ESCOMP / ESCOMP/CAM

Turning on spectralflux namelist option doesn't return spectral fluxes

Open
#1,032 3 comments 0 reactions 2 assignees Claimed by @brian-eaton View on GitHub
enhancement misc tag
Dominant language
No language data
Stars
91
Forks
183
Avg merge
6d 2m
Merged PRs (30d)
9

Description

### What happened?

I was running SCAM on my personal computer (cesm2.1.3, cam_cesm2_1_rel41). I found the bug when I set the `spectralflux` namelist option to true. This should have returned the spectral longwave and shortwave fluxes in addition to broadband fluxes and it didn't. When I added the corresponding fields to the fincl1 namelist option (LD, LU, SD, SU), SCAM stopped with and error because it couldn't find those fields.

### What are the steps to reproduce the bug?

To reproduce this bug, set 'spectralflux' to true and include the fields LD, LU, SD, SU in the fincl1 namelist option.

### What CAM tag were you using?

cam_cesm2_1_rel41

### What machine were you running CAM on?

Personal Computer

### What compiler were you using?

Other (please specify below)

### Path to a case directory, if applicable

_No response_

### Will you be addressing this bug yourself?

No

### Extra info

## Notes
I am running CAM through SCAM on my personal computer and I am not sure what compiler I am using. I have checked the most recent tag of CAM (cam6_3_160) and this bug is still in the CAM code. The buggy files `radconstants.F90` and `radiation.F90` are located at `/cam/src/physics/rrtmg`.

## Fixes
I have produced code that fixes the bug. The files are at `/glade/u/home/glydia/mods/fixedspectral/`.

In the radconstants file, I added the longwave and shortwave band midpoints as variables, as well as functions to retrieve them, similar to existing functions in the file to retrieve the band boundaries. I included the midpoints so that the history file has them and the boundaries as dimensions for the spectral fluxes.

The additions I made to the radconstants file are below:

```
real(r8),parameter :: wavenum_mid(nbndsw) = & ! in cm^-1
(/2906.89_r8, 3605.55_r8, 4312.77_r8, 4893.62_r8, 5627.83_r8, 6881.50_r8, 7873.06_r8, &
10170.67_r8,14338.76_r8,19036.81_r8,25629.08_r8,33196.39_r8,43588.99_r8, 1460.14_r8/)
```
```
real(r8), parameter :: wavenumber3_longwave(nlwbands) = &! Longwave spectral band mid-points (cm-1)
(/ 59.16_r8, 418.33_r8, 561.25_r8, 664.08_r8, 757.63_r8, 896.44_r8, 1028.79_r8, 1128.89_r8, &
1280.70_r8, 1434.29_r8, 1632.18_r8, 1934.94_r8, 2163.33_r8, 2318.94_r8, 2492.79_r8, 2906.89_r8 /)
```
```
public :: get_number_sw_bands, &
get_sw_spectral_boundaries, &
get_sw_spectral_midpts, &
get_lw_spectral_boundaries, &
get_lw_spectral_midpts, &
get_ref_solar_band_irrad, &
get_ref_total_solar_irrad, &
get_solar_band_fraction_irrad
```
```
subroutine get_lw_spectral_midpts(midpts, units)
! provide spectral midpts of each longwave band

real(r8), intent(out) :: midpts(nlwbands)
character(*), intent(in) :: units ! requested units

select case (units)
case ('inv_cm','cm^-1','cm-1')
midpts = wavenumber3_longwave
case('m','meter','meters')
midpts = 1.e-2_r8/wavenumber3_longwave
case('nm','nanometer','nanometers')
midpts = 1.e7_r8/wavenumber3_longwave
case('um','micrometer','micrometers','micron','microns')
midpts = 1.e4_r8/wavenumber3_longwave
case('cm','centimeter','centimeters')
midpts = 1._r8/wavenumber3_longwave
case default
call endrun('get_lw_spectral_midpts: spectral units not acceptable'//units)
end select

end subroutine get_lw_spectral_midpts
```
```
subroutine get_sw_spectral_midpts(midpts, units)
! provide spectral midpts of each shortwave band

real(r8), intent(out) :: midpts(nswbands)
character(*), intent(in) :: units ! requested units

select case (units)
case ('inv_cm','cm^-1','cm-1')
midpts = wavenum_mid
case('m','meter','meters')
midpts = 1.e-2_r8/wavenum_mid
case('nm','nanometer','nanometers')
midpts = 1.e7_r8/wavenum_mid
case('um','micrometer','micrometers','micron','microns')
midpts = 1.e4_r8/wavenum_mid
case('cm','centimeter','centimeters')
midpts = 1._r8/wavenum_mid
case default
call endrun('rad_constants.F90: spectral units not acceptable'//units)
end select

end subroutine get_sw_spectral_midpts
```

In the radiation file, I added the longwave and shortwave boundaries and midpoints as variables and added the LD, LU, SD, SU fields to the history master field list. I added the longwave and shortwave bands to the history file coordinates and added calls to the pbuf in the longwave and shortwave output functions to output the spectral fluxes to the history file.

The additions/changes I made to the radiation file are below:
```
use radconstants, only: nswbands, nlwbands, rrtmg_sw_cloudsim_band, rrtmg_lw_cloudsim_band, &
idx_sw_diag, get_lw_spectral_boundaries, get_sw_spectral_boundaries, &
get_lw_spectral_midpts, get_sw_spectral_midpts
use cam_history_support, only: fillvalue, add_hist_coord
```
```
! LW/SW band variables
real(r8) :: lw_bounds(2,nlwbands)
real(r8) :: sw_bounds(2,nswbands)
real(r8) :: lw_vals(nlwbands)
real(r8) :: sw_vals(nswbands)
```
```
subroutine radiation_register
use cam_history_support, only: add_hist_coord
...
if (spectralflux) then
! Add history coordinates for lw/sw bands
call get_sw_spectral_midpts(sw_vals,'cm^-1')
call get_lw_spectral_midpts(lw_vals,'cm^-1')
call get_sw_spectral_boundaries(sw_bounds(1,:),sw_bounds(2,:),'cm^-1')
call get_lw_spectral_boundaries(lw_bounds(1,:),lw_bounds(2,:),'cm^-1')

call add_hist_coord('sw_band', nswbands, 'Shortwave bands', &
'cm^-1', sw_vals, bounds_name='sw_band_bds', bounds=sw_bounds)
call add_hist_coord('lw_band', nlwbands, 'Longwave bands', &
'cm^-1', lw_vals, bounds_name='lw_band_bds', bounds=lw_bounds)
```
```
subroutine radiation_init(pbuf2d)
...
! Add shortwave radiation fields to history master field list.
do icall = 0, N_DIAG

if (active_calls(icall)) then
...
call addfld('SU'//diag(icall), (/ character(len=8) :: 'ilev','sw_band' /), 'A', 'W/m2', 'Shortwave spectral flux up', &
sampling_seq='rad_lwsw')
call addfld('SD'//diag(icall), (/ character(len=8) :: 'ilev','sw_band' /), 'A', 'W/m2', 'Shortwave spectral flux down', &
sampling_seq='rad_lwsw')
...
! Add longwave radiation fields to history master field list.

do icall = 0, N_DIAG

if (active_calls(icall)) then
...
call addfld('LU'//diag(icall), (/ character(len=8) :: 'ilev','lw_band' /), 'A', 'W/m2', 'Longwave spectral flux up', &
sampling_seq='rad_lwsw')
call addfld('LD'//diag(icall), (/ character(len=8) :: 'ilev','lw_band' /), 'A', 'W/m2', 'Longwave spectral flux down', &
sampling_seq='rad_lwsw')
```
```
subroutine radiation_output_sw(lchnk, ncol, icall, rd, pbuf, cam_out)
real(r8), pointer, dimension(:,:,:) :: su
real(r8), pointer, dimension(:,:,:) :: sd

call pbuf_get_field(pbuf, su_idx, su)
call pbuf_get_field(pbuf, sd_idx, sd)

call outfld('SU'//diag(icall), su, ncol, lchnk)
call outfld('SD'//diag(icall), sd, ncol, lchnk)
```
```
subroutine radiation_output_lw(lchnk, ncol, icall, rd, pbuf, cam_out, freqclr, flntclr)
real(r8), pointer, dimension(:,:,:) :: lu
real(r8), pointer, dimension(:,:,:) :: ld

call pbuf_get_field(pbuf, lu_idx, lu)
call pbuf_get_field(pbuf, ld_idx, ld)

call outfld('LU'//diag(icall), lu, ncol, lchnk)
call outfld('LD'//diag(icall), ld, ncol, lchnk)
```

Contributor guide

No contributing guide indexed for this repository

Assessment

This issue has not been assessed yet.

Get new issues in your inbox

A short digest of beginner-friendly GitHub issues.