cam6_4_156: fixes to prescribed volcanic aerosols

bld/namelist_files/use_cases/1850_cam_mt.xml

@@ -15,15 +15,15 @@
 
 <!-- ozone data :  -->
 <prescribed_ozone_datapath>   'atm/cam/ozone_strataero'  </prescribed_ozone_datapath>
-<prescribed_ozone_file    >   'ozone_strataero_cyclical_WACCM6_L70_CMIP6-piControl.001_y21-50avg_zm_5day_c180802.nc'       </prescribed_ozone_file>
+<prescribed_ozone_file    >   'ozone_strataero_b.e30_alpha07c_cesm.B1850C_MTt4s.ne30_t232_wgx3.251.001_1850climo_zm_5day_c251204.nc'</prescribed_ozone_file>
 <prescribed_ozone_name    >   'O3'  </prescribed_ozone_name>
 <prescribed_ozone_type    >   CYCLICAL   </prescribed_ozone_type>
 <prescribed_ozone_cycle_yr>   1850    </prescribed_ozone_cycle_yr>
 
 <!-- Prescribed stratospheric aerosols : -->
 <prescribed_strataero_cycle_yr>      1850   </prescribed_strataero_cycle_yr>
 <prescribed_strataero_datapath>      'atm/cam/ozone_strataero'   </prescribed_strataero_datapath>
-<prescribed_strataero_file>   'ozone_strataero_cyclical_WACCM6_L70_CMIP6-piControl.001_y21-50avg_zm_5day_c180802.nc' </prescribed_strataero_file>
+<prescribed_strataero_file>   'ozone_strataero_b.e30_alpha07c_cesm.B1850C_MTt4s.ne30_t232_wgx3.251.001_1850climo_zm_5day_c251204.nc'</prescribed_strataero_file>
 <prescribed_strataero_use_chemtrop>  .true.   </prescribed_strataero_use_chemtrop>
 <prescribed_strataero_type>          'CYCLICAL'    </prescribed_strataero_type>
 
@@ -87,6 +87,17 @@
 'SOAE -> 5.1004D0*$INPUTDATA_ROOT/atm/cam/chem/emis/cmip7/ne30/DRES-CMIP-BB4CMIP7-2-0_smoothed_20251102/MTERP_smoothed_input4MIPs_emissions_CMIP_DRES-CMIP-BB4CMIP7-2-0_gn_175001-202112_c20251102.nc'
 </srf_emis_specifier>
 
+<rad_climate>
+  'A:Q:H2O', 'N:O2:O2', 'A:CO2:CO2', 'N:ozone:O3', 'A:N2O:N2O', 'A:CH4:CH4', 'N:CFC11STAR:CFC11', 'A:CFC12:CFC12',
+  'M:mam4_mode1:$INPUTDATA_ROOT/atm/cam/physprops/mam4_mode1_rrtmg_aeronetdust_sig1.6_dgnh.48_c140304.nc',
+  'M:mam4_mode2:$INPUTDATA_ROOT/atm/cam/physprops/mam4_mode2_rrtmg_aitkendust_c141106.nc',
+  'M:mam4_mode3:$INPUTDATA_ROOT/atm/cam/physprops/mam4_mode3_rrtmg_aeronetdust_c141106.nc',
+  'M:mam4_mode4:$INPUTDATA_ROOT/atm/cam/physprops/mam4_mode4_rrtmg_c130628.nc',
+  'N:VOLC_MMR1:$INPUTDATA_ROOT/atm/cam/physprops/volc_camRRTMG_byradius_sigma1.6_mode1_c210211.nc',
+  'N:VOLC_MMR2:$INPUTDATA_ROOT/atm/cam/physprops/volc_camRRTMG_byradius_sigma1.6_mode2_c210211.nc',
+  'N:VOLC_MMR3:$INPUTDATA_ROOT/atm/cam/physprops/volc_camRRTMG_byradius_sigma1.2_mode3_c210211.nc',
+  'N:VOLC_MMR5:$INPUTDATA_ROOT/atm/cam/physprops/volc_camRRTMG_byradius_sigma1.2_mode5_c260228.nc'
+</rad_climate>
 
 <csw_time_type>CYCLICAL</csw_time_type>
 <csw_cycle_yr>1850 </csw_cycle_yr>

