Recent global nonhydrostatic modeling approach without using a cumulus parameterization to understand the mechanisms underlying cloud changes due to global warming

Clouds are the primary source of uncertainty in the prediction of climate change. To reduce the uncertainty of cloud simulations and overcome this difficulty in prediction, many climate modeling centers are now developing a new type of climate model, the global nonhydrostatic atmospheric model, which reduces the uncertainty arising from a cumulus parameterization by computing clouds explicitly using a cloud microphysics scheme. Among the global nonhydrostatic atmospheric models used in recent intercomparison studies, NICAM aims to project climate change by improving our understanding of cloud changes due to warming and related physical processes. NICAM is the first global nonhydrostatic model and was developed by our research team. This review summarizes the outcomes of a recent major five-year research program in Japan for studying climate using NICAM, as well as providing an overview of current issues regarding the use of global kilometer-scale simulations in high-resolution climate modeling.


Introduction
Clouds play an important role in the change in earth's surface temperature due to global warming; however, the large uncertainty in the projection of cloud changes makes it difficult to better predict climate change.The Coupled Model Intercomparison Project (CMIP) phase 6 data, the latest dataset of climate change simulations produced by general circulation models (GCMs) from world climate modeling centers, has updated and improved our understanding of cloud changes and their underlying mechanisms (Zelinka et al. 2020(Zelinka et al. , 2022)).In the CMIP6 dataset, climate sensitivity is the highest ever projected.
Climate scientists have attempted to reduce the wide uncertainty in climate sensitivity in GCM datasets; however, the CMIP6 dataset has a larger range of climate sensitivity than past evaluations (Meehl et al. 2020;Zelinka et al. 2020).In particular, 10 of the GCMs used in CMIP6 show climate sensitivity exceeding 4.5 K.The high climate sensitivity is attributed mainly to a larger positive feedback due to a reduced cloud albedo in middle and high latitudes (Zelinka et al. 2020), indicating the importance of clouds to climate.In fact, the CMIP6 dataset shows that cloud feedback is the biggest cause of the large uncertainty in climate sensitivity.
To reduce the uncertainty in cloud feedback, Sherwood et al. (2020) proposed an alternative evaluation method that takes into account past observational results, called the expert assessment.The latest results from the expert assessment show larger uncertainty in feedback caused by high clouds, because most GCMs fail to simulate the interannual variation in tropical clouds and associated upward radiation (Mauritsen and Stevens 2015;Williams and Pierrehumbert 2017;Zelinka et al. 2022).Zelinka et al. (2020) show that the large cloud feedback comes not from a single anomalously large component but from a systematic biased high.These findings suggest a need for efforts to improve the overall performance of cloud simulations.
Large uncertainty in the projection of clouds arises mainly from the fact that most GCMs simulate deep convection using a cumulus parameterization, as they cannot explicitly compute deep convection.It is difficult to model cloud processes sufficiently and to adequately obtain the relevant physical parameters.To improve cloud projections, several climate modeling centers have launched large inter-organizational or international research projects to develop global nonhydrostatic modeling with kilometer-grid spacing (e.g., Slingo et al. 2022;Mauritsen et al. 2022).In Japan, the world's first global nonhydrostatic model, Nonhydrostatic ICosahedral Atmospheric Model (NICAM), was developed in the early 2000s to improve a representation of clouds in a GCM (Satoh et al. 2014) A consistent result of those research projects is the increase in a coverage of high clouds and decrease in ice water path in the tropical atmosphere (Iga et al. 2007;Satoh et al. 2012;Noda et al. 2012).Recent studies also show that a larger number of optically thin and smallsize clouds contribute to the increase in high clouds in a warmer atmosphere (Noda et al. 2016;Chen et al. 2016;Satoh et al. 2018).However, the underlying mechanisms are not clear.Bretherton et al. (2005) emphasized the role of aggregation in deep convection due to warming.
For a deeper understanding of cloud changes due to warming, we conducted not only Atmospheric Model Intercomparison Project (AMIP)-type (Webb et al. 2017) simulations that consider land-ocean distributions and topography, but also idealized aqua-planet-type simulations, so-called radiative convective equilibrium (RCE) experiments.Such simpler planet conditions allow us to systematically examine the responses of physical parameters to cloud changes in order to study how the roles of physical processes, including cloud microphysics processes and subgrid-scale (SGS) turbulent mixing processes, in cloud formation would change in a warmer atmosphere.
In the TOUGOU program, we improved our understanding of cloud changes due to warming using NICAM and improved how physical processes are modeled in NICAM, including low-level mixed-phase clouds and SGS turbulent mixing.The results of this program are relevant to recent research activities developing highresolution GCMs (Stevens et al. 2019;Slingo et al. 2022), research into cloud changes due to warming.
The purpose of this review is to present the findings of the TOUGOU program and to discuss related issues.The review is organized as follows.Section 2 reviews research activities based on the results of NICAM that was conducted in the TOUGOU program; in particular, Sect.2.1 outlines findings regarding changes in clouds in response to global warming, and Sect.2.2 describes notable improvements in physical processes.Section 3 documents the future direction of cloud studies using a high-resolution GCM, current related research activities, and future issues of high-resolution numerical studies.The conclusion is provided in Sect. 4.

Projection of cloud changes due to warming 2.1.1 Convective aggregation
To evaluate the increase in earth's temperature, it is important to understand how anvil clouds respond to warming: the fractional coverage of clear-sky regions increases as tropical deep convection becomes more aggregated.In such an atmospheric state, the earth can release a larger amount of heat from the atmosphere into space without being shaded by anvil clouds, and vice versa (Fig. 1).To gain knowledge of potential future changes in the aggregation of tropical deep convection, we analyzed 14 km mesh NICAM simulation data for both present and future climate, assuming a world after a century (Kodama et al. 2015) to investigate changes in the fractional coverage of downdraft regions, which is often used as an index representing the degree of aggregation of deep convection (Coppin and Bony 2015).
We first divide the tropics into 10° × 10° subdomains and compute a mean of vertical velocity at an altitude of 5 km, < w > , over each subdomain.Then, the modi- fied subsidence fraction, SF' , is defined as the fractional coverage of negative w′ , where w′ = w− < w > in each 10° × 10° subdomain to exclude contributions of largescale circulations.Smaller SF' values correspond to the state that convection develops in narrower regions (i.e., convection is more aggregated).Figure 2 shows tropical distributions of < w > and SF' and their changes due to warming.The distribution of SF' corresponds well with < w > in both the present and future climate.The two regions with the greatest difference in < w > between the present and warmer conditions are the Indian Ocean (negative change) and Central Pacific (positive change) in the tropics, showing a weakened Walker circulation in a warmer world.The difference in SF' is mostly negative, except for some regions such as the Central Pacific.
Furthermore, we focus on the relationship between < w > and SF' in Fig. 3.The two values are correlated with each other, and they have the same correlation coefficient of 0.88, even in the different climate states; this means that the relationship between the intensity of large-scale circulations and convective aggregation is quite similar in the different climate states.In a warmer climate, the degree of convective aggregation becomes weaker and corresponds to reduced large-scale circulations.This result is consistent with past studies based on NICAM simulation data that showed larger numbers of smaller-scale clouds in a warmer world (Noda et al. 2014(Noda et al. , 2016)).It also suggests that high clouds that originate from deep convection cover larger parts of the tropical atmosphere because they develop to be more scattered (Fig. 1).

