SOCOLv3 bug fixes
Collection of bug fixes for various SOCOLv3 versions incl. SOCOL-AER after Nov 2020:
SOCOL/SOCOL-AER: bug fix in interactive wet deposition scheme (Mar 2023)
A wrong cloud cover variable is used in the convective scavenging of gas-phase species:
...
!!$ mcov(jl) = zoldcov(i)
mcov(jl) = zwashfrac(i)
...
New code: mo_socol_wetdep_main.f90
SOCOL/SOCOL-AER: bugs in interactive wet deposition scheme (Nov 2022)
SOCOL-AER did not run on 48 CPUs anymore
By change I realized that SOCOL-AER does not run on 48 CPUs (not sure if this problem is new and related to new compiler revisions etc. or if nobody run SOCOL-AER on 48 CPUs in the past...)
Problem: shape of argument SVI in call to socol_wet_deposition (src/socol_mezon.f90), replace SVI(1:kproma,:,:,krow) [or SVI(1:kproma,1:klev,1:n_aer_bin,krow) ] by SVI(:,:,:,krow)
New code:
...
if (linteractive_wetdep) CALL socol_wet_deposition(krow, kproma, pxtm1, pxtte, t, p, ph, paclc, &
plwc, piwc, pfrain_no, pfsnow_no, prate_r,prate_s, pr_cover, & !large scale properties
pmflxr, pmflxs, cv_precnew, cv_snownew, zmfu, kconbot, & !convective properties
SVI(:,:,:,krow))
...
Problems in calculation of convective cloud cover
For runs with interactive wet deposition, the model crashed from time to time with the following error:
Numerical Conversion Not Representable
Reason was a bug in the calculation of the convective cloud cover in mo_socol_wetdep_main.f90:
New code: mo_socol_wetdep_main.f90
SOCOL-AER: Reading volcanic SO2 emissions T${RES}_so2_emiss_volcYYYYMM (Nov 2022)
Error message src/socol_boundary_conditions.f90, array so2_erupt: Index 31 of dim 1 beyond upper bound of 29
Problem occurs at the beginning of February, as the model run is resumed on 31 Jan 23:30:00
Workaround: Add if-statement to src/socol_boundary_conditions.f90, subroutine set_socol_bcond_fluxes:
...
fac = 1._dp/(atomweight*amso2) ! [kg/m^2/s] -> [molec/m^2/s]
call get_date_components(current_date,yr,mo,dy)
if (dy .gt. size(so2_erupt,1)) goto 100
...
ENDDO
100 continue
! SO2 aircraft emissions:
...
Ozone decoupling: lo3_coupl = .FALSE., but lsrb_lya_heating = .TRUE. (Jun 2022)
The previous bug fix allowed to use an ozone climatology provided on pressure levels only. The current fix allows now to use all types of possible ozone climatologies for the calculation of the mesospheric heating rates.
New code [modules/mo_socol_sun.f90, subroutine srb_lya_heating]:
USE mo_radiation, ONLY: io3
USE mo_o3_lwb, ONLY: o3_lwb
USE mo_o3clim, ONLY: o3clim, o3clim_rawlev
...
REAL(dp) :: zdp(kbdim,klev)
...
! g/g -> mol/mol
oz: SELECT CASE (io3)
CASE (0)
vmr_o3(1:kproma,:) = EPSILON(1._dp)
CASE (2)
! layer pressure thickness
zdp(1:kproma,:)=paphm1(1:kproma,2:klev+1)-paphm1(1:kproma,1:klev)
vmr_o3(1:kproma,:) = o3_lwb(krow,zdp,paphm1)*amd/amo3
CASE (3)
vmr_o3(1:kproma,:) = o3clim(krow,kproma,kbdim,klev,paphm1,papm1)*amd/amo3
CASE (4)
vmr_o3(1:kproma,:) = o3clim(krow,kproma,kbdim,klev,paphm1,papm1)*amd/amo3
CASE (5)
vmr_o3(1:kproma,:) = o3clim_rawlev(krow,kproma,kbdim,klev)*amd/amo3
CASE (10) ! coupling with CTM
vmr_o3(:,:)=xtm1(:,:,idt_o3,krow)
CASE default
CALL finish('srb_lya_heating','o3: this "io3" is not supported')
END SELECT oz
Ozone decoupling: lo3_coupl = .FALSE., but lsrb_lya_heating = .TRUE. (Jul 2021)
The extra heating rates due to absorption in the Schumann Runge bands (175-200 nm) and by the Lyman alpha line (121.6 nm) depend on the ozone field. In case the coupling between chemistry and radiation via ozone is switched off (lo3_coupl = .FALSE.), the ozone climatology used for the radiative calculations has to be used in the
subroutine srb_lya_heating as well! For that, an additional if-statement has to added at the beginning of the subroutine. Details see below.
Old code [modules/mo_socol_sun.f90, subroutine srb_lya_heating]:
There are different code versions available, some include already a first bug fix.
Either:
vmr_o3(:,:)=xtm1(:,:,idt_o3,krow)
or:
USE mo_o3clim , ONLY: o3clim
...
if (lchem) then
vmr_o3(:,:)=xtm1(:,:,idt_o3,krow)
else
! g/g -> mol/mol
vmr_o3(1:kproma,:) = o3clim(krow,kproma,kbdim,klev,paphm1,papm1)*amd/amo3
end if
New code [modules/mo_socol_sun.f90, subroutine srb_lya_heating]:
USE mo_o3clim , ONLY: o3clim
...
if (lchem) then
if (lo3_coupl) then
vmr_o3(:,:)=xtm1(:,:,idt_o3,krow)
else
! g/g -> mol/mol
vmr_o3(1:kproma,:) = o3clim(krow,kproma,kbdim,klev,paphm1,papm1)*amd/amo3
end if
else
! g/g -> mol/mol
vmr_o3(1:kproma,:) = o3clim(krow,kproma,kbdim,klev,paphm1,papm1)*amd/amo3
end if
Problems with LSOLARVAR = .FALSE. (Dec 2020)
The SOCOL namelist parameter
LSOLARVAR determines, whether the solar input (irradiance, additional heating by Schumann-Runge-, Hartley-, Huggins- and Lyman-alpha-radiation, NO photolysis) varies with time or not. In the latter case monthly mean values averaged over two consecutive solar cycles (1977-1998) are used. These values are set in
mo_socol_namelist.f90.
For CCMI-1, the parameterization for extraheating as well as the solar input files have been updated, but not the settings in
mo_socol_namelist.f90 for
lsolarvar = .false., at least not in all model versions. In some model versions, the settings had been updated, but the values for Schumann-Runge- and Lyman-alpha-bands have been swapped.
Furthermore, the values for the two Huggings-parameters have to be converted from W/m2 to W/m2/angstrom. This was missing for the case
lsolarvar = .false..
Modified source code
And here the updated constant values for
mo_socol_namelist.f90:
REAL(dp):: sun_srb_const(12) = (/ &
1.51E-01_dp, 1.50E-01_dp, 1.48E-01_dp, 1.46E-01_dp, &
1.44E-01_dp, 1.42E-01_dp, 1.41E-01_dp, 1.42E-01_dp, &
1.44E-01_dp, 1.46E-01_dp, 1.48E-01_dp, 1.50E-01_dp /)
! Monthly values of Schumann-Runge
! parameter if lsolarvar=.FALSE.
! (Default: mean over 1977-1998)
REAL(dp):: sun_lya_const(12) = (/ &
7.03E-03_dp, 7.04E-03_dp, 6.92E-03_dp, 6.80E-03_dp, &
6.66E-03_dp, 6.57E-03_dp, 6.54E-03_dp, 6.60E-03_dp, &
6.68E-03_dp, 6.82E-03_dp, 6.91E-03_dp, 7.07E-03_dp /)
! Monthly values of Lyman alpha
! parameter if lsolarvar=.FALSE.
! (Default: mean over 1977-1998)
REAL(dp):: sun_har_const(12) = (/ &
3.28E-01_dp, 3.08E-01_dp, 2.52E-01_dp, 1.77E-01_dp, &
1.02E-01_dp, 4.68E-02_dp, 2.71E-02_dp, 4.77E-02_dp, &
1.02E-01_dp, 1.77E-01_dp, 2.52E-01_dp, 3.09E-01_dp /)
! Monthly values of Hartley
! parameter if lsolarvar=.FALSE.
! (Default: mean over 1977-1998)
REAL(dp):: sun_hug1_const(12) = (/ &
7.15E-01_dp, 6.67E-01_dp, 5.41E-01_dp, 3.69E-01_dp, &
1.97E-01_dp, 7.16E-02_dp, 2.68E-02_dp, 7.30E-02_dp, &
1.97E-01_dp, 3.68E-01_dp, 5.40E-01_dp, 6.69E-01_dp/)
! Monthly values of Huggins 1
! parameter if lsolarvar=.FALSE.
! (Default: mean over 1977-1998)
REAL(dp):: sun_hug2_const(12) = (/ &
3.11E+00_dp, 2.90E+00_dp, 2.34E+00_dp, 1.57E+00_dp, &
8.06E-01_dp, 2.49E-01_dp, 4.92E-02_dp, 2.53E-01_dp, &
8.06E-01_dp, 1.57E+00_dp, 2.33E+00_dp, 2.90E+00_dp/)
! Monthly values of Huggins 2
! parameter if lsolarvar=.FALSE.
! (Default: mean over 1977-1998)
REAL(dp):: sun_nopho_const(12) = (/ &
4.67E-06_dp, 4.66E-06_dp, 4.59E-06_dp, 4.52E-06_dp, &
4.44E-06_dp, 4.39E-06_dp, 4.37E-06_dp, 4.39E-06_dp, &
4.45E-06_dp, 4.52E-06_dp, 4.59E-06_dp, 4.66E-06_dp /)
! Monthly values of NO-photolysis
! parameter if lsolarvar=.FALSE.
! (Default: mean over 1977-1998)
SOCOL-AER: aerosol-radiation coupling (Dec 2020)
Due to an erroneous do loop over vertical levels, the extinction coefficients for AER aerosols had not been set to zero for the lower most 6 model levels. Furthermore, the extinctions had not been converted to optical depth. That leads to a double-counting of AOD from sulfate aerosols in the lowermost 6 levels:
