Evaluation of the ECHAM family radiation codes performance in the representation of the solar signal

Solar radiation is the main source of energy for the Earth’s atmosphere and in many respects defines its composition, photochemistry, temperature profile and dynamics. The magnitude of the solar irradiance variability strongly depends on the wavelength, making difficult its representation in climate models. Due to some deficiencies in the applied radiation codes, several models fail to show a clear response in middle stratospheric heating rates to solar spectral irradiance variability; therefore, it is important to evaluate model performance in this respect before doing multiple runs. In this work we evaluate the performance of three generations of ECHAM (4, 5 and 6) solar radiation schemes by a comparison with the reference high-resolution libRadtran code. We found that all original ECHAM radiation codes miss almost all solar signals in the heating rates in the mesosphere. In the stratosphere the two-band ECHAM4 code (E4) has an almost negligible radiative response to solar irradiance changes and the six-band ECHAM5 code (E5c) reproduces only about half of the reference signal, while representation in the ECHAM6 code (E6) is better – it misses a maximum of about 15 % in the upper stratosphere. On the basis of the comparison results we suggest necessary improvements to the ECHAM family codes by the inclusion of available parameterizations of the heating rate due to absorption by oxygen (O2) and ozone (O3). Improvement is presented for E5c and E6, and both codes, with the introduced parameterizations, represent the heating rate response to the spectral solar irradiance variability simulated with libRadtran much better without a substantial increase in computer time. The suggested parameterizations are recommended to be applied in the middle-atmosphere version of the ECHAM-5 and 6 models for the study of the solar irradiance influence on climate.


Introduction
Although solar ultraviolet radiation (SUV) comprises only a couple of percent of the total solar irradiance (TSI), it plays a crucial role, largely defining the structure of the middle atmosphere.While the radiation in the visible (VIS) and infrared spectral ranges of the solar spectrum propagates through the atmosphere without significant absorption, almost all solar ultraviolet irradiance below 300 nm is absorbed by ozone and oxygen above the troposphere and represents the main source of energy in these regions.Furthermore, the SUV is strongly modulated by the solar rotational and 11-year solar cycles.Whereas the variability of TSI during an 11-year solar activity cycle is around 0.1 %, SUV variations can be more than 10 times higher.Moreover, recent measurements by the SORCE (SOlar Radiation and Climate Experiment) suggest an SUV variability significantly higher than all previous estimates (Ermolli et al. (2013) and references therein).
Changes in SUV irradiance lead to significant ozone, temperature and zonal wind responses in the stratosphere and mesosphere, which has been shown in many modelling and observation data analysis studies (Hood and Soukharev, 2012;Austin et al., 2008;Gray et al., 2010;Haigh et al., 2010;Shapiro et al., 2013).The SUV is not considered as a direct radiative forcing for troposphere and surface, since it does not reach these altitudes, but there are indirect effects of solar irradiance variability, which are communicated downward in the so-called "top-down" mechanism: the modulation of stratospheric temperatures leads to dynamical feedbacks by affecting the Brewer-Dobson circulation and hence the stratosphere-troposphere exchange, resulting in decadal T. Sukhodolov et al.: Evaluation of the ECHAM family radiation codes performance climate changes in the lower atmosphere (Solomon et al., 2007;Gray et al., 2010;Ermolli et al., 2013).
A comprehensive study of the entangled possible effects of solar variability requires chemistry-climate models (CCMs), the main instruments which are capable of taking into account many atmospheric chemical, dynamical and temperature feedbacks.To this end, CCMs should contain a correct representation of the radiative transfer in the atmosphere.Accurate codes for radiative transfer solution exist, e.g.libRadtran (Mayer and Kylling, 2005), but they are too computationally expensive to be commonly used in global models.Therefore, different parameterizations have been designed to provide a compromise between accuracy and efficiency.Since most CCMs arise from global circulation models (GCMs), which are primarily tropospheric models, their radiation schemes carefully treat the longwave part of the spectrum, whereas the representation of the solar irradiance is coarse, approximating the entire UV/VIS spectral range by one or two spectral bands and not considering wavelengths shorter then ∼ 250 nm.The evaluation of the radiation codes performed in the framework of the Stratospheric Processes And their Role in Climate (SPARC) Chemistry-Climate Model Validation (CCMVal-2) project (Forster et al., 2011;SPARC CCMval, 2010) has shown that only a few CCM radiation codes are capable of reproducing the magnitude and vertical profile of heating rate differences between solar minimum and maximum, which in turn directly depend on the treatment of the spectral resolution in the codes.
As was pointed out by Forster et al. (2011), a good representation of the solar signal can be obtained by increasing the number of spectral intervals.However, such an approach implies an increase in computational costs, which is a sensitive issue for already numerically expensive global CCMs.Nissen et al. (2007) replaced the 4-band scheme of Fouquart and Bonnel (1980) above 70 hPa by a 49-band parameterization Freie Universität Berlin radiation scheme based on the Beer-Lambert law and allowing a good agreement with a reference model.They showed that the reduction of the FUBrad resolution to six bands results in a 20 % loss of the solarvariability-induced changes in heating rates.Another way is to apply parameterization only for the missed extra heating due to solar UV enhancement.It has been already used in Middle Atmosphere version of ECHAM 4 (MAECHAM-4) (Egorova et al., 2004) and Canadian Middle Atmosphere Model (CMAM) (Fomichev et al., 2004;Semeniuk et al., 2011) in order to parameterize the solar signal in missing and/or underrepresented spectral intervals.These parameterizations are also based on the Beer-Lambert law (Strobel, 1978;Nicolet, 1985;Zhu, 1994) but apply a smaller number of spectral bands (four to eight) compared to Nissen et al. (2007) and still demonstrate good accuracy and efficiency.The most recent way of obtaining satisfying results even with a relatively small number of spectral intervals is to use a completely different approach of incorporating non-grey gaseous absorption based on the so-called "correlated k-distribution" method (e.g.Fu and Liou, 1992).This method exploits the cumulative probability of the absorption coefficient in a spectral interval to replace wave number as an independent variable.Such a code is a part of ECHAM6 (Stevens et al., 2013), but its performance in respect to solar UV influence has not been checked, which limits its application for solar-climate studies.
In this paper we evaluate the performance of the ECHAM family radiation codes in reproducing the heating rate response to SUV variability through the detailed comparison with the reference libRadtran code.We demonstrate the weaknesses of the ECHAM family solar radiation codes and suggest possible ways of improving their performance.(Simmons et al., 1989).It covered only the lower part of the atmosphere up to the 25 hPa level.Therefore, its solar radiation scheme (Fouquart and Bonnel, 1980), inherited by ECHAM, was quite crude with respect to the shortwave part of spectrum: it had only one band covering the UV/VIS parts of the solar spectrum (250-680 nm) and one band covering near infrared (NIR), considered absorption by O 3 and H 2 O, and used TSI as input, i.e. change in the TSI was equally distributed among all spectral bands, and high shortwave variability was missed.This scheme (E4 hereafter) was used up to ECHAM4 until the NIR part of this scheme was extended to three bands (Table 1) in ECHAM5 (E5 hereafter).The weakness of both these versions in representing the solar signal was demonstrated several times in stand-alone form (Solomon et al., 2007;Forster et al., 2011) and within CCMs (Egorova et al., 2004;Cagnazzo et al., 2007;Nissen et al., 2007): basically, they have an almost negligible radiative response to solar irradiance changes due to the lack of wavelength dependence within the one broad UV/VIS band.E5 was further upgraded in Cagnazzo et al. (2007) by extending the number of spectral intervals from one in UV/VIS to three, with two covering the UV range and switching to spectral solar irradiance (SSI) as input (E5c hereafter).This allowed reproducing about half of the reference heating rate differences (Forster et al., 2011).However, this scheme still does not contain any O 2 absorption.One of the main improvements of ECHAM6 compared to previous versions was the adaptation of another solar radiation scheme, namely the rapid radiation transfer model optimized for general circulation modelling studies (E6 hereafter) (Stevens et al., 2013).This scheme is ∼ 10 times faster than previous schemes, it uses the correlated k distribution method, and solar irradiance is calculated over a prescribed number of pseudo wavelength or g-points regarding the absorbing features of certain wavelengths.Quadrature is performed over 112 g-points in the shortwave part of the spectrum, which then are grouped into 14 bands with 3 bands in UV (Table 1).The model has three UV spectral bands and considers oxygen absorption.However, the lowest wavelength boundary is 200 nm (Iacono et al., 2008) so that important features such as the solar Lyman-α (121.6 nm) line (LYA) and part of the Schumann-Runge oxygen absorption bands (SRB) are not taken into account.