Response of ice hydrometeors to surface warming
We also examined how cloud condensates respond to changes in sea surface temperature (SST) using 14 km mesh NICAM (Kodama et al. 2015), focusing on the link between cloud and ice hydrometeors, such as cloud ice, snow and graupel, which are not ordinarily resolved in conventional GGMs.The vertical distribution of the tropical mean cloud fraction is plotted against tropical mean SST (regression) in Fig. 4a, for both the present and warming climate in AMIP simulations with NICAM and GCM-oriented Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations (CALIPSO) Cloud Product.For the present climate, the cloud fraction and SST have a positive (negative) relationship at altitudes above (below) 14-16 km.This behavior is also found in observations (Zelinka and Hartmann 2011;Chen et al. 2022), although the relationship is stronger in NICAM simulations than in observations.For the warming climate, the relationship between the vertical profile of the tropical mean cloud fraction and the SST is similar to that in the present climate, but with an upward shift: the peak of the positive relationship appears at an altitude around 18.5 km under warmer conditions.This upwardly shifted profile corresponds to the increase in height of high clouds in the warming climate.Since the high cloud is formed by ice hydrometeors, such as cloud ice, snow and graupel, the distribution of ice water content contributed by each ice hydrometeor in relation to SST is also illustrated in Fig. 4b.For the present climate, the vertical profiles of cloud ice and snow are consistent with the cloud fraction shown in Fig. 4a, although the peak of the positive regression is lower than that for the cloud fraction.In contrast, the contribution of graupel dominates mainly in the lower part of the high cloud region and does not have a large influence on the distribution of the cloud fraction.In the warming climate, the vertical profiles of ice water content contributed by cloud ice, snow, and graupel shift upward.The amplitude of the annual variation in the SST has a similar scale in the present and warming climates, while convective activities become stronger in the Central/Eastern Pacific and the cloud fraction doubles in the warming climate (data not shown).The fact that the variation in the cloud fraction in the warming climate is double that in the present climate may be strongly related to the difference in convection activities (Chen et al. 2022).

Changes in the responses of cloud microphysics processes based on radiative convective equilibrium experiments
One of the advantages of using global nonhydrostatic models for climate projection science is that they explicitly calculate processes in cloud layers and can relate cloud feedbacks to changes in the processes associated with global warming.Here, we show how the processes controlling high clouds change with global warming and how these changes relate to changes in the high-cloud coverage, focusing on the cloud microphysics.In order to simplify the problem, high clouds and processes were analyzed using RCE simulations.The RCE simulations were conducted with fixed SST, no diurnal cycle of radiation, and no rotation, using an earth-size spherical computational domain.The qualitative similarity of highcloud cover response to SST changes in RCE and in more realistic simulations using NICAM have been examined by Ohno and Satoh (2018) and Ohno et al. (2019Ohno et al. ( , 2020Ohno et al. ( , 2021)).First, we examined the dominant cloud microphysical processes for high clouds.Figure 5 shows domain-averaged vertical profiles of the tendencies of cloud ice due to cloud microphysical processes at the equilibrium states at an SST of 300 K and 304 K.It is evident that the tendencies of cloud ice were dominated by the sedimentation, diffusional growth, sublimation and collection processes near the cloud top layers (~ 220 K).The sedimentation and collection rates depend on the terminal velocity of ice particles (Pruppacher and Klett 2010).The terminal velocity of particles generally increases with the decrease in pressure due to the reduction of the aerodynamic effects (e.g., Heymsfield 2007), which suggests that the sedimentation and collection are enhanced at lower pressures.The deposition and sublimation rates are proportional to the diffusivity of water vapor in air (Pruppacher and Klett 2010), which decreases with an increase in air pressure (Montgomery 1947;Hall and Pruppacher 1976).These facts indicate that the reduction in air pressure in high-cloud layers due to the elevation of high clouds associated with global warming enhances the dominant cloud microphysical processes for high clouds.This can lead to the shortening of the lifetime of high clouds and the reduction of high-cloud cover in warmer climates.
The effects of pressure dependencies of these processes on the high-cloud cover and its response to changes in SST were examined using sensitivity simulations.We conducted simulations using a value of 440 hPa as the lower pressure limit for evaluations of the terminal velocity of ice particles (VTp440) and the water vapor diffusivity in air (KVp440).The value of 440 hPa was chosen based on the upper limit of the cloud-top pressure as per the International Satellite Cloud Climatology Project definition of clouds.Figure 6a and b shows the high-cloud cover and its response to changes in the SST in the control and sensitivity simulations.The high-cloud covers also increase in the sensitivity experiments.High-cloud cover is reduced by suppressing the pressure dependencies of the terminal velocity of ice particles and the water vapor diffusivity in air (Fig. 6c).The estimated effects of the pressure dependencies of the terminal velocity of ice particles and of the water vapor diffusivity in air on the high-cloud response to SST change were both negative (Fig. 6d).These results suggest that the reduction in the time scales of the dominant processes due to the upward shift of high clouds in warmer climates reduces high cloud cover (Fig. 7).Note that Fig. 6d also shows complicated responses in different optical thickness clouds, e.g., the result of VTp400 shows a remarkable decrease of thin clouds, and small decrease and increase of medium and thick clouds, respectively, while that of KVp400 indicates a small increase of thin clouds, and remarkable increases of medium and thick clouds.Further research is needed to clarify the reason.

Projection of changes of tropical cyclones
It is also important to evaluate how tropical cyclones (TCs) change in a future climate, as they are major meteorological phenomena causing severe damage in a large part of the world.In order to deepen our understanding of the response of TCs to global warming, we investigated future changes in TC structure and TC seeds (Yamada et al. 2017(Yamada et al. , 2021)).The projected change in TC frequency in a future climate differs among different models (Knutson et al. 2020).Some studies have noted the response of TC seeds (incipient vortices of TCs) to global warming and the relationship of this response to the frequency of TC genesis (Vecchi et al. 2019;Hsieh et al. 2020;Lee et al. 2020;Sugi et al. 2020).Following Vecchi et al. (2019), Yamada et al. (2021) used the outputs of six models in CMIP6 HighResMIP (Haarsma et al. 2016; Roberts et al.  2021), © The American Geophysical Union.Used with permission 2020a, 2020b) to decompose tropical cyclogenesis into the contributions of TC seeds and survival rate, which is the rate of TC seeds developing into TCs.They showed that tropical cyclogenesis frequency decreased significantly from 1990 to 2049 in the multi-model ensemble, which is attributed to a decrease in the number of TC seeds.However, the main contributor varied between models and their horizontal resolutions.This indicates that decomposing tropical cyclogenesis into TC seed and survival rate likely addresses the cause of the uncertainty in the projected frequencies of tropical cyclogenesis in previous studies.
TCs cause severe disasters, the magnitude of which is modulated by TC intensity, path and size.How these features will respond to global warming is not fully understood or whether (Knutson et al. 2020).Although GCMs are a useful tool for assessing future changes in TC activity, a GCM with a lower resolution tends to produce TCs with a spatial scale larger than that of observed TCs (Camargo et al. 2005).Due to improvements in computing power and modeling, we were able to run a climate simulation with a finer horizontal grid spacing of 14 km (Kodama et al. 2015).While the 14 km mesh is not likely to completely reproduce the finer structures of TCs, it did reproduce broad structures like the primary and secondary circulations as well as the warm core (Yamada et al. 2017).As horizontal scales of TCs differ from TC to TC (e.g., Wu et al. 2015), Yamada et al. (2017) evaluated the horizontal scale of TCs as a function of their lifetime maximum intensities.The result indicated that the radial-averaged tangential wind of TCs with a central From Ohno et al. (2021) pressure less than 980 hPa will be enhanced outside the eye wall cloud under warmer climate conditions, even if their minimum central pressures are in the same category.This enhancement possibly expands the horizontal scale of TCs to a value such as the radius of gale force wind under warmer climate conditions.
A future change in the spatial distribution of mean TC size is also important in projections of socio-economic damage.We regarded the TC size as a radius of 8 m s −1 wind (R08), as per Schenkel et al. (2022).R08 was calculated from outputs of NICAM AMIP-type simulations (Kodama et al. 2015;Satoh et al. 2015;Yamada et al. 2017) and was defined as the outermost radius exceeding 8 m s −1 in the radial profile of azimuthal-averaged tangential wind speed at a height of 10 m.The simulation for 1979-2008 was regarded as the control simulation, and the simulation for 2075-2104 was regarded as the global warming simulation.Figure 8 shows the spatial distributions of mean R08 and its future change.R08 varied by region and increased with latitude, which is consistent with previous studies (e.g., Merrill 1984;Chavas and Emanuel 2010;Schenkel et al. 2022).The dependence on latitude differed between the northern and southern hemispheres (Fig. 8d).Zonal means of R08 increased between 10° S and 20° S in the southern hemisphere than in the northern hemisphere due to warming.In terms of the future change, although the pattern of change was complicated on the global ocean area, R08 became larger over the Arabian Sea, the tropical western North Pacific, the east coast region of the USA, and the north coast of Australia (Fig. 8a-c).