bld/namelist_files/use_cases/1850_trop_strat_t4s_cam7.xml

@@ -3,7 +3,7 @@
 <namelist_defaults>
 
 <!-- Initial conditions -->
-<ncdata hgrid="ne30np4" nlev="93">atm/cam/inic/se/1850C_T4S_ne30pg3_spinup01.cam.i.0002-01-01_c241114.nc</ncdata>
+<ncdata hgrid="ne30np4" nlev="93">atm/cam/inic/se/b.e30_alpha07c_cesm.B1850C_MTt4s.ne30_t232_wgx3.251.001.cam.i.0032-01-01_c260304.nc</ncdata>
 
 <!-- Solar data -->
 <solar_irrad_data_file>atm/cam/solar/SolarForcingCMIP7-4.6_18491230-20240101_sumEPP_c20250630.nc</solar_irrad_data_file>

cime_config/testdefs/testmods_dirs/cam/outfrq9s_cam7_presc_volc/user_nl_cam

@@ -18,7 +18,7 @@ rad_climate =
  'N:VOLC_MMR1:$DIN_LOC_ROOT/atm/cam/physprops/volc_camRRTMG_byradius_sigma1.6_mode1_c210211.nc',
  'N:VOLC_MMR2:$DIN_LOC_ROOT/atm/cam/physprops/volc_camRRTMG_byradius_sigma1.6_mode2_c210211.nc',
  'N:VOLC_MMR3:$DIN_LOC_ROOT/atm/cam/physprops/volc_camRRTMG_byradius_sigma1.2_mode3_c210211.nc',
- 'N:VOLC_MMR5:$DIN_LOC_ROOT/atm/cam/physprops/volc_camRRTMG_byradius_sigma1.2_mode3_c210211.nc'
+ 'N:VOLC_MMR5:$DIN_LOC_ROOT/atm/cam/physprops/volc_camRRTMG_byradius_sigma1.2_mode5_c260228.nc'
 
  fincl2 = 'VOLC_MMR1','VOLC_MMR2','VOLC_MMR3','VOLC_MMR5','VOLC_RAD_GEOM1','VOLC_RAD_GEOM2','VOLC_RAD_GEOM3','VOLC_RAD_GEOM5',
-          'VOLC_SAD'
+          'VOLC_SAD'

doc/ChangeLog

@@ -1,5 +1,124 @@
 ===============================================================
 
