Improved simulation of clouds over the Southern Ocean in a General Circulation Model

Abstract. The present generation of global climate models is characterized by insufficient reflection of short-wave radiation over the Southern Ocean due to a misrepresentation of clouds. This is a significant concern as it leads to excessive heating of the ocean surface, sea surface temperature biases, and subsequent problems with atmospheric dynamics. In this study we modify cloud micro-physics in a recent version of the Met Office's Unified Model and show that choosing a more realistic value for the shape parameter of atmospheric ice-crystals, in better agreement with theory and observations, benefits the simulation of short-wave radiation. In the model, for calculating the growth rate of ice crystals through deposition, the default assumption is that all ice particles are spherical in shape. We modify this assumption to effectively allow for oblique shapes or aggregates of ice crystals. Along with modified ice nucleation temperatures, we achieve a reduction in the annual-mean short-wave cloud radiative effect over the Southern Ocean by up to 4 W/m2, and seasonally much larger reductions. By slowing the growth of the ice phase, the model simulates substantially more supercooled liquid cloud. We hypothesize that such abundant supercooled liquid cloud is the result of a paucity of ice nucleating particles in this part of the atmosphere.



Introduction
One of the major known drawbacks in the present-day global climate models is an excess in the absorbed short-wave (SW) radiation over the Southern Ocean (SO) (Trenberth and Fasullo, 2010;Ceppi et al., 2012;Hwang and Frierson, 2013;Hyder et al., 2018). In Chapter 9 of the 5 th Assessment Report of the Intergovernmental Panel on Climate Change (IPCC AR5) (Flato 15 et al., 2014), it points out that most of the 5 th Coupled Model Intercomparison Project (CMIP5) models (Taylor et al., 2012) have a positive SW cloud radiative bias of magnitude of up to 20 Wm −2 over the SO, suggesting that inadequately simulated clouds allow substantially too much sunlight to reach the ocean surface.
Several studies have focused on the relation between various aspects of cloud representation in the model and observed radiation biases. Bodas-Salcedo et al. (2012, using cyclone compositing cluster analyses, suggest the need to increase 20 the optical depth of the low-level clouds and improve the simulation of mid-level cloud regime, to help reduce the biases in the model. By modifying the shallow convection detrainment in their global climate model, Kay et al. (2016) showed that the resultant increase in the supercooled liquid clouds (SLC) enable large reductions in long-standing climate model SW radiation biases. By implementing a new parametrisation that includes the turbulent production of mixed-phase clouds, Furtado et al. (2016) show that the radiation biases can be substantially improved, especially over the SO. In another study by Furtado and Field (2017), the importance of ice micro-physics parametrisation in determining the phase composition, and thus the liquid water content of the SO clouds is highlighted.
Discrepancies in the response of clouds to anthropogenic forcings are recognized as a leading reason for a persistent, large spread in the climate sensitivity throughout various generations of climate models . We thus conjecture that 5 this model problem contributes to this large spread, and thus solving it would increase confidence in projections of anthropogenic climate change (Tan et al., 2016).
In the present study, we investigate the role of parameters involved in atmospheric ice formation within a global climate model in causing the above mentioned SW radiation bias. Here, we define a SO region as the latitudinal band between 50°S and 70°S. The control climate model used in this study is an accrual of the most recent version of the Met Office's Unified Model, GA7.1 (Walters et al., 2019) with modified micro-physics scheme for riming process and several other scientific changes. Appendix A summarizes the scientific set-up for this model version. The resolution used here is N96L85 (i.e. a horizontal resolution of 1.875 • ×1.25 • and 85 terrain-following hybrid-height levels extending to 85 km of altitude). It uses the "ENDGAME" dynam- 15 ical core with a semi-implicit semi-Lagrangian formulation to solve the non-hydrostatic, fully compressible deep-atmosphere equations of motion (Wood et al., 2014)

Model set-up
In the present study, control run follows the Atmospheric Model Intercomparison Project (AMIP) climate model development protocol (Gates et al., 1999;Schuddeboom et al., 2019), using prescribed sea-surface temperature climatology. Excess atmo- 20 spheric ice has been a persistent concern in the control version of the model ( fig. 1), which is especially pronounced over the SO region. Ice clouds have a significant influence on the global climate through their effects on the Earth's radiation budget e.g. (Hartmann and Doelling, 1991;Waliser et al., 2009). Hence, sensitivity set-ups in our study are aimed at modifications to the micro-physics scheme such that the ice growth in the model is controlled. We achieve this by modifying those parameters that control the growth of existing ice by vapor deposition and heterogeneous nucleation of new ice. The classical theory of ice 25 crystal growth uses an electrostatic analogy due to the similarity between the equations governing the water vapor distribution around an ice crystal and the electrostatic potential distribution around an electric conductor of the same shape as the ice crystal (Chiruta and Wang, 2003). Thus, the growth rate of ice crystals by diffusion depends on a shape (also known as capacitance) parameter C, which is a function of both ice crystal size and habit. To determine the ice crystal growth rates in models, it is necessary to know the value of C (Chiruta and Wang, 2003;Hobbs, 1976). The standard equation that is used for calculating 30 the growth rate of ice crystals in the model is, where D v is the diffusivity of water vapor in air, C is the capacitance, ρ α and ρ s are distributions of vapor densities at and away from crystal's surface.
From eq.1, it is evident that once the value of capacitance C is known, the growth rate of ice crystals can be determined. All other quantities on the right-hand side of eq.1 are independent of the shape (Chiruta and Wang, 2003;Hobbs, 1976). Thus, C in the model effectively defines the shape of ice crystals, which in turn is fed through to the ice processes of deposition/sublimation 5 and melting without affecting any other ice processes.
Technically the capacitance, C, is defined as 1.0 x d in the model, where d is the particle maximum size. In our sensitivity studies, we modified the value to 0.5 x d (corresponding to any oblate ellipsoid with two unequal axes, thought to be more appropriate for aggregates and plate-like crystals rather than the assumption of spherical crystals alone). Our value of 1.0 or 0.5 is a non-dimensional capacitance (Field et al., 2008). The effect of this change in the shape parameter is tested independently 10 as well as in combination with changing the temperatures at which heterogeneous and homogeneous freezing start in the cloud micro-physics scheme. The ice nucleation temperature is the temperature at which heterogeneous nucleation of ice first starts to occur in the model. The default value of −10 • C was changed to −40 • C and −20 • C, to investigate the effect of delaying the heterogeneous ice nucleation in the model. Two further parameters that were modified in the parametrised convection scheme that control ice formation in the model are the temperature at which detraining condensate as ice begins in the model (start- 15 ice temperature) and the temperature at which all condensate is detrained as ice (all-ice temperature). The values used in our numerical simulations are summarised in Table 1.
The ice nucleation temperature is meant to be similar both in large-scale and convective cloud schemes. Hence, the experiment where the nucleation temperature is reduced to -40°C (i.e. exp2) is physically unrealistic. However, it is still a useful sensitivity scenario to study the importance of detrained ice vs. large-scale freezing. All simulations were run for twenty years 20 under steady-state present-day conditions.

Observational data
We use the National Aeronautics and Space Administration (NASA) Clouds and the Earth's Radiant Energy System -Energy Balanced And Filled (CERES EBAF Ed4.0, Terra-Aqua) surface and top-of-the-atmosphere (TOA) data set, covering the period 2000 to 2018 as an observational reference for radiative fluxes. This data set, in an earlier version (Loeb et al., 2009), was also 25 used in AR5. The overall uncertainty in the monthly all-sky TOA flux for the CERES EBAF Ed4.0 data set is estimated to be 2.5 Wm −2 (for for both SW and LW fluxes). For clear-sky TOA, uncertainties in SW and LW fluxes are 5 Wm −2 and 4.5 Wm −2 respectively (Loeb et al., 2018). We also use the European Centre for Medium-Range Weather Forecasts (ECMWF) Re-Analysis 5 (ERA5) monthly mean data for comparison of cloud-ice content (ERA5, 2017).
3 Results where both capacitance and nucleation temperature are modified will have an added impact on the IWP. 10 Zonally averaged distribution of IWP and LWP, over both hemispheres, for all boundary layer types for annual and seasonal means are provided in the supplementary material (figs. S1 to S4).  fig. 3a), the decrease in LW radiation at the TOA, in absolute terms, is larger than the increase in SW TOA over the SO. This is visible in the distribution of net radiation at the TOA (i.e. LW plus SW at TOA) as well (solid mustard lines in figs. 3a to 3c). For exp1, there is an increase in the net outgoing TOA radiation whereas for exp2 and exp3, it shows a decrease over the SO region.   that of the CERES EBAF observational data. As evident, there has been significant improvements of the SO radiation biases in GA7 compared to its predecessor GA6. One of the major reasons behind this improvement was a better representation of mixed-phase clouds and supercooled liquid in the cloud micro-physics scheme (Furtado et al., 2016). Fig. 6c  An overestimation of ice in clouds is a known shortcoming of many of the present-day global climate models ( fig. 1). It is coupled to an underestimation of SLC. This problem is of particular importance in the SO region characterized by abundant SLCs (Kay et al., 2016;Bodas-Salcedo et al., 2016;Huang et al., 2012;Hu et al., 2010). When ice and supercooled liquid coexist, the ice grows at the expense of the liquid by the Wegener-Bergeron-Findeisen (WBF) mechanism (Wegener, 1911;5 Bergeron, 1935;Findeisen, 1938). Acknowledging the complexities in representing the many possible background microphysical processes that are responsible for this in a global climate model, the primary idea of modifying the shape parameter of ice-crystals is to reduce the rate of depositional growth of ice particles. This reduction essentially slows down the deposition growth of ice crystals, which leaves more water vapor to be available for condensation into liquid phase particles. At the scale represented in global climate models, for conditions of very low ice-nucleating particle (INP) concentrations and temperatures Our choice of 0.5 * d (d being the particle maximum size) for capacitance is based on theory and observational studies (Field et al., 2008). The atmosphere-only model studied here does perform better with this value than with the default value of 1.0 * d. The SW radiation over SO is improved but results are more mixed for the other fluxes (figs. 3 and 4). The uncertainty in the surface radiation budget observations also needs to be considered. As already noted earlier, the experiment where the 20 nucleation temperature is reduced to -40°C (i.e exp2) is physically unrealistic and is intended to be a useful sensitivity scenario to study the importance of detrained ice vs. large-scale freezing.
Even though there is noticeable reduction in the SW radiation bias over the SO in all the experimental scenarios ( fig. 5), we recognize persisting shortcomings in this regard in other parts of the world (figs. 6c. Similar to the control version in this study, certain regions have not shown much of an improvement in terms of the SW CRE bias (e.g. the Bay of Bengal, areas around 25 southeast Asia, eastern south Pacific etc). Previous studies have suggested that some cloud micro-physics parameterisations produce unrealistically bright clouds, especially over the Northern Hemisphere (NH) (Furtado et al., 2016). Since the SW biases over the NH were smaller than those over the SO, a significant brightening of modelled NH clouds is undesirable (Furtado et al., 2016). While the changes to nucleation temperature has significant impact on the tropics as well, the capacitance changes are more localized to the high latitudes.

Conclusions
In this study we improve the SW radiation biases in a recent version of the UK Met Office's Unified Model. This and other contemporary climate models are characterized by excess cloud ice causing biases in SW radiation biases which are especially 5 pronounced over the SO. Here, we modify the capacitance or shape parameter which represents ice crystal shape and habit. In our sensitivity studies, we reduce this parameter from 1 to 0.5 (corresponding to any oblate sphere shape in general, where the horizontal axes are longer than the vertical axis and more representative of an aggregate or flat ice crystal). We also examine the impact of changing other temperature thresholds in the cloud micro-physics scheme for the onset of heterogeneous ice production. Our analysis shows that the SW radiation bias has significantly reduced over the SO after the modification of these parameters. However, disparities still exist in other regions. INPs that are currently not represented in the cloud micro-physics scheme might be a factor in this model behavior. The fact that nucleation temperature changes currently is associated with the same effects globally is undesirable, it further motivates the future work to couple the nucleation temperature to a prognostic or the least a regionally specified INP concentration. Several changes were introduced in the GA7.1* control model used in this study relative to its predecessor GA7.1 (Walters et al., 2019;Brown et al., 2012). These changes range from minor bug fixes and optimisation techniques to major science 20 changes. As far as our study is concerned, the main modification to GA7.1 is the inclusion of the modified micro-physics scheme which includes a shape dependence of riming rates using the parameterization by Heymsfield and Miloshevich (2003), as a measure to prevent small liquid droplets from riming (Furtado and Field, 2017