Extratropical cyclones in a warmer climate
Extratropical cyclones play a key role in day-to-day weather at midlatitudes.Their large-scale dynamical features may be represented satisfactorily by coarse-resolution GCMs.However, extratropical cyclones consist of fine-scale dynamical and microphysical features such as cold and warm fronts and cloud precipitation systems, and a much higher resolution model could better reproduce the precipitation and intense wind associated with such features.Catto et al. (2019) explained that projected future changes in the precipitation area, fronts, and wind are highly uncertain and suggested the need for higher-resolution models.Kodama et al. (2019) (hereafter, K19) investigated the response of extratropical cyclones to global warming using 14 km mesh NICAM.In this review, we perform a preliminary analysis of High-ResMIP NICAM data using a similar approach as K19 and show that cyclone intensity and precipitation in a future climate depend on the model resolution.
K19 showed that the change in intensity of oceanic extratropical cyclones (those staying over the ocean for more than half of their lifetime) due to warming is not evident.This result is also confirmed in this study using 14-56 km mesh NICAM (not shown).Close inspection reveals that the number of moderately intense extratropical cyclones (e.g., the wind speed between 30 and 50 m s −1 at the 850 hPa level) decreases with time in the southern hemisphere (Fig. 9b-d).Note that the extratropical cyclones with stronger wind occur more frequently as horizontal resolution increases in both southern hemisphere (Fig. 9a) and the northern hemisphere (not shown), implying the need for 14 km (and maybe even finer) mesh models to assess the extreme wind speed, such as higher than 50 m s −1 , associated with extratropical cyclones.
Figure 10 shows changes in precipitation around the center of the oceanic extratropical cyclones due to global warming.Precipitation increases mainly at the poleward side of the direction of cyclone movement; this increase is larger around the intense extratropical cyclones than around all the extratropical cyclones.These preliminary results, consistent with those of K19, do not depend on the horizontal resolution of the model.Quantitatively, the simulated precipitation increases with the model resolution, and this trend is more pronounced when averaged over only the intense extratropical cyclones.Such a dependence of precipitation change on resolution seems to correspond to the change in surface air temperature (not shown), and further in-depth analysis is needed to understand the cause of such a dependency.

Evaluation of effective climate sensitivity
In the present project, we also attempted to evaluate effective climate sensitivity (ECS) using data from global nonhydrostatic climate simulations, for the first time in climate research.According to Shiogama et al. (2014), we can calculate ECS as, where and R and T denote the global means of radiation at the top of the atmosphere and of the earth's surface temperature, respectively.CTL, 4 × CO 2 , and SST + 4 K refer to experiments on the present climate, quadruple atmospheric CO 2 concentration, and the increase in SST at 4 K homogeneously over the globe in the present climate, ( 1) (3) respectively.Each simulation was conducted for 5 years.Figure 11 shows change over time in the ECS calculated from the data for year one to year five.Over that time, ECS gradually becomes closer to 3.6 ~ 3.7°, values that are in the uncertainty range of the CMIP5 model results (Andrews et al. 2015).Those conventional GCMs compute deep convection based on a cumulus parameterization, while NICAM does so explicitly based on a cloud microphysics scheme.It is notable that the results of two such different types of GCMs show similar ECS values.