+Tag name: cam6_4_156
+Originator(s): fvitt, tilmes
+Date: 6 Mar 2026
+One-line Summary: Fix issue with prescribed volcanic aerosols
+Github PR URL: https://github.com/ESCOMP/CAM/pull/1495
+
+Purpose of changes (include the issue number and title text for each relevant GitHub issue):
+
+ Fix issue with 5th mode of the prescribed volcanic sulfate.
+ Issue #1494
+
+ Update IC file for B1850C_MTt4s compset.
+ Issue #1493
+
+Describe any changes made to build system: N/A
+
+Describe any changes made to the namelist: N/A
+
+List any changes to the defaults for the boundary datasets: N/A
+
+Describe any substantial timing or memory changes: N/A
+
+Code reviewed by: cacraigucar
+
+List all files eliminated: N/A
+
+List all files added and what they do: N/A
+
+List all existing files that have been modified, and describe the changes:
+
+M       bld/namelist_files/use_cases/1850_cam_mt.xml
+ - update prescribed ozone and strataero inputs
+ - set rad_climate to include VOLC_MMR5 (prescrided 5th mode of volcanic sulfate)
+   and use corrected physprop file for VOLC_MMR5
+
+M       bld/namelist_files/use_cases/1850_trop_strat_t4s_cam7.xml
+ - update IC file (to include O2)
+
+M       cime_config/testdefs/testmods_dirs/cam/outfrq9s_cam7_presc_volc/user_nl_cam
+ - set rad_climate to include VOLC_MMR5 (prescrided 5th mode of volcanic sulfate)
+ - use corrected physprop file for VOLC_MMR5
+
+M       src/chemistry/modal_aero/aero_model.F90
+ - output diameters of all 5 modes when history_cesm_forcing is TRUE
+
+M       src/chemistry/mozart/mo_chm_diags.F90
+ - output so4 of the 5th mode when history_cesm_forcing is TRUE
+
+M       src/physics/cam/aerosol_optics_cam.F90
+ - use correct radius (5th mode) for VOLC_MMR5
+ - call output_bin_diags for BAM aerosols (more diagnistcis)
+
+M       src/physics/cam/phys_prop.F90
+ - add volcanic_radius5 to the list of valid optics methods
+
+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:
+  FAIL ERI_D_Ln18.ne16pg3_ne16pg3_mt232.FHIST_C4.derecho_intel.cam-outfrq3s_eri
+   - pre-existing failure
+      ERI bug in CICE -- See: https://github.com/ESCOMP/CESM_CICE/issues/34
+
+  DIFF ERP_D_Ln9.ne30pg3_ne30pg3_mt232.F1850C_MTso.derecho_intel.cam-outfrq9s
+   - expected difference due to corrections to prescribed volcanic aerosols
+
+  DIFF ERP_Ln9.f19_f19_mg17.FWsc1850.derecho_intel.cam-outfrq9s
+   - difference only in AODVISstdn (now includes prescribed volcanic sulfate)
+     otherwise bit-for-bit unchanged
+
+  DIFF ERP_Ln9.ne30pg3_ne30pg3_mg17.FHISTC_WAma.derecho_intel.cam-outfrq9s
+  DIFF ERP_Ld3.ne16pg3_ne16pg3_mg17.FHISTC_WAt1ma.derecho_intel.cam-reduced_hist1d
+  DIFF ERS_Ln9.ne30pg3_ne30pg3_mg17.FHISTC_WXma.derecho_intel.cam-outfrq9s_ctem
+  DIFF SMS_D_Ln9.f19_f19_mg17.FWma2000climo.derecho_intel.cam-outfrq9s
+  DIFF SMS_D_Ln9.f19_f19_mg17.FWma2000climo.derecho_intel.cam-outfrq9s_waccm_ma_mam4
+  DIFF SMS_Ld1.f09_f09_mg17.FW2000climo.derecho_intel.cam-outfrq1d
+  DIFF SMS_Ln9.f09_f09_mg17.FW1850.derecho_intel.cam-reduced_hist3s
+  DIFF SMS_C2_D_Ln9.ne16pg3_ne16pg3_mg17.FHISTC_WXma.derecho_intel.cam-outfrq9s
+  DIFF ERS_Ln9.f19_f19_mg17.FXSD.derecho_intel.cam-outfrq9s
+   - default output now includes diameters all 5 modes (for history_cesm_forcing=.true.)
+     otherwise bit-for-bit unchanged
+
+  DIFF ERC_D_Ln9.f19_f19_mg17.QPMOZ.derecho_intel.cam-outfrq3s
+  DIFF ERP_Lh12.f19_f19_mg17.FW4madSD.derecho_intel.cam-outfrq3h
+   - differences in fill patterns of AOD diags, otherwise bit-for-bit unchanged
+
+derecho/nvhpc/aux_cam: All PASS
+
+izumi/nag/aux_cam:
+  DIFF ERC_D_Ln9.ne3pg3_ne3pg3_mt232.FHISTC_LTso.izumi_nag.cam-outfrq9s_nochem  
+  DIFF SMS_D_Ln6.ne5_ne5_mg37.QPWmaC4.izumi_nag.cam-outfrq3s_physgrid_tem
+   - differences in fill patterns of AOD diags, otherwise bit-for-bit unchanged
+
+izumi/gnu/aux_cam:
+  DIFF ERC_D_Ln9.f10_f10_mg37.QPC4.izumi_gnu.cam-outfrq3s_diags
+   - differences in fill patterns of AOD diags and other aerosol optics now includes
+     BAM aerosols, otherwise bit-for-bit unchanged
+
+  DIFF ERC_D_Ln9.f10_f10_mt232.FHIST_C5.izumi_gnu.cam-rad_diag
+  DIFF ERC_D_Ln9.ne3pg3_ne3pg3_mt232.F1850_C4.izumi_gnu.cam-co2rmp
+  DIFF ERC_D_Ln9.ne3pg3_ne3pg3_mt232.FHIST.izumi_gnu.cam-nochem_clubbmf
+  DIFF SMS_D_Ln3.f10_f10_mg37.QPMOZ.izumi_gnu.cam-outfrq3s_chemproc
+   - differences in fill patterns of AOD diags, otherwise bit-for-bit unchanged
+
+  DIFF SMS_D_Ln9.f10_f10_mg37.FWmaHIST.izumi_gnu.cam-outfrq9s_mee_fluxes
+  DIFF SMS_D_Ln9.f10_f10_mg37.QPWmaC4.izumi_gnu.cam-outfrq9s_apmee
+  DIFF SMS_D_Ln9.f19_f19_mt232.FW4madSD.izumi_gnu.cam-outfrq9s
+   - default output now includes diameters all 5 modes (for history_cesm_forcing=.true.)
+     otherwise bit-for-bit unchanged
+
+Summarize any changes to answers:
+
+ larger than roundoff for F1850C_MTso compset, otherwise bit-for-bit unchanged
+
+===============================================================
+===============================================================
+
 Tag name: cam6_4_155
 Originator(s): cacraig, jimmielin
 Date: Feb 24, 2026