Validation
To demonstrate the capabilities of the original codes we performed calculations with stand-alone versions of E4, E5c and E6 for the tropical standard atmosphere, with a solar zenith angle equal to 10 • and for solar minimum and maximum conditions.We have not analysed E5 separately since it has the same single UV/VIS band as E4.To validate the original schemes, we compare all our calculations to the reference code libRadtran (Mayer and Kylling, 2005), which has shown high accuracy in a number of intercomparison studies.
For the 120-440 nm range libRadtran considers more than 16 000 wavelengths, resolving in detail all relevant spectral features.Figure 1 shows the input information that we used to simulate solar variability: the solar irradiance changes, i.e. the relative difference between the irradiances during solar maximum and minimum conditions, and resulting solarinduced ozone changes.The irradiance spectrum for solar minimum and maximum conditions was calculated with the calculated with the COSI (COde for Solar Irradiance) code (Shapiro et al., 2010), following the approach presented in Shapiro et al. (2011).The solar minimum and maximum conditions correspond to sunspot numbers equal 0 and 120, respectively.We note that the spectral profile of the solar irradiance variability on the 11-year timescale yielded by the approach presented in Shapiro et al. (2011) agrees well with other reconstructions (Ermolli et al., 2013).Figure 1 shows that the solar irradiance variability is a very sophisticated function of wavelength.The ozone changes during an 11- year solar activity cycle were estimated from a composite of observational data (Soukharev and Hood, 2006;Austin et al., 2008;SPARC CCMVal, 2010).
Figure 2 illustrates the heating rates calculated by original E4, E5c and E6 schemes and by libRadtran for solar minimum conditions and heating rate differences between solar maximum and minimum caused only by the solar irradiance changes.In terms of absolute values, E5c and E6 overestimate heating rates compared to libRadtran by up to 2 and 3.5 K day −1 , respectively.This overestimation arises from 250-440 nm (E5c) and 263-345 nm (E6) model bands, i.e. from Hartley (HAR) and Huggins (HUG) ozone absorption bands.In the mesosphere, E5c underestimates absolute values by up to 5 K day −1 since it does not take into account any oxygen absorption.E6 considers absorption by oxygen and shows adequate absolute values in the mesosphere although its lowest wavelength bound is 200 nm.
The E4 and E5c absolute values comparison is also shown in Fig. 3  libRadtran to 690 nm.For this analysis we have calculated the daily averaged shortwave heating rates for the tropical atmosphere following the same approach as in Cagnazzo et al. (2007).Although E4 starts from 250 nm, it shows an almost perfect agreement with libRadtran with a slight overestimation around 40 km, which is fully consistent with Nissen et al. (2007).The fact that E5c shows higher heating rates than E4 is consistent with Cagnazzo et al. (2007) and Forster et al. (2011); however, the value of this difference is higher in Cagnazzo et al. (2007), and libRadtran results are positioned between E4 and E5c in Forster et al. (2011).In these two comparisons, NIR was also included, producing additional heating (Fomichev, 2009) and additional distinctions between the models that can probably explain this inconsistency.Cagnazzo et al. (2007) also used another reference model that was more consistent with E5c in the upper stratosphere; this means that deviations found using libRadtran are comparable to the uncertainty range between high-resolution models.
In terms of heating rate responses to SUV changes (Fig. 2), all schemes highly underestimate the solar signal in the mesosphere.At these altitudes, heating rates are significantly defined by oxygen absorption in a highly variable LYA and SRB, which is completely missed in E4 and E5c and only covered slightly in E6.In the upper stratosphere, E5c and E6 first bands covering the Herzberg continuum and part of HAR are reproduced well.However, the contribution from the second bands containing HAR and HUG is noticeably underestimated, causing the main deviation from the reference model to result in a total maximum of 45 and 15 % deviation at 49 km for E5c and E6, respectively.E4 is able to reproduce only 10 % of the signal at 49 km.The results of E4 and E5c are in agreement with previous comparison studies (Forster et al., 2011;SPARC CCMval, 2010).The underestimation of all schemes in HAR-HUG bands can be explained by a high spectral inhomogeneity of the solar irradiance variability in these regions (see Fig. 1), which is smoothed in integrated fluxes.Since the main disagreement appears in this wavelength region, it should be paid more attention in the future evolution of heating rate parameterizations.If the higher UV variability suggested by SORCE (Ermolli et al., 2013) is correct, the absolute values of the missed solar signal in heating rates would, correspondingly, be higher, resulting in more discrepancy in all feedbacks related to solar irradiance changes.