Model improvement 2.2.1 Low-level mixed-phase clouds
This subsection describes the improvement in supercooled liquid water simulations of low-level mixed-phase clouds over the Southern Ocean.The NICAM singlemoment bulk scheme with six water categories (NSW6) (Roh and Satoh 2014;Roh et al. 2017), which was used for the CMIP6 project (Kodama et al. 2021), was revised in reference to the NICAM double-moment bulk scheme with six water categories (NDW6) (Seiki and Nakajima 2014;Seiki et al. 2014Seiki et al. , 2015) ) and satellite observations.For details of the cloud microphysics schemes, refer to the original publications or a review of the cloud microphysics schemes used for GCMs (Seiki et al. 2022).
The shortwave cloud radiative forcing biases over the Southern Ocean in the GCMs had been a longstanding issue in the past few decades (e.g., Bodas-Salcedo et al. 2012;Williams et al. 2013).In addition, high-resolution weather forecasting models with more detailed cloud microphysics schemes had also suffered from the underestimation biases in low-level mixed-phase clouds (e.g., Field et al. 2014).Analyses of the three-dimensional liquid-ice partitioning using CALIPSO indicated that the shortwave biases originate mainly from the underestimation of supercooled liquid water in low-level clouds over the polar regions (e.g., Forbes and Ahlgrimm 2014;Tan and Storelvmo 2016).Recently, GCM communities have revealed that an increase in supercooled liquid water as a result of modifying liquid-ice partitioning efficiently increases the lifetime of low-level clouds over the Southern Ocean (Forbes and Ahlgrimm 2014;Tan and Storelvmo 2016;Kawai et al. 2019).These studies indicate that the ice growth speed in low-level clouds over the Southern Ocean was overestimated in cloud microphysics schemes used for conventional GCMs.Roh et al. (2020) found that low-level clouds derived from the NICAM with the NDW6 scheme effectively represented the characteristics of low-level mixedphase clouds from CALIPSO satellite observations, whereas those from the NICAM with NSW6 remained biased.Seiki and Roh (2020) confirmed that the biases clearly appeared within 10 min of the numerical integration, as was shown in past studies.Therefore, Seiki and Roh (2020) demonstrated that the longstanding biases in ice cloud microphysics schemes can be solved using a single-column model with no physical processes other than cloud microphysics and no external forcing.As it takes only one second to integrate six hours with the single-column model, it was easy to comprehensively test sensitivity experiments.Thanks to the single-column model, all the production and reduction terms in the NSW6 scheme were compared to those in the NDW6 scheme.The initial condition was prepared by simplifying the vertical profiles of typical low-level mixed-phase clouds over the Southern Ocean from the NICAM global simulations with the NDW6 scheme (see Seiki and Roh (2020) for detail).
The single-column model simulation with the NSW6 scheme successfully reproduced the rapid reduction of a liquid cloud layer under the supercooled condition, whereas the simulation with the NDW6 scheme sustained the liquid cloud layer (Fig. 12a-j).The budget analyses indicated that the Bergeron-Findeisen process excessively worked in the NSW6 scheme, as was shown in previous studies (e.g., Tan and Storelvmo 2016;Kawai et al. 2019).In addition, the subsequent occurrence of riming of cloud water by snow and graupel rapidly consumed supercooled liquid water.To alleviate the bias in the NSW6, Seiki and Roh (2020) revised four processes: the initiation of cloud ice was suppressed by changing the ice nucleation scheme; the initiation of graupel through the freezing of rain was suppressed by changing the auto-conversion and accretion schemes for cloud water; the initiation of vapor deposition was delayed with the mixing-ratio threshold of snow and graupel; and the initiation of riming was delayed with the cut-off diameter of snow and graupel.Finally, the lifetime of the low-level mixed-phase clouds simulated using the revised NSW6 scheme was comparable to that obtained using the NDW6 scheme (Fig. 12a and k).
The strategy behind this revision was to suppress the growth of cloud ice, snow, and graupel in the cloud layer under supercooled conditions.This approach can be applied to other atmospheric models that use single-moment bulk cloud microphysics schemes.In fact, the Unified Model and the numerical weather prediction model HARMONIE-AROME applied similar approaches to improve their cloud microphysics schemes (Furtado and Field 2017;Engdahl et al. 2020).By contrast, a fundamental issue remains in NSW6, as in most of the single-moment schemes currently available: the diagnosis of number concentration can be affected by the cloud system.We solve this issue by using the double-moment approach (which predicts the number concentration) only for cloud ice categories.
We also applied this new microphysics scheme to the global domain.We evaluated the new microphysics scheme using 14 km mesh NICAM temporally integrated for one year.We first compare the difference in the fractional mixing ratio of water clouds to the sum of water and ice clouds, F (= q c /(q c + q i )), as a function of air temperature in Fig. 13.The data were computed as follows: first, we smoothed the 14 km mesh data down to 2.5° mesh data, and then averaged them temporally to produce monthly data.Finally, we determined the mixing ratios of cloud water and cloud ice and the temperature at each 2.5° grid value to calculate F(= q c /(q c + q i )).
The new scheme shows larger fractions of water clouds at lower temperatures.In the old scheme, for example, the formation of major ice clouds occurs in a range between − 19 and 0 °C.In the new scheme, however, water clouds develop by − 30 °C.The past observation shows the existence of liquid droplets at an air temperature near − 40 °C, so the result of the new scheme is closer than that of the old scheme to observations.Figure 13 also shows a larger range of F than in the old scheme at the same temperature, implying that more room in cloud feedback associated with mixed-phase clouds is allowed to occur by improving the processes underlying mixedphase clouds.
Next, we compare shortwave radiation at the top of the atmosphere and the liquid water path in Figs.14.The new scheme improves not only clouds over the Southern Ocean, but also over the tropics, and reduces the bias in overestimating the reflection of incident solar radiation over the globe, as was found in the old scheme.The values of the liquid water path increase almost over the entire globe and are closer to the observation, showing that the improved simulation of optically thick water clouds leads to a better simulation of shortwave radiation at the top of the atmosphere.We also confirmed that negative biases in the shortwave radiative field, and relevant positive biases in the shortwave cloud radiative effect, are also reduced greatly throughout the year (Fig. 6 in Noda et al. 2021).

Cirrus clouds
This subsection describes the improvement in cirrus cloud modeling mainly through the revision of ice terminal velocity and collisional growth in the NDW6 scheme (Seiki and Ohno 2022).Cirrus clouds broadly extend over the tropics (e.g., Sassen et al. 2008), and their longwave radiative forcing dominates the earth's energy budget, even with their small optical thickness (e.g., Liou 1986).The dominant cloud microphysical processes in tropical cirrus clouds are ice nucleation, aggregation, vapor deposition/sublimation, and gravitational sedimentation (e.g., Seeley et al. 2019;Ohno et al. 2021).Double-moment bulk cloud microphysics schemes potentially capture the rapid growth of ice particles by aggregation and the corresponding increase in ice terminal velocity, which depends mainly on the maximum dimension.
The present study focused on the uncertainties related to the simplification of the ice terminal velocity formulation, because modification of this formulation results in a strong change in the global radiation budget (e.g., Mitchell et al. 2008;Hourdin et al. 2017).In general, the ice terminal velocity υ t in bulk cloud microphysics schemes is approximated by a power law relationship to the maximum dimension D, as follows: Here, the coefficient a υ and exponent b υ are generally derived by fitting to observations (e.g., Locatelli and Hobbs 1974;Heymsfield and Kajikawa 1987).Now, two issues are involved with this formulation: a υ and b υ are given as a global constant; and a υ and b υ are derived by fitting in a narrow size range.Figure 15 shows the dependence of the ice terminal velocity υ t on the maximum dimension D based on the theoretical formulation (Böhm 1989;Mitchell 1996;Seiki and Ohno 2022).It is clear that the exponent varies by the size range: b υ approximates 2 at smaller sizes and approaches 0.5 at larger sizes.Therefore, use of a single pair of a υ and b υ results in systematic biases in cloud microphysical processes in global simulations.The issues in gravitational sedimentation in the NDW6 scheme were already dealt with by Seiki et al. (2014) by using additional pairs of a υ and b υ for different size ranges.However, the issues in collisional growth in the NDW6 scheme have not yet been solved.Seifert et al. (2014) showed that use of Eq. ( 4) causes non-negligible (4) errors in the curve of collisional growth.Therefore, Seiki and Ohno (2022) revised how collisional growth is determined and then examined the impact of this revision on the simulation of tropical cirrus clouds.Note that heterogeneous and homogeneous ice nucleation were also revised in the new NDW6 scheme.
The collisional growth terms are evaluated by integrating the collection kernel, which consists of the product of the collisional cross-section and the difference in terminal velocity.However, it is difficult to evaluate the difference in terminal velocity in bulk cloud microphysics schemes (e.g., Seifert et al. 2014;Karrer et al. 2021).As a result, in the original NDW6 scheme, the collisional cross-section and the difference in terminal velocity were integrated separately.In addition, the power law relationship was used for the ice terminal velocities.In the revised scheme, we numerically integrated the collection kernel with the Gauss-Legendre quadrature.Thus, the difference in terminal velocity can be directly integrated with the theoretical formulation of ice terminal velocities.
We found that the original scheme overestimated aggregation between cloud ice to form snow (autoconversion of cloud ice) by approximately 300-400% at sizes smaller than 20 μm (not shown).This error corresponds to the bias in the terminal velocity due to the power law relationship (cf. Figure 15).Similarly, accretion of rain on graupel was overestimated in the original scheme (not shown).Therefore, the original scheme causes severe errors in thin cirrus clouds and intense rainfall systems.
The lifetime of tropical cirrus clouds was examined in reference to the CloudSat satellite observations (Fig. 16).Here, this study only analyzed cirrus clouds over the ocean, since cirrus clouds over orographic features are strongly affected by atmospheric disturbances (Seiki et al. 2019).With the new NDW6 scheme, aggregation of cloud ice and snow significantly decreases, particularly at altitudes above 8 km (not shown).Correspondingly, the frequency of radar echoes larger than − 20 dBZ near the tropopause (14-16 km) decreases in the global simulations with the new NDW6 scheme (Fig. 16b-c).As a result, slowly growing ice crystals as a result of weak collision and vapor deposition in the new NDW6 scheme maintain thin cirrus clouds in the upper troposphere.This signal is represented as the distinct mode value of the radar echoes of -30 to -20 dBZ at altitudes above 12 km.In addition, increasing radar echoes toward the cirrus cloud base (8-12 km) become more distinct when using the new NDW6 scheme.All the revisions increase the cirrus cloud fraction over the tropics (Fig. 16d), and consequently, shortwave cloud radiative forcing and longwave cloud radiative forcing over the tropics improve by approximately 4.1 and 7.5 W m −2 , respectively, compared to the Clouds and Earth's Radiant Energy System (CERES) satellite observations.These findings indicate that the lifetime of cirrus clouds is reasonably represented by the revised NDW6 scheme.on the maximum dimension D (m).In the database compiled by Mitchell (1996), cloud ice (blue) is assumed to form hexagonal columns, snow (green) is assumed to comprise assemblages of planar polycrystals in cirrus clouds, and graupel (red) is assumed to be lump graupel.From Seiki and Ohno (2022), © American Meteorological Society.Used with permission