src/chemistry/modal_aero/aero_model.F90

@@ -444,7 +444,7 @@ subroutine aero_model_init( pbuf2d )
           call add_default( dgnum_name(n), 1, ' ' )
           call add_default( dgnumwet_name(n), 1, ' ' )
        endif
-       if ( history_cesm_forcing .and. n<4 ) then
+       if ( history_cesm_forcing ) then
           call add_default( dgnumwet_name(n), 8, ' ' )
        endif
 

src/chemistry/mozart/mo_chm_diags.F90

@@ -28,7 +28,7 @@ module mo_chm_diags
   integer :: id_hf,id_f,id_cof2,id_cofcl,id_ch3br
   integer :: id_br,id_bro,id_hbr,id_hobr,id_ch4,id_h2o,id_h2
   integer :: id_o,id_o2,id_h, id_h2o2, id_n2o
-  integer :: id_co2,id_o3,id_oh,id_ho2,id_so4_a1,id_so4_a2,id_so4_a3
+  integer :: id_co2,id_o3,id_oh,id_ho2,id_so4_a1,id_so4_a2,id_so4_a3,id_so4_a5
   integer :: id_num_a2,id_num_a3,id_dst_a3,id_ncl_a3
   integer :: id_ndep,id_nhdep
   integer :: id_clno2,id_brno2,id_br2
@@ -160,6 +160,7 @@ subroutine chm_diags_inti
     id_so4_a1  = get_spc_ndx( 'so4_a1' )
     id_so4_a2  = get_spc_ndx( 'so4_a2' )
     id_so4_a3  = get_spc_ndx( 'so4_a3' )
+    id_so4_a5  = get_spc_ndx( 'so4_a5' )
     id_num_a2  = get_spc_ndx( 'num_a2' )
     id_num_a3  = get_spc_ndx( 'num_a3' )
     id_dst_a3  = get_spc_ndx( 'dst_a3' )
@@ -568,6 +569,7 @@ subroutine chm_diags_inti
           if (m==id_so4_a1) call add_default( spc_name, 8, ' ')
           if (m==id_so4_a2) call add_default( spc_name, 8, ' ')
           if (m==id_so4_a3) call add_default( spc_name, 8, ' ')
+          if (m==id_so4_a5) call add_default( spc_name, 8, ' ')
 
           if (m==id_num_a2) call add_default( spc_name, 8, ' ')
           if (m==id_num_a3) call add_default( spc_name, 8, ' ')

src/physics/cam/aerosol_optics_cam.F90

@@ -178,6 +178,9 @@ subroutine aerosol_optics_cam_init
     modal_active = nmodes>0
     carma_active = nbins>0
     bulk_active = nbulk_aerosols>0