Implementation of the parameterizations
We do not consider E4 further because its upgraded version was already discussed in Egorova et al. (2004) and Forster et al. (2011) and currently it is not used as widely as E5c and E6 anymore.To improve the representation of the solar signal, we implemented the parameterizations of the heating rates in the spectral regions, where we found problems in the previous section.All parameterizations use the same approach based on Strobel (1978), deriving heating rates H from the  1).
From Zhu (1994), we used the following for SRB: where σ srb = 2.07 × 10 −24 m 2 , x srb = N 2,top /N 2 0.3 σ srb and y srb = 0.0152.And for HAR and HUG from Zhu (1994) we used where M = 0.01273 Å −1 , λ short , λ long = (2805, 3015) Å, σ har , σ hug = 8.7 × 10 −22 , 1.15 × 10 −6 m 2 , and F 1,hug and F 2,hug are the integrated solar fluxes in the 280.5-305.5 and 305.5-360 nm ranges.First, we performed separate tests of these parameterizations, which showed that the parameterizations for HAR and HUG are in a good agreement with libRadtran.However, for LYA and SRB, according to the test results, we changed σ lya and added altitude-dependent x srb .The outcome of these tests is presented of these tests are presented in Fig. 4. Because the original ECHAM schemes can partly reproduce the response of the heating rate to the solar UV variability obtained with the reference scheme, we apply these parameterizations to cover only the missing part of the signal.The scaling coefficients for each of the four applied parameterizations were calculated from the following system of equations: where m is the number of levels in vertical direction, n is a number of unknown k coefficients, A ij is an n column m row array containing heating rate difference between solar maximum and minimum calculated with LYA, SRB, HAR and HUG parameterizations, and B i and C i are an m element vectors containing the same difference calculated with the reference model and original ECHAM codes.In our case, m = 42 and n = 4, meaning that the system of equations is overdetermined and has no exact solution.The approximate solution of the system is then calculated by the standard leastsquares procedure from the IDL (Interactive Data Language) linear algebra package library.The set of scaling coefficients was calculated separately for E5c and E6 and is presented in Table 2. Since E5c does not have original absorption by oxygen and therefore underestimates the absolute values in the mesosphere, the heating parameterizations for LYA and SRB have been added to the original scheme using the full flux integrated within a specific band in order to improve the scheme in respect to the calculation of the absolute heating rates.However, to avoid an overestimation in the upper stratosphere, related to the fact that the original codes partially treat O 3 absorption in the Hartley and Huggins bands, we recommend not using the full radiative flux but the difference between solar minimum and maximum.The same should be done for LYA and SRB in E6 to  the mesosphere since the absolute values in the mesosphere are already reproduced well.In global models this can be done choosing the month or day with the lowest SSI in which all extra heating will be equal to 0, and then, for calculations at all other dates, one should use the SSI difference from this "grand minimum" value.The implementation of the proposed parameterizations does not require any retuning of the original codes, and another important advantage is that these parameterizations take negligible computer time compared to the time taken by radiation schemes.

Changing UV
Figure 5 shows the improvement of the original schemes' performance due to the implemented parameterizations of O 2 and O 3 absorption calculated under changing UV and constant ozone conditions for tropical standard atmosphere and a solar zenith angle equal to 10 • .The implemented parameterizations of O 2 and O 3 absorption allowed us to get very good agreement in solar-variability-induced heating rate changes with the reference model in the mesosphere and the stratosphere.The only notable difference appears in the lower mesosphere around 67 km, but we suspect that this is the artefact of the vertical resolution used.For E5c, since we used the full radiative flux for LYA and SRB, we also improved the representation of the absolute values in the mesosphere.Both radiation schemes, E6 and E5c, overestimate the total heating rate in certain regions in absolute values compared to libRadtran.This overestimation is a feature of original schemes and is larger during is larger during solarminimum conditions since E5c and E6 underestimate the additional heating through spectral irradiance variability over the 11-year solar cycle.By the inclusion of the additional parameterizations, the extra heating rate maximally reaches 0.06 K day −1 for E6 and 0.21 K day −1 for E5c around 46 km during the solar maximum; it, therefore, decreases the discrepancy of E5c+ and E6+ with libRadtran, which is now constant in time.In a transient simulation such deviation will be always equal to the difference during the "grand minimum".
Results of calculations with four other different atmosphere models (midlatitude summer, midlatitude winter, subarctic summer, subarctic winter (McClatchey et al., 1972)) and three solar zenith angles (10, 40, 70 • ) presented in Fig. 6 show that the parameterizations work well for all conditions, and the applied scaling coefficients do not depend strongly on the position of the Sun and latitude and can be used in models with high confidence.It should be noted that, for other radiation schemes and other SSI data sets, these coefficients will differ and have to be calculated carefully with regard to the specific features of each particular scheme.