Consideration of SGS ice condensation in a turbulent closure scheme
Turbulent processes play important roles in the life cycle of various types of clouds (e.g., Squires 1958;Klaassen and Clark 1985;Grabowski 1993Grabowski , 2007;;Grabowski and Clark 1993;Gasparini et al. 2019); as such, the representation of moist processes in turbulent schemes is crucial for their performance.In general, turbulent schemes assume that the time scale of SGS cloud condensation is sufficiently shorter than that of turbulent mixing, and they employ saturation adjustment-type approaches (e.g., Olson et al. 2019).However, the time scale of ice condensation is several orders of magnitude longer than that of liquid water, and the phase relaxation time for ice clouds under typical conditions is several orders of magnitude longer than the time-step length commonly used in the past relevant numerical studies (Heymsfield and Miloshevich 1995;Khvorostyanov and Curry 2014;Gryspeerdt et al. 2018).These studies speculate that the use of a saturation adjustment-type approach for representing SGS ice clouds in turbulent schemes overestimates the effects of the phase change on the turbulent mixing.This background knowledge motivated us to reconsider the representation of ice phase clouds in a turbulent closure scheme and to evaluate the effects of this scheme on the high-cloud cover in response to global warming.
The effects of the saturation adjustment-type approach for simulating SGS ice clouds in turbulent schemes on high clouds and their response to a warmer climate were investigated by a sensitivity study based on RCE simulations (Ohno et al. 2020).Figure 17 shows the domainaveraged high-cloud cover at the equilibrium state in the RCE simulations with and without an SGS condensation scheme for ice water condensate and with SSTs of 300 and 304 K.The suppression of SGS ice condensation reduced the high-cloud cover and altered the sign of the cloud cover response to the SST change, results that are similar to the effects of reducing the turbulent diffusivity K by reducing the turbulent mixing length (Ohno et al. 2019).
Next, the effects of SGS ice clouds on K were investigated.Figure 18a and b shows binned vertical profiles of K calculated with and without the SGS ice condensation  (Rossow and Schiffer 1999).From Ohno et al. (2020), © The Meteorological Society of Japan.Used with permission scheme sorted by the ice water path, respectively.K and ice water path were calculated using a snapshot dataset from the ICE simulation with an SST of 300 K.The values of K were large at the convective core region and in the vicinity of the cloud top in both cases due to the frequently occurring static instability (not shown).The static instability was enhanced by the SGS ice condensation scheme.Consequently, the magnitude of K in the upper troposphere was generally larger with the SGS ice condensation scheme than without it.Similar results can be seen with an SST of 304 K, as shown in Fig. 18c and d.The impacts of suppressing the SGS ice condensation scheme for the turbulence scheme on K were consistent with those of reducing the turbulent mixing length in the study of Ohno et al. (2019).These results indicate that the application of the SGS ice condensation scheme for the turbulence scheme changed high cloud covers and their response to a change in SST by altering the static stability in the cloud layers.
Since the phase relaxation time of ice clouds is much longer than that of liquid clouds, the effect of the phase change in ice clouds on the dynamical fields through the buoyancy should be smaller than those of liquid phase clouds.This suggests that the application of a saturation adjustment-type approach for representing SGS ice clouds in the turbulent scheme overestimates the turbulent diffusivity and causes model biases in the high-cloud fields.

Impacts of ice hydrometeors on a radiative field
Optical characteristics differ depending on the types of ice hydrometeors.However, most GCMs simplify the treatment of those ice species, which often causes errors in modeled radiative fields.It is also an interesting application of a high-resolution GCM to evaluate the extent to which the simplified cloud modeling leads to model biases, knowledge of which can hint at the source of errors in recent GCM simulations.Here, we show how Fig. 18 Binned vertical profiles of the turbulent diffusivity K (color and white lines) and the ice condensate (black lines) calculated a with and b without the SGS ice condensation scheme in the turbulent closure scheme, sorted by the ice water path.K and ice water path were calculated using a snapshot dataset from the ICE simulation with an SST of 300 K. c and d are the same as (a) and (b), respectively, for the simulation with an SST of 304 K. From Ohno et al. (2020) the ice hydrometeors, such as cloud ice, snow and graupel, that are simulated in NICAM, may affect the simulated field of longwave radiation.We use 3-month-long simulations during boreal summer from 1 June to 31 August 2004 with 14-km mesh NICAM using a doublemoment bulk cloud microphysics scheme (Seki and Nakajima 2014).To estimate the radiative effect attributed to each ice hydrometeor, the radiation transfer model MSTRNX (Sekiguchi and Nakajima 2008), which is the same radiative code as implemented in NICAM, is run offline.We designed four experiments, as follows: 'ctrl' represents the longwave radiative effects attributed to all ice hydrometeors (i.e., cloud ice, snow, and graupel); 'no_s' represents all ice hydrometeors except snow; 'no_gs' represents all ice hydrometeors except snow and graupel; and in 'allice' , snow is replaced by cloud ice to estimate the maximum longwave cloud radiative forcing.
The effect on the outgoing longwave radiation (OLR) for each experiment is plotted in Fig. 19a-c, and the effect on the cloud optical depth (COD) is plotted in Fig. 19d-f.COD is underestimated and OLR is overestimated in 'no_s' (Fig. 19a and d) and 'no_sg' (Fig. 19b  and e).This effect is prominent in the tropics and in the storm track regions at mid-latitudes, where cloud ice and snow are abundant.The positive bias in OLR due to the removal of snow leads to an average bias of 1 W m -2 and reaches a maximum bias of 2 W m −2 in the Indian Ocean region (Fig. 19a).The effect on the OLR field attributed to graupel is negligible, as seen by the effect on the OLR fields in 'no_s' and 'no_sg' , although the effect on the COD field attributed to graupel is non-negligible (Fig. 19d and e).By contrast, OLR is underestimated over the tropics and mid-latitudes when snow is replaced by cloud ice (Fig. 19c), due to the overestimation of COD (Fig. 19f ).The underestimation of OLR reaches − 0.4 W m −2 over the Indian Ocean (-0.2 W m −2 on average over the intertropical convergence zone).The horizontal distribution of these effects of snow on OLR is similar to that estimated from CloudSat observations (Waliser et al. 2011) and GCMs (Li et al. 2016).However, the magnitude of the effect is about half that reported in previous studies.This difference is due mainly to the vertical distribution of ice hydrometeors (Chen et al. 2018).