+    if (masterproc) then
+       write(iulog,*) prefix,'nmodes,nbins,nbulk_aerosols: ',nmodes,nbins,nbulk_aerosols
+    end if
 
     ! count aerosol models
     if (modal_active) then
@@ -863,15 +866,19 @@ subroutine aerosol_optics_cam_sw(list_idx, state, pbuf, nnite, idxnite, tauxar,
           case('nonhygro', 'insoluble')
              aero_optics=>insoluble_aerosol_optics(aeroprops, aerostate, list_idx, ibin)
 
-          case('volcanic_radius','volcanic_radius1','volcanic_radius2','volcanic_radius3')
+          case('volcanic_radius','volcanic_radius1','volcanic_radius2','volcanic_radius3','volcanic_radius5')
+
+             ! construct name of radius physics buffer field
              pbuf_fld = 'VOLC_RAD_GEOM '
              if (len_trim(opticstype)>15) then
                 pbuf_fld = trim(pbuf_fld)//opticstype(16:16)
              endif
+
              ! get microphysical properties for volcanic aerosols
              idx = pbuf_get_index(pbuf_fld)
-             call pbuf_get_field(pbuf, idx, geometric_radius )
+             call pbuf_get_field(pbuf, idx, geometric_radius)
 
+             ! construct aerosol optics object
              aero_optics=>volcrad_aerosol_optics(aeroprops, aerostate, list_idx, &
                   ibin, ncol, pver, geometric_radius(:ncol,:))
 
@@ -904,9 +911,7 @@ subroutine aerosol_optics_cam_sw(list_idx, state, pbuf, nnite, idxnite, tauxar,
                       ! ref: Fig. 1d of Jasper F. Kok et al. (2017),
                       ! Smaller desert dust cooling effect estimated from analysis of dust size and abundance
 
-                      if (.not.aeroprops%model_is('BAM')) then
-                         call update_diags( is_coarse_dust=coarse_dust_mode )  ! dopaer is updated in update_diags.
-                      end if
+                      call update_diags( is_coarse_dust=coarse_dust_mode )  ! dopaer is updated in update_diags.
 
                       ! dmleung: update_diags updated dopaer(icol) as a diagnostic.
                       ! Aerosol optical and radiative properties are subsequently modified given dopaer update in update_diags.
@@ -941,9 +946,9 @@ subroutine aerosol_optics_cam_sw(list_idx, state, pbuf, nnite, idxnite, tauxar,
              bam_cnt = bam_cnt+1
              call aer_vis_diag_out(lchnk, ncol, nnite, idxnite, bam_cnt, taubam, &
                   list_idx, troplev)
-          else
-             call output_bin_diags()
-          end if
+          endif
+
+          call output_bin_diags()
 
        end do binloop
     end do aeromodel
@@ -1018,33 +1023,45 @@ subroutine update_diags( is_coarse_dust )
             case('dust')
                dustvol(icol) = vol(icol)
                burdendust(icol) = burdendust(icol) + specmmr(icol,ilev)*mass(icol,ilev)
-               scatdust(icol) = vol(icol) * specrefindex(iwav)%re
-               absdust(icol)  =-vol(icol) * specrefindex(iwav)%im
+               if (associated(specrefindex)) then
+                  scatdust(icol) = vol(icol) * specrefindex(iwav)%re
+                  absdust(icol)  =-vol(icol) * specrefindex(iwav)%im
+               end if
                hygrodust(icol)= vol(icol)*hygro_aer
             case('black-c')
                burdenbc(icol) = burdenbc(icol) + specmmr(icol,ilev)*mass(icol,ilev)
-               scatbc(icol) = vol(icol) * specrefindex(iwav)%re
-               absbc(icol)  =-vol(icol) * specrefindex(iwav)%im
+               if (associated(specrefindex)) then
+                  scatbc(icol) = vol(icol) * specrefindex(iwav)%re
+                  absbc(icol)  =-vol(icol) * specrefindex(iwav)%im
+               end if
                hygrobc(icol)= vol(icol)*hygro_aer
             case('sulfate')
                burdenso4(icol) = burdenso4(icol) + specmmr(icol,ilev)*mass(icol,ilev)
-               scatsulf(icol) = vol(icol) * specrefindex(iwav)%re
-               abssulf(icol)  =-vol(icol) * specrefindex(iwav)%im
+               if (associated(specrefindex)) then
+                  scatsulf(icol) = vol(icol) * specrefindex(iwav)%re
+                  abssulf(icol)  =-vol(icol) * specrefindex(iwav)%im
+               end if
                hygrosulf(icol)= vol(icol)*hygro_aer
             case('p-organic')
                burdenpom(icol) = burdenpom(icol) + specmmr(icol,ilev)*mass(icol,ilev)
-               scatpom(icol) = vol(icol) * specrefindex(iwav)%re
-               abspom(icol)  =-vol(icol) * specrefindex(iwav)%im
+               if (associated(specrefindex)) then
+                  scatpom(icol) = vol(icol) * specrefindex(iwav)%re
+                  abspom(icol)  =-vol(icol) * specrefindex(iwav)%im
+               end if
                hygropom(icol)= vol(icol)*hygro_aer
             case('s-organic')
                burdensoa(icol) = burdensoa(icol) + specmmr(icol,ilev)*mass(icol,ilev)
-               scatsoa(icol) = vol(icol) * specrefindex(iwav)%re
-               abssoa(icol) = -vol(icol) * specrefindex(iwav)%im
+               if (associated(specrefindex)) then
+                  scatsoa(icol) = vol(icol) * specrefindex(iwav)%re
+                  abssoa(icol) = -vol(icol) * specrefindex(iwav)%im
+               end if
                hygrosoa(icol)= vol(icol)*hygro_aer
             case('seasalt')
                burdenseasalt(icol) = burdenseasalt(icol) + specmmr(icol,ilev)*mass(icol,ilev)
-               scatsslt(icol) = vol(icol) * specrefindex(iwav)%re
-               abssslt(icol) = -vol(icol) * specrefindex(iwav)%im
+               if (associated(specrefindex)) then
+                  scatsslt(icol) = vol(icol) * specrefindex(iwav)%re
+                  abssslt(icol) = -vol(icol) * specrefindex(iwav)%im
+               end if
                hygrosslt(icol)= vol(icol)*hygro_aer
             end select
          end do
@@ -1432,15 +1449,19 @@ subroutine aerosol_optics_cam_lw(list_idx, state, pbuf, tauxar)
           case('nonhygro', 'insoluble')
              aero_optics=>insoluble_aerosol_optics(aeroprops, aerostate, list_idx, ibin)
 
-          case('volcanic_radius','volcanic_radius1','volcanic_radius2','volcanic_radius3')
+          case('volcanic_radius','volcanic_radius1','volcanic_radius2','volcanic_radius3','volcanic_radius5')
+
+             ! construct name of radius physics buffer field
              pbuf_fld = 'VOLC_RAD_GEOM '
              if (len_trim(opticstype)>15) then
                 pbuf_fld = trim(pbuf_fld)//opticstype(16:16)
              endif
+
              ! get microphysical properties for volcanic aerosols
              idx = pbuf_get_index(pbuf_fld)
-             call pbuf_get_field(pbuf, idx, geometric_radius )
+             call pbuf_get_field(pbuf, idx, geometric_radius)
 
+             ! construct aerosol optics object
              aero_optics=>volcrad_aerosol_optics(aeroprops, aerostate, list_idx, &
                   ibin, ncol, pver, geometric_radius(:ncol,:))
 

src/physics/cam/phys_prop.F90

@@ -537,7 +537,7 @@ subroutine aerosol_optics_init(phys_prop, nc_id)
    case ('insoluble')
       call insoluble_optics_init(phys_prop, nc_id)
 
-   case ('volcanic_radius','volcanic_radius1','volcanic_radius2','volcanic_radius3')
+   case ('volcanic_radius','volcanic_radius1','volcanic_radius2','volcanic_radius3','volcanic_radius5')
       call volcanic_radius_optics_init(phys_prop, nc_id)
 
    case ('volcanic')