Changing ozone
For the previous calculations we used only changing UV fluxes with a constant ozone profile, but the ozone profile is modulated by solar irradiance changes, affecting the irradiance propagation.To check the parameterization applicability by taking into account the ozone feedback we also calculated the heating rate response to the solar-induced ozone changes, keeping the UV fluxes unchanged.Results of these calculations are shown in Fig. 7.In this case the original codes work well, and since we use irradiance differences to calculate extra heating, we do not affect heating rates by ozone changes because extra-heating rates in this case are equal to 0. The total heating rate (UV + ozone) also looks good compared to the reference model.

Conclusions
We evaluated the performance of the ECHAM4, six-band ECHAM5 and ECHAM6 radiation codes in the representation of the solar-UV-variability-induced changes in the heating rates.All schemes showed high underestimation in the mesosphere.In the stratosphere, ECHAM4 code is able to reproduce only 10 % of the reference solar signal, while six-band ECHAM5 code misses 45 % and ECHAM6 code misses about 15 %.We suggested an accurate method to correct the problems revealed: the implementation of parameterizations of extra heating due to oxygen and ozone absorption.This approach was implemented in the six-band ECHAM5 and ECHAM6 schemes and allowed us to get very good agreement with the reference model in the representation of the solar signal in the mesosphere and stratosphere without a significant increase in computational time.This method does not require tuning of the original codes, but it only provides the solar-induced addition to the original heating rates.Therefore, this method is suitable for any other radiation scheme to correct the solar signal in heating rates due to missing or underrepresented spectral intervals.

Figure 1 .
Figure 1.Variability of solar irradiance in the 120-440 nm wavelength range calculated by COSI (left) and resulting ozone response from a composite of observational data from Soukharev and Hood (2006) and Austin et al. (2008) (right).

Figure 2 .
Figure 2. Shortwave heating rates in K day −1 for tropical standard atmosphere and solar zenith angle equal to 10 • calculated by E5c and E4 (left panels) and E6 (right panels).Top panels: absolute values during solar minimum.Bottom panels: differences between minimum and maximum (max-min) of the 11-year solar cycle.Solid lines: ECHAM results.Dotted lines: libRadtran results for the same spectral intervals.Different spectral intervals are designated by colours; yellow line -E4 250-680 band; black dashed line -libRadtran results for 120-440 nm (i.e.including shortest wavelengths > 120 nm).

Figure 4 .
Figure 4. Shortwave heating rate differences of the 11-year solar cycle (solar max minus solar min) in K day −1 for tropical standard atmosphere and solar zenith angle equal to 10 • calculated by extra heating parameterizations and libRadtran.Solid lines: results of parameterizations.Dotted lines: libRadtran results for the same spectral intervals (Table1).

Figure 5 .
Figure 5. Shortwave heating rates in K day −1 for tropical standard atmosphere and solar zenith angle equal to 10 • .Left panel: differences between minimum and maximum (max-min) of the 11year solar cycle in the case of UV only variability and a constant ozone profile.Right panel: absolute values during solar maximum.Coloured solid lines: results from original E5c and E6 codes.Black solid line: libRadtran results for reference.Dashed lines: results from improved parameterizations.

Figure 7 .
Figure 7. Shortwave heating rate differences (solar max minus solar min) of the 11-year solar cycle in K day −1 for tropical standard atmosphere and solar zenith angle equal to 10 • .Left panel: including only ozone changes.Right panel: UV + ozone changes.Original codes results are denoted by solid lines, improved codes results by dashed lines.

Table 1 .
ECHAM solar radiation scheme spectral intervals and main absorbers in the UV part of spectrum.

Table 2 .
Wavelength intervals and scaling coefficients of the extra heating parameterizations.