Improvement of bulk microphysical scheme by bin scheme
We also improved the two-moment bulk microphysics scheme based on a bin microphysics scheme to improve the simulation of water clouds in NICAM.Kuba et al. ( 2020) compared the two-moment bulk scheme (NDW6, based on Seifert and Beheng (2006) and modified slightly by Seiki and Nakajima (2014)) and the two-moment bin scheme (developed by Kuba and Fujiyoshi (2006) and modified by Kuba and Murakami (2010)).Their aim was to constrain the parameters in cloud bulk schemes by using observational data from satellite remote sensing.They used a dynamic-kinematic model to avoid the interactions between cloud microphysics processes and dynamics.Kuba et al. (2020) also used the Joint Simulator for Satellite Sensors (hereafter referred to as the Joint-Simulator; Hashino et al. 2013;Satoh et al. 2016) to calculate the horizontally averaged radar reflectivity Z m and the optical depth from the cloud top τ d .They studied the conversion processes from cloud droplets to raindrops in shallow cumulus clouds.In their study, they conducted sensitivity experiments on vertical velocity, the concentration of cloud condensate nuclei, and size distribution parameters and studied the relationships between τ d and Z m .
In Kuba et al. (2020), the size distributions of cloud droplets and raindrops for the bulk scheme are represented by generalized gamma distributions, as follows: where a = c for cloud droplets, a = r for raindrops, and x is the mass of a cloud droplet or a raindrop.
Figure 20 shows the relationship between the domainaveraged accumulated surface rainfall at 120 min and the cloud droplet number concentration near the cloud base at about 15 min using the bulk and bin schemes for thin, medium-thickness, and thick clouds.For the bulk scheme, the values of (v c , v r ) in the generalized gamma distributions were varied in Fig. 20 (a: (1, − 1/3), b: (5) f a (x) = α a x υ a exp − a x µ a , (− 1/3, − 1/3), c: (1, 1), d: (− 1/3, 1).The values of (1, − 1/3) for (v c , v r ) are based on Seifert and Beheng (2006) and Seiki and Nakajima (2014).The rainfall amount at 120 min was similar between the bulk and bin schemes for the thick clouds (green circles in Fig. 20a − d) and the smaller numbers of cloud droplets.Conversely, the rainfall amount was smaller in the bulk scheme than in the bin scheme for clouds with a large number of cloud droplets, particularly for the thin clouds (shown as black circles whose cloud droplet number concentrations are 160 cm −3 and larger).In the case of (v c , v r ) values of (− 1/3, 1), the difference between the bulk and bin schemes is small (Fig. 20d).Decreasing the shape parameter v c from 1 to − 1/3 leads to an increase in the auto-conversion rate, and an increase in ν r from − 1/3 to 1 decreases the width of the raindrop size distribution (which means a decrease in the falling velocity of raindrops).In the bulk scheme, a decrease in the falling velocity of raindrops increases the falling time (which means an increase in rain production).Kuba et al. (2020) also compared bin and bulk schemes with (v c , v r ) of (1, − 1/3) and (− 1/3, 1) using the relationships between horizontally averaged radar reflectivity Z m and optical depth from the cloud top τ d .Figure 21 shows that the optical thicknesses of droplets with a radius larger than 8 μm were similar to those of all droplets in the case of the bin scheme (see Fig. 21a and d).This means that almost all droplets simulated using the bin scheme had a radius larger than 8 μm near the cloud top.However, with the bulk scheme, the optical thicknesses of droplets with a radius larger than 8 μm were smaller  Kuba et al. (2020) than those of all droplets (see Fig. 21b and e, c and f ).The results of the bulk scheme using (v c , v r ) of (-1/3, 1) and excluding raindrops with a radius less than 8 μm (Fig. 21f ) are most similar to those of the bin scheme (Fig. 21a and d).Above, we suggested a set of improved parameters in a bulk microphysics scheme to evaluate the growth of cloud droplets to rain according to the degree of atmospheric pollution.In future, the bulk scheme will need to be improved for other meteorological conditions using the method described in Kuba et al. (2020).

CMIP6 HighResMIP simulations
Before the TOUGOU program, NICAM simulation data had been analyzed mostly by developers and users of the NICAM code.Recent advances in computer technology now allow us to perform much longer-term NICAM simulations.In this context, Haarsma et al. (2016) proposed the HighResMIP to endorse model intercomparison projects in CMIP6, which we consider to be a good opportunity not only to contribute to the CMIP community but also to try NICAM simulations at scales of more than half a century.
Practical climate simulations using NICAM originated from NICAM AMIP-type simulations (Kodama et al. 2015), in which NICAM.12 with a mesh size of 14 km was used.Here, NICAM adopts a release number (e.g., "12") based on the last 2 digits of the year when major updates from the previous release version are completed.NICAM.12 successfully simulated a wide variety of phenomena, particularly the distribution and seasonal march of tropical cyclogenesis, and the simulation data were intensively used for an analysis of tropical cyclones (Satoh et al. 2015;Yamada et al. 2017;Matsuoka et al. 2018;Sugi et al. 2020), among others.However, major issues remained in the simulation of basic-state climatology such as surface air temperature, cloud, and precipitation, which motivated us to develop NICAM.16 and its CMIP6 version, NICAM.16-S("-S" represents the use of a single-moment cloud microphysics scheme), as reported in Kodama et al. (2021).In summary, the update of the single-moment bulk cloud microphysics scheme (Roh and Satoh 2014;Roh et al. 2017) affected the vertical distributions of snow and cloud ice, leading to better simulations of the amount of high thin cloud.In addition, the improved treatment of natural and anthropogenic aerosols; updates to the land surface model, surface albedo, and sea ice thickness; and the introduction of an orographic gravity wave drag scheme together contributed to an improvement in the simulation of climatology.
NICAM.16-S was used to simulate timeframes of a century with 56 and 28 km mesh and a decade with 14 km mesh to produce a dataset for CMIP6 HighResMIP.The simulated year per day (SYPD) of NICAM.16-S on the Earth Simulator 3 (NEC SX-ACE system) was 0.22-0.63 for 14-56 km mesh using 10-160 nodes and 40-640 MPI processes (Table 7 in Kodama et al. 2021), and it often took comparable time to wait for the job to start on the Earth Simulator 3. The time needed to perform postprocessing-such as remapping and adding metadata to meet the CMIP6 standard format-was comparable to that required for the simulation itself.The total data size of the final product was 115 TB for the 11-year 14 km mesh dataset and 261 TB for the 101-year 28 km mesh dataset.The HighResMIP dataset is now widely used in climate research (e.g., Roberts et al. 2020b;Yamada et al. 2021;Liang-Liang et al. 2022;Priestley and Catto 2022).

Future directions
Future changes in weather and climate extreme events in a changing climate with global warming are a great concern, and they are intensively assessed in Chapter 11 of the Sixth Assessment Report of the Intergovernmental Panel on Climate Change (IPCC 2021; Seneviratne et al. 2021).Among extreme events, the torrential rain and violent wind associated with mesoscale convective systems or TCs are represented more accurately by high-resolution models.Global storm-resolving models, or global kilometer-scale models, with a horizontal mesh size of O(km), are becoming popular.Their use for modeling future changes in extreme events is prospective (Satoh et al. 2019;Slingo et al. 2022).Studies of future changes in TCs require global high-resolution simulations that simultaneously reproduce large-scale circulations and the inner structure of TCs that consider the interaction of the atmosphere-ocean coupled system (Knutson et al. 2020).We studied future changes in TCs in Sect.2.1.4,in which the atmospheric 14-km mesh model was used.Ideally, the simulations for the TC projection should be with a global atmospheric-ocean coupled kilometer-scale model.In particular, to assess near-future changes in TCs in a world with a 1.5 °C or 2.0 °C warming level, a large number of decadal ensemble simulations with such kilometer-scale-type models are eagerly awaited.A highresolution version of Database for Policy Decision-Making for Future Climate Change or d4PDF (Ishii and Mori 2020) is also in high demand for more accurately assessing projected changes in TCs.Such experiments using global kilometer-scale models will allow us to project not only TCs (Sect.2.1.4)but also related processes of mesoscale convective systems, such as convective aggregation (Sect.and upper clouds (Sect.2.1.2).
Several scientific and technical challenges must be overcome to establish NICAM as a kilometer-scale climate model.How such a model should be configured for clouds in kilometer-scale climate simulations remains under debate.Although the kilometer-scale model resolves deep convection, albeit partially, it is obvious that its spatial resolution is completely insufficient for representing shallow convection and boundary-layer clouds.By and large, there are two possible approaches to achieve this: the explicit cloud microphysics (i.e., cloudresolving) approach, and the cloud parameterization (i.e., GCM) approach.In the former approach, vertical resolution and turbulent mixing may be key to dealing with the issue of lower-level clouds.For example, refining the vertical resolution only for a specific physics scheme (Yamaguchi et al. 2017) may be a promising approach.Flux adjustment, albeit classical, may also be an option for representing some scientific targets such as TCs.In the latter approach, additional empirical schemes such as a shallow convection scheme may be introduced.Another popular scheme is eddy diffusivity/mass flux or EDMF parameterization (e.g., Suselj et al. 2022), in which eddy diffusivity and mass flux formulations are used to represent downdraft and updraft regions, respectively.In both cases, it is necessary to consider the applicability of the schemes to a kilometer-scale grid box and also to consider the relative role of and consistency between parameterization schemes and microphysics schemes.
In terms of modeling technology, accelerating the computation is necessary to run the model at the climate scale.As previously stated in this review, the SYPD of the 14 km mesh NICAM on the Earth Simulator 3 is around 0.22.The model speeds up by a factor of around 1.8 if most of the double-precision floating-point arithmetic is replaced with single-precision (Nakano et al. 2018).Owing to the excellent weak-scaling performance of NICAM (Yashiro et al. 2016), the wall-clock time per model time-step is almost constant when the horizontal grid spacings are halved and the number of computation nodes is quadrupled.Note that the time-step interval should be halved to halve the horizontal grid spacings.According to a benchmark, our model runs a few times faster on the Fugaku supercomputer than on the Earth Simulator 3. Therefore, we expect an SYPD of around 0.3 for the 3.5 km mesh NICAM on the supercomputer Fugaku, which is comparable with that for the 5 km mesh coupled ICON model (Hohenegger et al. 2022).Other possible approaches for speeding up computation include reducing the required byte-per-flops with even lower floating-point arithmetic (e.g., Paxton et al. 2022), performing calculations faster with the aid of accelerators, and replacing physics schemes with AI-based surrogate models (Arakawa et al. 2022).In addition, implementing an analysis platform on the supercomputer or cloudbased system is urgently needed to reduce the time needed for data transfer.Future technological trends in climate simulation are discussed in detail in WMO (2021).

Conclusions
We have reviewed the results of global nonhydrostatic simulations using NICAM to examine future changes in clouds, which were conducted in the core Japanese research program for climate change, the TOUGOU program, in the five fiscal years since 2017.This review also describes current research activities around the world and future issues in high-resolution climate modeling.
For the projection of clouds, previous NICAM studies have predicted a larger number of high clouds, which is characterized by a larger number of smaller high clouds in a warmer atmosphere and then an increase in the high-cloud amount.However, it was not clear why high clouds change in such a way.The lack of knowledge about the underlying mechanisms has been a large motivation behind our climate research.Noda et al. (2019) reveals that a key factor was a reduced degree of disorganized tropical convection due to warming.In a warmer atmosphere, the smaller degree of cloud organization leads to a more scattered development of high clouds, which then change to cover a larger area of the tropical atmosphere.In fact, our results show that the large-scale circulation strongly correlates with the degree of cloud organization, suggesting that tropical convection obtains more vapor from a boundary layer.In addition, large-scale ascendant flows act to generate convection.Further efforts are needed to gather observational data in order to evaluate the modeled result.We also evaluated changes in the vertical structure of ice species in high clouds due to warming and their relationship to changes in surface temperature.
Using RCEs also enabled us to systematically examine processes regarding clouds by using a sensitivity analysis of the cloud processes, and we found a change in behavior of cloud microphysics processes to cloud ice in a warmer atmosphere.Conventionally, most previous studies focus on a change in cloud-top temperature due to change (e.g., Hartman and Larson 2002;Zelinka and Hartmann 2010); by contrast, less attention has been paid to the role of a change of cloud-top pressure.It is well recognized that the mean altitude of high clouds increases in a warmer atmosphere (e.g., Zelinka et al. 2013), which also means that cloud-top pressure decreases (Zelinka and Hartmann 2010).Our study first pointed out that reduced cloud-top pressure plays an important role in reducing the lifetime of anvil clouds by enhancing the deposition growth of cloud ice from vapor, thereby leading to an increase in sedimentation.
We also investigated changes in well-organized cloud disturbances, such as TCs and extratropical cyclones, due to warming.The modeled projection of a reduced number of TCs is attributed to a reduced number of TC seed; by contrast, the survival rate of TCs (defined as a ratio of the number of TCs to TC seeds) is almost identical.We also showed that TC size increases with latitude, on average, and its regional characteristics.We also identified a larger number of moderately intense extratropical cyclones in warming climate.
We also evaluated an effective climate sensitivity, which is the first attempt based on data from global nonhydrostatic climate simulations.In practical terms, the detailed value of climate sensitivity can depend on cloud schemes and their parameters (e.g., Kodama et al. 2012).Therefore, in the next step, we need to improve our knowledge of the extent to which microphysics parameters affect cloud feedback processes and eventually climate sensitivity, as well as how those parameters can be constrained based on observation data in future.Those efforts will lead to a better estimation of the climate sensitivity of our planet.
To improve the performance of future global nonhydrostatic simulations, we have also proposed methods to improve the physical processes in such global high-resolution models, including low-level mixed-phase clouds, cirrus, warm rain, and SGS turbulent mixing processes, along with evaluating the importance of taking precipitating categories of ice species into account.
Efforts to develop a new type of climate model that is capable of resolving storms, like NICAM, have been enthusiastically pursued in world climate modeling centers.A global nonhydrostatic model is a useful tool for climate research; however, it has many aspects that could be improved, as we have reviewed here.For example, further efforts are needed to improve the physical schemes regarding clouds and the spatial resolution to better simulate the behavior of clouds.For advanced climate research, it is important that these efforts continue and that research groups around the world collaborate and cooperate in sharing knowledge.

Fig. 1
Fig. 1 Schematic of the influence of convective aggregation on outgoing longwave radiation in cases where convection becomes (left) more aggregated and (right) less aggregated.Wavy upward arrows indicate longwave radiation.Longwave radiations with thick arrows in darker color and thin arrows in lighter color show those emitted from the earth's surface and from high clouds, respectively

Fig. 3
Fig.3The relationship between < w > and SF' a for present (CTL; red dots) and global warming (GW; black dots) simulations, and b their difference (i.e., red minus black).The correlation coefficients of < w > and SF' are shown above the panels.FromNoda et al. (2019)

Fig. 4 a
Fig. 4 a Vertical distribution of cloud fraction (January 1989-December 2008) regressed onto the SST for a NICAM AMIP-type experiment under present (aqua) and warming (pink) conditions.b Vertical distribution of ice hydrometeors regressed onto the SST in a NICAM AMIP-type experiment, where yellow, green, and blue denote ice water content from cloud ice, snow, and graupel, respectively.In b solid and dashed lines denote present and warming conditions, respectively.Error bars show the standard deviation in the tropical region (30° S-30° N) for each variable.From Chen et al. 2022, © American Meteorological Society.Used with permission

Fig. 5
Fig. 5 Domain-averaged vertical profiles of a the net, b positive, and c negative tendencies of ice clouds due to cloud microphysical processes at equilibrium states in radiative-convective equilibrium simulations with a SST of 300 or 304 K. From Ohno et al. (2021), © The American Geophysical Union.Used with permission

Fig. 6 aFig. 7
Fig. 6 a Domain-averaged high-cloud cover and b its response to a SST change in control simulations (CTRL) and simulations with a lower limit of pressure for the evaluation of the terminal velocity of ice particles (VTp440) and diffusivity of water vapour in air (KVp440).Differences in c high-cloud cover and d the responses between the control (CTRL) and sensitivity (SENS) simulations.From Ohno et al. (2021)

Fig. 8
Fig. 8 Horizontal distribution of the mean 8 m s −1 radius (R08) on a 5º × 5º grid box for a the control experiment (CTL) and b the global warming experiment (GW), and c the difference between them (GW minus CTL).Panel d indicates the zonal mean of R08 in CTL (solid line) and GW (dashed line) for each ocean basin: the North Atlantic (70ºW-20ºW), the eastern North Pacific (180º-120ºW), the western North Pacific (140ºE-180º), the South Pacific (160ºW-120ºE), and the South Indian Ocean (50ºE-110ºE)

Fig. 9
Fig. 9 Monthly frequency of all the simulated oceanic extratropical cyclones in the northern hemisphere binned by lifetime-maximum 850 hPa wind speed (5 m s −1 bin size).a Results for 1951-1960 (past) simulated by the 56 km (blue), 28 km (green), and 14 km (red) mesh model.b-d Changes in 2001-2010 (present-day; green line) and 2041-2050 (future; red line) with reference to the past.The mesh sizes of the model 56 km (b), 28 km (c), and 14 km (d), respectively.The gray shadings in (b)-(d) show the results for the past divided by 10

Fig. 10
Fig. 10 Composite of future minus past precipitation amount averaged over all the oceanic extratropical cyclones (top row) and the intense oceanic extratropical cyclones (bottom row) for (left to right): 56 km, 28 km, and 14 km mesh models.Gray circles indicate the distance from the center of the extratropical cyclone in intervals of 500 km

Fig. 12
Fig. 12Time series of mixing ratios of cloud water (qc), rain (qr), cloud ice (qi), snow (qs), and graupel (qg) from the single-column simulation with the NDW6 scheme (a-e), the original (OLD) NSW6 scheme (f-j), and the revised (NEW) NSW6 scheme (k-o).The freezing level (z ~ 220 m) is indicated by the solid black line.The units are g kg -1 .FromSeiki and Roh (2020), © American Meteorological Society.Used with permission

Fig. 14
Fig. 14 Fields of shortwave radiation at the top of the atmosphere for a observation, and one-year means of results from experiments using b the old scheme, and c the new scheme.d-f are the same as (a), (b), and (c), but for the liquid water path, respectively.Climatology of the CERES data and that of the Observation data in (a) and (d) are Multisensor Advanced Climatology Mean Liquid Water Path (MAC-LWP) data are shown in (a) and (d), respectively.From Noda et al. (2021)

Fig. 15
Fig.15Dependence of the terminal velocity υ t (m s −1 ) on the maximum dimension D (m).In the database compiled byMitchell (1996), cloud ice (blue) is assumed to form hexagonal columns, snow (green) is assumed to comprise assemblages of planar polycrystals in cirrus clouds, and graupel (red) is assumed to be lump graupel.FromSeiki and Ohno (2022), © American Meteorological Society.Used with permission

Fig. 16
Fig. 16 Diagram of the contoured frequency of the 94 GHz radar echo by altitude from a CloudSat satellite observations, b global simulations using the original version of the NDW6 scheme, and c those using the revised version of the NDW6 scheme.The vertical profiles of the cirrus cloud fraction are also indicated (d).From Seiki and Ohno (2022), © American Meteorological Society.Used with permission

Fig. 17
Fig. 17 Domain-averaged cloud cover for a total, b thin, c medium, and d thick high clouds for the simulations with and without an SGS condensation scheme for ice water condensate (labeled as ICE and NOICE, respectively) using SSTs of 300 (black) and 304 (red) K. e Cloud cover response to increasing SST for the total (purple), thin (green), medium (blue), and thick (orange) high clouds.High clouds were defined based on the International Satellite Cloud Climatology Project definition of cloud types (Rossow and Schiffer 1999).From Ohno et al. (2020), © The Meteorological Society of Japan.Used with permission

Fig. 21
Fig. 21 Relationships between horizontally averaged radar reflectivity (Z m ) and optical depth from the cloud top (τ d ) for the case of thin clouds with polluted cloud condensate nuclei.Results are shown for a the bin scheme, b the bulk scheme with (ν c , ν r ) = (1, − 1/3), c the bulk scheme with (ν c , ν r ) = (− 1/3, 1).d-f are the same as a-c, respectively, but excluding droplets with a radius less than 8 μm from the calculation performed by the Joint Simulator.The flow of time is drawn in black, red, green and blue.FromKuba et al. (2020) This study was supported by the Integrated Research Program for Advancing Climate Models (TOUGOU) Grants JPMXD0717935457 and the Advanced Studies of Climate Change projection (SENTAN) JPMXD0722680395 from the Japanese Ministry of Education, Culture, Sports, Science and Technology.CK was also supported by the Environment Research and Technology Development Fund (JP-MEERF20172R01) of the Environmental Restoration and Conservation Agency of Japan (ERCA), JSPS KAKENHI Grant Number JP17H04856, and JSPS KAKENHI Grant Number JP20H05728.
. Since then, NICAM has been used in several large Japanese climate research projects: the Innovative Program of Climate Change Projection for the 21st Century, or KAKUSHIN program in Japanese, from 2007 to 2012; the Program for Risk Information on Climate Change, or SOUSEI in Japanese, from 2012 to 2017; and the Integrated Research Program for Advancing Climate Models, or TOUGOU in Japanese, from 2017 to 2022.