The impacts of climate change on the winter water cycle of the western Himalaya

Abstract

Some 180 million people depend on the Indus River as a key water resource, fed largely by precipitation falling over the western Himalaya. However, the projected response of western Himalayan precipitation to climate change is currently not well constrained: CMIP5 GCMs project a reduced frequency and vorticity of synoptic-scale systems impacting the area, but such systems would exist in a considerably moister atmosphere. In this study, a convection-permitting (4 km horizontal resolution) setup of the Weather Research and Forecasting (WRF) model is used to examine 40 cases of these synoptic-scale systems, known as western disturbances (WDs), as they interact with the western Himalaya. In addition to a present-day control run, three experiments are performed by perturbing the boundary and initial conditions to reflect pre-industrial, RCP4.5 and RCP8.5 background climates respectively. It is found that in spite of the weakening intensity of WDs, net precipitation associated with them in future climate scenarios increases significantly; conversely there is no net change in precipitation between the pre-industrial and control experiments despite a significant conversion of snowfall in the pre-industrial experiment to rainfall in the control experiment, consistent with the changes seen in historical observations. This shift from snowfall to rainfall has profound consequences on water resource management in the Indus Valley, where irrigation is dependent on spring meltwater. Flux decomposition shows that the increase in future precipitation follows directly from the projected moistening of the tropical atmosphere (which increases the moisture flux incident on the western Himalaya by 28%) overpowering the weakened dynamics (which decreases it by 20%). Changes to extreme rainfall events are also examined: it is found that such events may increase significantly in frequency in both future scenarios examined. Two-hour maxima rainfall events that currently occur in 1-in-8 WDs are projected to increase tenfold in frequency in the RCP8.5 scenario; more prolonged (1-week maxima) events are projected to increase fiftyfold.

1 Introduction

1.1 Western disturbances

Western disturbances (WDs; Pisharoty and Desai 1956; Rao and Srinivasan 1969) are mid-latitude cyclonic storms, typically embedded in the subtropical westerly jet stream, that are capable of bringing heavy precipitation to the mountain ranges of South Asia. They are most frequent during the boreal winter (Dec–Mar) during which they are responsible for the majority of precipitation over Pakistan and north India (Dimri 2013), and almost all instances of extreme precipitation (Hunt et al. 2018b). WDs are extratropical in nature (Dhar et al. 1984) and usually originate as dynamical instabilities in the jet over the Mediterranean or western Asia (Cannon et al. 2016; Hunt et al. 2018a) before intensifying in its baroclinically unstable environment (Puranik and Karekar 2009). On subsequent approach to the subcontinent, WDs draw moisture from the Arabian Sea (from isotope modelling: Krishan et al. 2012; from dynamical considerations: Rangachary and Bandyopadhyay 1987; Hatwar et al. 2005; Hunt et al. 2018b) and advect it towards the nearby orography: the Hindu Kush, the Karakoram, and the western Himalaya. Given that about two-thirds of the annual discharge of the Indus River comes from snowmelt in this region (Bookhagen and Burbank 2010) and that snow deposition here is principally a consequence of WD activity, it is clearly important to investigate how that relationship changes under future climate forcing.

1.2 Observed trends of WDs

A number of studies have assessed trends of western disturbance frequency in station-based observational records, all finding either no significant trend or a significant negative trend. Das et al. (2002) derived the frequency of Mar–May WDs in the western Himalaya from 1971 to 2000, finding a significant decline only in May. Shekhar et al. (2010), also looking at the western Himalaya, found no significant trend in either Nov–Apr or Jan–Mar WD frequency over the period 1984–2007, and noted that this contrasted with the significant decline in snowfall days over the same periods. However, Kumar et al. (2015) noted a significant decrease in Dec–Mar WD frequency over the Central Himalaya from 1977 to 2007. Midhuna et al. (2020) collated seasonal WD frequency from daily weather reports issued by the India Meteorological Department from 1987 to 2016, for which there was a significant, declining trend. In each of Shekhar et al. (2010), Kumar et al. (2015), and Midhuna et al. (2020), the observed trend in frequency was a decline of about 1% per year. Other authors have used tracking algorithms to identify WDs in reanalysis datasets; both Cannon et al. (2016), using CFSR, and Hunt et al. (2018a), using ERA-Interim, identified insignificant (though declining) trends in winter western disturbance frequency.

Interestingly, authors that have tried to quantify WD frequency trends using proxies tend to disagree with those that use a more direct method. Madhura et al. (2015) used a high-pass filter on 500 hPa geopotential height before using EOF analysis to compute winter WD-like activity in 1948–2011 NCEP/NCAR reanalysis data, and found a significant increasing trend over this period. Cannon et al. (2015) performed wavelet analysis on 200 hPa geopotential height to investigate the interannual variability of Dec–Mar WDs. They found that the 5–15 day bandpassed synoptic power of this field had a significant positive trend over the western Himalaya, but a significant negative trend over the Central Himalaya, in CFSR reanalysis (1979–2010). Krishnan et al. (2018) used the standard deviation of 4–15 day bandpassed 200 hPa geopotential height over upstream Iran as a proxy for Dec–Apr WD activity. They found that this field exhibits a significant positive trend in both ERA-40 (1958–2002) and ERA-Interim (1979–2015) reanalyses. These disagreements with more direct counting methods are in part explained by Midhuna et al. (2020) who tried to derive an index of WD activity based on 850 hPa and 200 hPa geopotential height EOFs, but could only achieve a correlation coefficient of 0.18 between that and observed counts, suggesting that even well thought out proxy techniques are currently inadequate.

1.3 Observed trends of precipitation

One key motivation for investigating western disturbances is their central role in winter precipitation over these mountain ranges, as well as over Pakistan and north India in general. Many studies have attempted to determine trends in observed precipitation over this region, and are generally sensitive to the location and methodology used, as well as the altitudes of the stations in question. Indeed, stations in this region are not only sparse but have a low-altitude bias (Palazzi et al. 2013), and snow measurements at all altitudes are typically fairly error-prone (Winiger et al. 2005).

Nevertheless, some insights can be gained as long as we acknowledge their shortcomings. For example, Gusain et al. (2014) stratified 1976–2011 station data by altitude and showed that 19 of the 26 stations below 4000 m recorded declining trends in seasonal snowfall (though only 4 significantly so) compared to only 5 (only 2 significant) of the 11 above. Similarly, Singh et al. (2015) examined data from 11 stations over varying periods (most commonly 1987–2009), with all bar one below 3250 m; eight exhibited a decreasing trend in seasonal snowfall, four significantly so; in contrast, eight had an increasing trend in seasonal rainfall.

1.4 Projected trends of temperature and precipitation

Many studies have looked at the projected impacts of climate change on the Himalaya and Hindu Kush regions. Palazzi et al. (2015) found that CMIP5 models generally projected wetter summers over the region in both RCP4.5 and RCP8.5 scenarios, but could not agree on the fate of winter precipitation, further noting that there was appreciable intermodel variability in the representation of precipitation seasonality. In contrast, Panday et al. (2015) found a small but significant increase in winter precipitation over the western Himalaya in the CMIP5 RCP8.5 scenario, finding in addition that extreme winter precipitation events there are projected to become less frequent but more intense. Su et al. (2016) noted substantial spatial variability in the precipitation response when they looked at the output from 21 CMIP5 GCMs over the Indus Basin, but found that in both RCP4.5 and RCP8.5 scenarios, precipitation in the south of the basin (which typically falls as rain) was projected to significantly decrease whereas precipitation in the north of the basin (typically snow) was projected to significantly increase. Generally speaking, CMIP5 models agree on future changes to temperature and humidity, but still retain substantial uncertainty in precipitation; this is a core motivation of this study.

Some authors have also examined the specific response of the cryosphere in more detail. Terzago et al. (2014) found that CMIP5 models projected a robust and significant decline in winter snowpack depth over the Hindu-Kush Karakoram-Himalaya region, but hypothesised that warming, rather than reduced snowfall, was the likely cause. Immerzeel et al. (2013) forced a very high resolution (90 m) glacio-hydrological model with output from four CMIP5-class GCMs, finding robust increases in total annual (but not seasonal) runoff feeding both the Indus and Ganges, largely as a result of increased precipitation in both regions.

1.5 Projected trends of western disturbances

However, only a handful of studies have tried to assess the fate of WDs in future climate projections. Ridley et al. (2013) used a clustering technique to identify four weather regimes that were associated with anomalous snowfall over the Himalaya in a regional climate model. They performed two RCP4.5 experiments, differing only in the GCM used to force the boundary conditions, finding that one showed a significant increase in WD-like regimes, and the other did not. Krishnan et al. (2018) found, using the same bandpassed geopotential variance proxy discussed earlier, that this activity metric was projected to increase in an RCP4.5-forced regional projection. They attributed this to an enhanced upper-level quasi-stationary trough over central Eurasia which occurs as a result of differential heating over the Tibetan Plateau. Since we have seen (Sect. 1.2) that proxy techniques are not reliable when it comes to counting WDs, we cannot be certain of the significance of these results in the context of changing WD frequency in a future climate. However, Hunt et al. (2019) used a feature-tracking algorithm to detect WDs in the output of 24 CMIP5 models. Compared to a 1980–2005 baseline, they found that 2080–2100 WD frequency decreased significantly by an average 11% and 17% in RCP4.5 and RCP8.5 experiments respectively. In the latter, the mean WD intensity also fell by around 10%, which they attributed to a decrease in baroclinic instability in the environment of the subtropical jet. They also noted that models could not agree on the sign of change of extreme precipitation events associated with WDs.

1.6 Previous WRF simulations

How well are the Himalayan region and western disturbances actually simulated in WRF, the model to be used here? Norris et al. (2017) used a 6.7 km, convection-permitting WRF setup to generate a one-year dataset over the whole Himalaya. They found that precipitation patterns were well simulated compared to observations, including the diurnal and seasonal cycles and elevation dependence. Further marginal gains were found on increasing the resolution to 2.2 km; both simulations performed markedly better than a 20 km simulation and the 10 km High Asia Refined analysis (Maussion et al. 2014), both of which had parameterised convection.

The performance of WRF in simulating WDs has been the subject of a large number of studies, including a number (Rakesh et al. 2009a, b; Raju et al. 2015, and others) that have assessed the value of data assimilation in WRF simulations of WDs, though for the sake of brevity we will not expand upon those studies here. Thomas et al. (2014) used a 15 km resolution experiment with five cases to explore sensitivity to three land surface schemes. They found that mean precipitation was generally well captured even at this relatively coarse resolution, and that the dynamics of the WD were largely insensitive to the choice of land surface scheme. Dimri and Chevuturi (2014) used a 9 km resolution experiment with three cases to explore sensitivity to five microphysics parameterisations. They found that the the so-called Eta scheme (Rogers et al. 2001) produced results closest to observations for both dynamics and precipitation, but that the forecasts did not differ significantly. Patil and Kumar (2016) also used a 9 km resolution experiment to explore sensitivity to five microphysics parameterisations, though they used four cases. They found all five schemes resulted in accurate precipitation footprints, and concluded that a simple, single-moment scheme performed best with their setup. Norris et al. (2015) found that WD-induced heavy snow events were well simulated in a convection-permitting 6.7 km setup; and also found that this resolution was generally sufficient to capture the complexity of orographically-driven processes involved in Himalayan precipitation.

More generally, authors have argued that high-resolution models are required for accurate water cycle simulations over orography. Higher resolution leads to steeper features and hence greater uplift, enhancing precipitation; it also increases the height of the topography, increasing the portion of the surface above the freezing level and hence increasing the snowfall (e.g. Boberg et al. 2010; Rauscher et al. 2010). Garvert et al. (2007) reasoned through a combined observation and modelling study that a horizontal resolution of 4 km or less was needed to capture and simulate precipitating processes over complex orography due to the interaction with valley/ridge updrafts and downdrafts; this was further corroborated by Roberts et al. (2009). However, Kendon et al. (2012) showed that for certain cases over the United Kingdom, even a horizontal resolution of 1.5 km was too coarse to simulate realistic orographic updrafts.

1.7 Summary

So, we have seen that GCMs almost invariably have a horizontal resolution too coarse to capture the key orographic processes over the western Himalaya, resulting in poor representation of both climatological and extreme precipitation. This, in turn, engenders divergent estimates of precipitation trends, despite observational records showing a decline in snowfall. The fate of western disturbances—which bring the majority of winter precipitation to the region—in a changing climate is also unclear, though largely due to inconsistent methodologies; though feature-tracking results from CMIP5 models suggest a projected decline in both frequency and intensity. A crucial, unresolved question is how the tug-of-war between weaker synoptic scale forcing and the increased moisture in a future climate will impact on orographic precipitation—either on the mean climatology or the extremes. Since we cannot satisfactorily constrain this by interrogating CMIP5 output, we will use a suite of high-resolution, convection-permitting WRF experiments to examine the relationship between WDs and future climate scenarios. We will take 40 existing WD case studies from previous literature and simulate their impact on the water cycle over the western Himalaya lying within the Indus River basin. We will then perturb the boundary conditions to reflect both future and past climate scenarios and determine the resulting changes to precipitation before isolating the large- and small-scale processes responsible for the differences.

This study is laid out as follows: data and methods are described in Sect. 2, with data sources outlined in Sect. 2.1, case study selection in Sect. 2.2.1, and model setup in Sect. 2.2.2. Then we discuss the results in Sect. 3, with changes to mean precipitation in Sect. 3.1, the role of elevation in Sect. 3.3, the changes to extreme precipitation in Sect. 3.4, and analysis of the moisture flux budget in Sect. 3.5.1. Finally we discuss the results and future avenues of research in Sect. 4.

Table 1 Details of the 40 western disturbance case studies used in this paper: the date given is when the vortex crossed $70^\circ \,\hbox {E}$ for the last time

2 Data and methodology

2.1 Data

2.1.1 CMIP5

For this study, we use the 32 freely-accessible CMIP5 models (Taylor et al. 2012) for which monthly pressure level data were available in the four required experiments (see following paragraph). Where possible, the r1p1i1 ensemble member was chosen as the representative of each model, so as not to unfairly weight the results towards any particular model. The exception was EC-EARTH, for which, due to data availability reasons, member r9p1i1 was used.

In this study, we use data from four of the CMIP5 experiments: historical, pre-industrial, RCP4.5, and RCP8.5. The historical experiments of all models used here are forced with observed natural and anthropogenic contributions, usually over the period 1850–2005, from which we take a representative period of 1981–2005, against which all perturbations are computed. The pre-industrial experiment comprises longer simulations with no variation in anthropogenic forcings; these have varying baseline periods depending on the model, so we take the representative period as being the last 25 years of the run. The two future scenarios used here, RCP4.5 and RCP8.5, correspond to effective net changes in radiative forcing in 2100 of $4.5\,\hbox {W}\,\hbox {m}^{-2}$ and $8.5\,\hbox {W}\,\hbox {m}^{-2}$, equivalent to roughly 650 ppm and 1370 ppm $\hbox {CO}_2$ respectively (Van Vuuren et al. 2011), or increases in global surface temperature of about $2\,^\circ \hbox {C}$ and $4\,^\circ \hbox {C}$. In both cases we choose the final 25 years (2076–2100) as our representative period.

2.1.2 ERA-Interim reanalysis

For the initial and lateral boundary conditions of our model, we use the European Centre for Medium-Range Weather Forecasts Interim reanalysis (ERA-I; Dee et al. 2011). The surface fields, as well as soil temperature and moisture at selected depths are used only to produce the initial conditions; atmospheric variables, which include wind, moisture and temperature fields are used to construct both initial and boundary conditions. All fields are available at 6-h intervals with a horizontal resolution of T255 ($\sim 78$ km at the equator), with the three-dimensional fields further distributed over 37 vertical levels spanning from the surface to 1 hPa. Data are assimilated into the forecasting system from a variety of sources, including satellites, ships, buoys, radiosondes, aircraft, and scatterometers. ERA-I outputs that derive purely from the model (i.e. are not analysed), for example precipitation and cloud cover, are not used in this study.

2.1.3 GPM-IMERG

For model verification, we use the gridded surface rainfall product IMERG (Integrated Multi-Satellite Retrievals for GPM; Huffman et al. 2015). This has global coverage at a half-hourly, $0.1^\circ$ resolution, starting June 2000 and continuing to the present day. Over the tropics, IMERG primarily ingests retrievals from (for 2000–2014) the now-defunct Tropical Rainfall Measuring Mission (TRMM; Kummerow et al. 1998, 2000) 13.8 GHz precipitation radar and microwave imager (Kozu et al. 2001) and (for 2014-) the Global Precipitation Measurement (GPM; Hou et al. 2014) Ka/Ku-band dual-frequency precipitation radar. Where an overpass is not available, precipitation is estimated by calibrating infrared measurements from geostationary satellites. While GPM-IMERG performs well when compared against gauge-based products, performance falls significantly at higher altitudes or when quantifying extreme rainfall events (Prakash et al. 2018).

2.2 Experimental setup

2.2.1 Case study selection

The nature of this study makes the use of a climatology impractical and the use of a small number of case studies unsuitable, thus we start by compiling a large number of case studies from previous literature, which are tabulated in Table 1. These were subject to two criteria: that the WD must have a clear signal in relative vorticity in the mid-troposphere, and that it must occur during boreal winter. The former ensures that the system is indeed a canonical WD, and the latter that the large-scale atmospheric conditions are roughly consistent. These conditions remove only a handful of systems, leaving us with 40 WDs that comprise a broad, representative sample.

2.2.2 WRF model

We use version 4.0 of the Advanced Research Weather Research and Forecasting (WRF) model (Skamarock et al. 2008) to simulate our 40 case studies. For each case study, the model is initialised 120 h before the tracked western disturbance (i.e. the centre of mean vorticity in the 450–300 hPa layer) crosses $70^\circ \,\hbox {E}$ for the final time. This ensures a consistent framework as well as a suitable spinup period. Two domains are employed (see Fig. 1): an outer domain with 20 km grid spacing and parameterised convection, and an inner, nested domain with 4 km grid spacing in which convection is treated explicitly. The inner domain was chosen to contain the region of highest climatological winter precipitation, along the foothills of the northwest Himalaya, as well as to capture as much of the upper Indus catchment area as possible, within computational constraints. It also avoids the dry region of Ladakh to the northeast of the basin. The outer domain was chosen to contain a majority of the moisture trajectories associated with extreme winter precipitation events in the inner domain (see Fig. 13(b) of Hunt et al. 2018b), and to appropriately pad the inner domain from the reanalysis boundary conditions. Together, these criteria place the inner domain east of the centre of the outer domain, since most moisture trajectories originate from the west or southwest; this is additionally useful because climatological winds in this region are typically strong westerlies, thus allowing more time for spinup of smaller-scale ‘granular’ features.

For the choice of parameterisation schemes (see Table 2), we follow Norris et al. (2017) who used WRF to generate a year-long high resolution dataset of precipitation over north India, the Himalaya, and the Tibetan Plateau; and also note the similarity with both the NCAR convection-permitting forecast setup and that used by Patil and Kumar (2016) to assess the capabilities of WRF in simulating western disturbances. Vertically, we use 35 eta levels with a model lid at 50 hPa. Initial and lateral boundary conditions are supplied every 6 h from ERA-Interim reanalysis data.

Table 2 Parameterisations and other physics schemes used in the WRF setup

Fig. 1

Map showing the orography of the region with rivers (blue) and glaciers (grey) from the Natural Earth database marked. The outer (20 km resolution) and inner (4 km resolution) domains are marked in red. The boundary of the Indus drainage basin is marked in black

2.2.3 Climate perturbation

The main purpose of this study is to probe the effects of climate change over the northwest Himalaya; dominated by steep orography, the use of coarse resolution climate models from the CMIP5 experiments would be unsuitable for a complete assessment particularly given their biases in this region (Palazzi et al. 2015). To circumvent this problem, we utilise a technique known in its general form as pseudo global warming (PGW, e.g. Kimura and Kitoh 2007; Prein et al. 2017). These simulations are run in exactly the same way as the control simulations, except for applying a climate perturbation signal to the initial and boundary conditions. This technique is demonstrated schematically in Fig. 2.

Fig. 2

Idealised example of the boundary condition climate perturbation methodology, here shown for a hypothetical eastern boundary. The central panel, which represents the climate perturbation signal, is added to the left panel, representing a direct interpolation from analysis. The sum makes the right panel, which is the perturbed boundary condition used in the experiment. If the experiment were an RCP8.5 re-run of the January 1997 case, the central panel would be computed by subtracting the January 1980–2005 CMIP5 historical MMM field from that of the January 2075–2100 CMIP5 RCP8.5 MMM

To compute the climate perturbation signal, we extract monthly (multi-model) means for 2075–2100 for future climate scenarios and subtract the equivalent for the historical experiment over 1980–2005. This method is used for both RCP8.5 and RCP4.5 experiments, but for the pre-industrial we use the final 25 years of CMIP5 model output. This is done for all variables (at all levels, including the surface) present in the reanalysis-derived boundary files: winds, temperature, geopotential, humidity. Perturbing all variables in the same manner is important to ensure a physically consistent setup. Due to the comparatively small size of our domain, spectral nudging is not employed.

The advantages of this technique are several. Firstly, it allows us to directly interrogate the change in impact of a known event were it to happen in a future climate (or “counterfactually”, had anthropogenic climate change not occurred); secondly, it removes the longer-term impacts of changes to the large-scale circulation (e.g. the projected poleward movement of the jet acting to reduce the frequency of WDs in the future), which have been explored previously, allowing us to focus explicitly on how WDs interact with the changes to local moist thermodynamics. Additionally, there are several benefits of using this technique over simpler, statistical methods, such as downscaling the output of GCMs. Firstly, the poor representation of orography in GCMs results in significant local errors since important features such as updrafts—as well as non-linear mechanisms such as vortex–orography interaction—are not simulated (McAfee et al. 2011). Secondly, GCM biases are mitigated because in PGW perturbations, we use the difference between two periods rather than the absolute values.

3 Results

3.1 Validation

We start with a short validation of our model setup. Figure 3 compares mean rainfall from the GPM-IMERG blended satellite product (see Sect. 2.1.3) against the outer domain of our WRF setup for the 26 cases in Table 1 that occurred since June 2000 (the start of the GPM-IMERG dataset). We must be aware of the caveat that there is some uncertainty in IMERG precipitation over high orography (i.e. snowfall). This is because much of the data ingested over the tropical region by the IMERG algorithm comes from the TRMM satellite; TRMM microwave sensors cannot make accurate retrievals over snow-covered ground, and thus typically do not detect near-surface snowfall¹ (Huffman and Bolvin 2018). Even so, in the winter Himalaya, GPM-IMERG performs about as well as any gauge-based dataset as the lower direct measurement uncertainty in the latter is offset by the relative sparseness of observations (Andermann et al. 2011). The attentive reader will also have noticed an apparent discontinuity in panel (b) on crossing the boundary between the inner and outer domains, this is because we have a direct feedback operating in our model setup: where gridpoints are coincident between the outer and inner domains, values of fields in the former are overwritten by those of the latter at each timestep. The broad structures of mean rainfall are the same in both panel (a) and (b) of Fig. 3: a relatively thin band of high rainfall at the base of the Himalaya, extending through the centre of our inner domain and much of the northern part of the Indus Basin; and some less spatially coherent rainfall attached to the orography of the Karakoram and Hindu Kush ranges further to the west. The model confines rainfall to a thinner, more intense band than exists in the satellite observations (in agreement with Norris et al. 2015), but the area-mean values over this region are comparable in both. Given this, and the results of similar WRF setups in this region (e.g. Bhatt et al. 2014; Norris et al. 2015, 2017; Kedia et al. 2018), we can be confident that the model is performing satisfactorily over our domain.

Fig. 3

Comparison of mean precipitation rate [$\hbox {mm} \,\hbox {day}^{-1}$] between GPM-IMERG and outer domain of the WRF setup used here, computed over the 26 cases that occur over the operational period (2000-). The inner domain is marked in red, and the boundary of the Indus catchment region in black. Smoothed elevation contours at intervals of 1.5 km are marked in magenta

3.2 Precipitation footprint

Now, as the focus of this study is constraining and understanding the effects of climate change on WD-induced precipitation in the western Himalaya, the appropriate starting point is to determine the magnitude (and sign) of those changes. The panels on the left-hand side of Fig. 4 show the mean rainfall and snowfall² over all 40 WD case studies for the inner domain. Therein, rainfall has a mean value of $2.17\,\hbox {mm}\,\hbox {day}^{-1}$ and is non-zero over most of the domain, though somewhat concentrated in the foothills towards the northeast. Snowfall, on the other hand has a mean value of $2.42\,\hbox {mm}\,\hbox {day}^{-1}$, falling in more localised areas among the higher altitude regions towards the northeast. These localised areas tend to be in between the rivers, confirming their higher elevation.

Fig. 4

Average precipitation rates [$\hbox {mm}\,\hbox {day}^{-1}$] in the inner domain separated into rainfall (top row) and snowfall (middle row). Left: the mean over all control experiment case studies; centre: the difference between control and pre-industrial means; right: the difference between RCP8.5 and control means. Note that colour scales differ between columns. Glaciers and major rivers are marked in black and blue respectively. Glacier data are from the GLIMS database (Raup et al. 2007)

The panels in the middle column of Fig. 4 show how these values would change had the pre-industrial climate persisted to the present day. Averaged over the domain shown, the absolute change—that is, control minus pre-industrial—is $0.08\,\hbox {mm}\,\hbox {day}^{-1}$ for both rainfall (where it is an increase of about 3.7%) and snowfall (where it is a decrease of about 3.2%). The pattern of the change is not identical across species, but the similarity in the magnitude of change suggests that much of it might be explained by warming-induced phase change.

The panels on the right hand side of Fig. 4 show how the mean precipitation would change in the period 2075–2100, assuming an RCP8.5 future scenario. Rainfall increases by an average of $0.66\,\hbox {mm}\,\hbox {day}^{-1}$, about 31%; this is concentrated in the foothills towards the northeast, notably at a higher elevation than the peak rainfall in the left-hand panel. Snowfall decreases, with a very similar footprint to that of rainfall, by an average of $0.41\,\hbox {mm}\,\hbox {day}^{-1}$ (about 17%). Thus, unlike the counterfactual pre-industrial experiment, there is a net change in regional total precipitation which amounts to an increase of 6%. We note that snowfall, in both of the perturbation experiments presented here, tends to increase (though only slightly) over the glaciated regions in the far northeast of the domain, in contrast to other parts of the domain where increasing environmental temperature causes a reduction.

The higher-resolution inner domain discussed here covers about a half of the area of climatologically high orographic winter precipitation across the foothills of the West and Central Himalaya, and inspection of analogous results for the outer domain show that this selection is representative. In all, these results suggest that from pre-industrial to the present day, and from the present day to future climate scenarios, there is a reduction in snowfall which is generally compensated by a rise in rainfall. This has implications for water management in the region, for which spring snowmelt is an important hydrological resource.

3.3 Elevation dependency

We saw in Fig. 4 that there is a relationship between surface elevation, mean precipitation, and precipitation response to environmental temperature change. We now seek to quantify that explicitly. Figure 5 shows a vertical cross section of temperature, perpendicular to the orography, passing through part of the inner domain, averaged over all case studies. The mean freezing levels ($0\,^\circ \hbox {C}$ isotherms) are shown for each experiment and have average heights of 2.60 km, 2.71 km, 3.06 km, and 3.37 km for pre-industrial, control, RCP4.5 and RCP8.5 respectively. Over the control experiment inner domain, this amounts to an area of terrain slightly less than $27{,}500\,\hbox {km}^2$ above the freezing level. This is a reduction of 3.2% compared to the pre-industrial period and shrinks by a further 12.0% in the RCP4.5 scenario and 22.8% in the RCP8.5 scenario. These values are comparable to the fractional species changes described above, supporting our hypothesis that conversion of snowfall to rainfall in a warming climate is caused largely by the increasing altitude of the freezing level.

Fig. 5

Mean thermal structure of the lower troposphere, averaged over all case studies along the transect (two gridpoints wide) shown in the inset. Filled contours show temperature for the control experiment as a function of altitude along the cross-section indicated in the inset map. Coloured lines denote the mean height of the freezing level ($0\,^\circ \hbox {C}$ isotherm) for each climate experiment. Orography is shown in grey

Fig. 6

Contoured frequency by altitude diagrams (CFADs) for accumulated rain and snow. For each pixel in each case study, the accumulated surface precipitation [separated into rain (left) and snow (right)] is binned by its magnitude and altitude. These histograms are then normalised by altitude. For the control run case studies, the histogram is shown by the coloured grid, with the $3.5\times 10^{-2}$ contour indicated in black. For the other experiments, the full histogram is not shown, but the $3.5\times 10^{-2}$ contour is given in blue (pre-industrial), yellow (RCP4.5) and red (RCP8.5). There are twelve bins per precipitation decade

To develop this further, we look at the elevation-specific responses to warming by using contoured frequency by altitude diagrams (CFADs; e.g. Yuter and Houze 1995; Houze et al. 2007; Liu et al. 2010; Hence and Houze 2011). These are two-dimensional histograms in which we bin the accumulated precipitation for each pixel in each case study by surface elevation and quantity, then normalise the distribution at each elevation bin. Fig. 6 shows these CFADs for the control experiment for both rainfall and snowfall, as well as a representative contour for each of the other experiments.

Rainfall in the control simulation reaches its highest accumulations between elevations of 1 and 2 km, i.e. where the foothills start to trigger orographic precipitation. There is a steep decline in mean accumulation above the freezing level ($\sim 2.5\,\hbox {km}$) where rainfall ceases to comprise a significant fraction of precipitation. Snowfall behaves in a similar but opposite fashion, albeit with a more clearly defined peak; consistent large accumulations persist at elevations above about 2 km, with a sharp reduction below this. A satellite-based study by Puranik and Karekar (2009) showed that the perhaps surprising existence of some snow at very low altitudes is permitted by the cold air intrusions sometimes associated with WDs.

For each experiment in both panels of Fig. 6 we have further drawn a single representative contour at $3.5\times 10^{-2}$, chosen arbitrarily to capture the general structure of the peak. For a given altitude, the contour captures highly likely values of accumulated precipitation during the passage of a WD; where there is no contour, one can conclude that the accumulated precipitation is of small magnitude, high uncertainty, or both. Rainfall, on the left, undergoes two distinct changes as the climate warms: an increase in altitude and an increase in magnitude. Between the pre-industrial and present-day control, there is little change, except for a slight advance in altitude at the expense of retreating snowfall; as future climate scenarios develop, however, the difference becomes pronounced. At altitudes above 1 km, very heavy rainfall becomes significantly more likely, for example: at 1.5 km, rainfall accumulations exceeding 200 mm become 68% and 115% more frequent in the RCP4.5 and RCP8.5 experiments respectively (cf. 16% less frequent in the pre-industrial); at 2 km, they become twice and four times more likely in the RCP4.5 and RCP8.5 scenarios respectively (39% less frequent in the pre-industrial).

In contrast, snowfall does not seem to undergo the frequency shift experienced by rainfall: orography above 3 km appears to be insensitive to changes in environmental temperature. At altitudes below this, frequency is controlled in the expected manner by climate, with increasing environmental temperatures resulting in increasingly rare snowfall at altitudes between 1 and 2.5 km. For example, at 1.5 km, compared to the control experiment, accumulations exceeding 50 mm are 53% and 85% less likely in the RCP4.5 and RCP8.5 experiments respectively (and 46% more likely in the pre-industrial). In summary, the climatological conversion of snowfall to rainfall in this region will result in challenges for water resource management: from a system based on an accumulating reserve and gradual springtime melt, to one of more instantaneous rains that pass through the hydrological cycle more quickly and thus are less readily available for crops later in the year.

3.4 Extreme rainfall

We have seen that rainfall increases over the inner domain as the environment warms, either at the expense of snowfall, or as part of increased overall precipitation. Either mechanism increases the risk of flooding, and the former has a secondary impact in that it also replaces the meltwater used for irrigation in the pre-monsoon months. One of the significant uncertainties arising from previous assessments of WDs and their impacts in CMIP5 model output was the nature of extreme rainfall, which is not well represented at typical GCM resolutions. With our high-resolution experiments, we are in a position to explore how such events, when related to the passage of western disturbances, respond to climate change. Figure 7 shows how extreme values of selected integration periods vary across the four experiments. Each line therein represents a distribution containing all raining gridpoints across all cases for a given experiment.

Fig. 7

Percentiles for rainfall in the inner domain integrated over selected integration periods (2, 24, and 168 h), for each of the four experiments. The distributions comprise values for each individual raining gridpoint. Daily accumulations are calculated over consecutive blocks, starting from the initialisation; weekly accumulations are taken as the single, continuous 168-h period with the highest rainfall for a given gridpoint

The 2-h accumulations represent the rainfall from either short storms or the most intense parts of longer storms, and even at this short timescale are relevant for landslides and flash flooding (e.g. Caine 1980; Crosta and Frattini 2003). At all percentiles, raising environmental temperature acts to increase the rainfall threshold; or, conversely, a given extreme 2-h event will become more likely under the future climate scenarios presented here. For example, at the far-right end of the spectrum, a $65\,\hbox {mm}\,(2\,\hbox {hr})^{-1}$ event—the 99.999th percentile threshold in the control run—becomes three times more likely to occur in the RCP4.5 scenario, and ten times more likely under RCP8.5. Furthermore, the increasing spread with percentile of the four experiments indicates the underlying distributions have increasing variance as well as mean with warming, which has added implications for extreme events (Schär et al. 2004). The tail of one-day accumulation distribution, capturing the effects of the longer-lived storms, behaves similarly to the 2-h case, with extreme events increasing markedly in frequency with increasingly severe future climate projections.

For the 1-week accumulations, which represent the effect of multiple storms embedded within a single synoptic-scale event, the effect is even more pronounced: 99.999th percentile events in the control ($\sim \,530\,\hbox {mm}\,\hbox {week}^{-1}$) become 15 and 50 times more likely in the RCP4.5 and RCP8.5 scenarios respectively. We note that given the size (and quantity of gridpoints) of the domain used, this percentile corresponds to an occurrence somewhere in the western Himalaya once every eight active WD events. Over this accumulation period, we also see a divergence with increasing percentile between the control and pre-industrial extremes; the 99.999th percentile pre-industrial threshold ($\sim \, 450\,\hbox {mm}\,\hbox {week}^{-1}$) has become about six times more likely in the present day. Unfortunately, the observational network in this region is too sparse to support this claim. To contextualise this further, we note that the RCP8.5 1-day accumulation distribution almost overtakes the pre-industrial 1-week distribution at the most extreme values.

3.5 Moisture budget

3.5.1 Current climate

As we saw in Sect. 3.1, only about half of the projected increase in rainfall from the control climate to RCP8.5 can be accounted for by the corresponding decrease in local snowfall (and at any rate, this is an upper bound). The rest must be explained by external sources, and the most practical place to start is an exploration of the moisture budget. We start by taking a single case study—for this we choose the January 1997 WD for its unexceptional nature (it was not associated with any extreme hazards) and citability (see Table 1)—and rerunning the control experiment to extract instantaneous output at 10-min intervals, as opposed to the 2-h frequency we have used in the rest of this study. In Fig. 8, we take a rectangular box within the inner domain, aligned with the orographic gradient and encapsulating the region of highest precipitation, which is shown in the inset.

Fig. 8

Components of the moisture flux budget over the foothills of the western Himalaya during the western disturbance of late January 1997, computed using 10-min output. The top panel shows the budget separated into moisture transport mechanisms, the middle separates net moisture flux into the domain into contributions at each of the domain edges, and the bottom separates net moisture flux into the box into vapour/hydrometeor categories. The box used, and boundary labels for the middle panel, are given in the inset panel. Also shown in the inset panel are the mean 500 hPa winds for this case

In the top panel, we frame moisture³ flux in the context of net sources and sinks of atmospheric moisture in the box; though for ease of comparison, precipitation fields have their sign changed and are presented as positive fluxes (i.e sources) when they are strictly negative (i.e. sinks). To check for budget closure, we compare the difference of explicit sources (net flux into the box and evaporation) and sinks (precipitation)—the so-called residual term—with the directly computed storage term (the integral of column water vapour over the volume of the box), in other words:

$$\begin{aligned} \oint _{\varvec{\ell }}\frac{\partial (\tilde{q{\mathbf {u}}})}{\partial t}\cdot \mathrm {d}\varvec{\ell } = P_{w} + P_{i} + E + \iint _{A}\frac{\partial {\tilde{q}}}{\partial t}\;\mathrm {d}A \, , \end{aligned}$$

(1)

where $\varvec{\ell }$ is the perimeter of the box and A is its area, $P_w$ and $P_i$ are surface rain and snow rates respectively, E is evaporation,

$$\begin{aligned} q = q_\mathrm {vapour} + q_\mathrm {cloud} + q_\mathrm {ice} + q_\mathrm {snow} + q_\mathrm {graupel} \end{aligned}$$

(2)

and the tilde denotes a vertical integral:

$$\begin{aligned} {\tilde{q}} = \frac{1}{g}\int ^{0}_{p_s}q\;\mathrm {d}p \, . \end{aligned}$$

(3)

If the moisture budget is exactly closed, then the storage and residual terms, denoted by the yellow and dotted green lines, will be equal. In the top panel of Fig. 8, we see that the two fields do not match exactly; this is due to the finite size of the box in question, and the fact that it is not aligned with the predominantly southerly winds. This can result in the depositing of moisture in the box by fluxes that are not otherwise accounted for because they enter and leave the box in the same timestep. We note that this effect is negligible and disappears altogether with increasing box dimensions. Now, considering the moisture budget itself, we see that the contribution of evaporation to the total moisture in the box is minuscule: only about 5% of that brought in from outside. It is reasonable to assume, therefore, that evaporation does not play a large part in the changes seen in the previous sections, and that—at least for these winter cases—we should pursue a moisture flux-orientated approach.

The middle panel of Fig. 8 decomposes the net flux (black line) in contributions through each of the walls of the box, labelled according to the most appropriate cardinal point. Southwesterlies prevail in this region during periods of WD-induced precipitation, and the January 1997 case presented here is no different (see inset). Flux through the southern boundary (which is actually rotated $35^\circ$ past due south) is comfortably the major contributor, even allowing for the fact that the eastern boundary—the only other net positive contributor—is five eighths of its length. In this case, the northern and western boundaries represent moisture sinks, opposite in sign and almost equal in magnitude (minus the water left behind) to their parallel counterparts. It is interesting to note that it is the cross-barrier flow, rather than the along-barrier flow that dominates precipitation generation in the foothills.

The bottom panel of Fig. 8 completes the budget by decomposing net flux into its ingredient phases: vapour, cloud (liquid), ice, snow, and graupel. Vapour flux, perhaps unsurprisingly, is easily the most significant contributor. In fact, it exceeds net flux by about 10%, which is accounted for, to within a small margin of error, by the expulsion of snow. This process is best illuminated in the day preceding the WD event (i.e. before midday UTC on Jan 18)—the net flux is close to zero because positive vapour flux is cancelled out by negative snow flux, indicating the presence of an important process whereby some southerly moisture flux is converted into snow on meeting the orography, but advected upslope before accumulating on the surface. The remaining species, ice, liquid cloud, and graupel, constitute only a fraction of a percent towards the total flux between them in this case.

3.5.2 Future climate

So, how does this budget adjust under the conditions of future climate? To address this, we composite selected budget terms for all 40 case studies, centering each about the timestep (2-h period) of heaviest total precipitation averaged over the inner domain. This is done for both the control and RCP8.5 experiments in Fig. 9, and we do not separate the flux by type here, given the result of the bottom panel of Fig. 8. Starting with the control experiment (left panel), we note that the general apportionment of the budget is as in the case study presented above, though some fields now behave a little more tidily. The two precipitation fields peak together and simultaneously with the net flux into the box—plausible given the error term discussed regarding the top panel of Fig. 8. The directional flux components make up similar fractions of the total as in the case study, with southerly (cross barrier) moisture flux dominating. One feature emerging from the composite that was not obvious from the case study is the role of atmospheric storage. Here, a positive value indicates that mean column water vapour in the box is increasing (and vice versa), or in other words, that there is area-integrated moisture flux convergence. This field, averaged over the case studies, is positive and increasing until the peak precipitation event, whereupon it suddenly changes sign, indicating that the occurring precipitation is not only feeding off local moisture flux, but also the basic-state humidity.

Fig. 9

Components of the moisture flux budget computed for the red box in the inset of Fig. 8, averaged over all case studies for a the control experiment and b the RCP8.5 perturbation experiment. The budget for each case study is centred about the timestep of maximum area-averaged precipitation in that box; all other definitions are as in Fig. 8

Figure 9b shows how these fields change when the case studies are rerun subject to an RCP8.5 climate perturbation signal. As we expect from the earlier discussions, total precipitation increases, but solid precipitation decreases, which is approximately accounted for by an increase in the net flux entering the box. In terms of peak fluxes (not time-integrated accumulations), liquid precipitation increases by 41% (whereas solid precipitation decreases by 28%), net moisture entering the box increases by 13%, and the atmospheric storage term increases by 34%. So, the increase in total precipitation is not entirely caused by changes to the local moisture flux, but also the change in the ability of the atmosphere to store moisture during periods of large-scale convergence.

Fig. 10

Net moisture flux into the red box in the inset of Fig. 8, averaged over all case studies. Values against time are shown for two pairs of experiments: a RCP8.5 and control, and b control and pre-industrial. In each case, the mean net fluxes are given in grey and black respectively, with the terms explaining the difference given in blue ($q'$), red (${\mathbf {u}}'$), and magenta ($q'{\mathbf {u}}'$)

However, not only is boundary flux still a significantly larger contributor to changes in the moisture balance than atmospheric storage in terms of absolute values, but also directly supports it. We now turn our attention to understanding how and why these large-scale fluxes change. We can decompose the change in moisture flux into contributions from humidity and dynamics thus:

$$\begin{aligned} (q{\mathbf {u}})' = q{\mathbf {u}}-\overline{(q{\mathbf {u}})} = q'\bar{{\mathbf {u}}}+ {\bar{q}}{\mathbf {u}}'+q'{\mathbf {u}}' \, , \end{aligned}$$

(4)

where the symbols all have their usual meanings, ${\bar{x}}$ is the value of x in the control experiment, and $x'=x-{\bar{x}}$ is the difference between the perturbation and control values of x. In Fig. 10, we use these fields to compare the changes in net moisture flux between the control and RCP8.5 experiments (a) and control and pre-industrial experiments (b). Unlike Fig. 9, these composites are constructed relative to the initialisation time, rather than the time of peak precipitation. Though this is perhaps slightly less informative, there is no guarantee that two parallel case studies from different experiments have peak precipitation at the same time, and choosing one over the other would unfairly bias the results and would cloud their interpretation.

Figure 10a compares the net flux into the box used for Fig. 8 between the control experiment cases and the RCP8.5 experiment cases. The mean flux in the RCP8.5 experiment ($2.84\times 10^6\,\hbox {kg}\,\hbox {s}^{-1}$) is not strikingly different from that of the control ($2.59\times 10^6\,\hbox {kg}\,\hbox {s}^{-1}$), in line with earlier analysis that suggested an increase of several percent. This increase of 9% in the net flux is much smaller than the apparent contribution from increased humidity alone ($q'\bar{{\mathbf {u}}}$: 28%), which is due to mitigating effects from the dynamics (${\bar{q}}{\mathbf {u}}'$: $-20\%$). For comparison, we note that Hunt et al. (2019) found that WDs were likely to weaken by about 12% in an end-of-century RCP8.5 climate.

For completeness, we also look at this decomposition for the change in flux between the pre-industrial and control experiments in Fig. 10b. As expected, given the very weak change in total precipitation seen in Fig. 4, there is correspondingly weak ($\sim \,1\%$) and insignificant change in the net flux; the respective contributions from moisture and dynamics are small and of opposite sign. We can thus be sure that the observed changes to precipitation species in this region are not due to changes in large-scale fluxes over South Asia. We finish by noting that in both experiments discussed here, the humidity-only changes to moisture flux are close to those predicted by the Clausius-Clapeyron relationship of 7% $\hbox {K}^{-1}$: the mean lower-tropospheric temperature over our domain increases by about 0.8 K from the pre-industrial to the control experiment, with a corresponding increase of 5% to the humidity contribution. From the control to RCP8.5 scenario experiment, the mean lower-tropospheric temperature increases by 4.5 K, corresponding to an increase of 28% in the humidity contribution.

Fig. 11

Vertically integrated moisture flux for all case studies over the outer domain (with inner domain marked in red). The left panel shows the mean vector field and its magnitude for the control cases. The middle panel shows the changes to those fields in the RCP8.5 cases considering only changes to specific humidity (analogous to the blue line in Fig. 10). The right panel is as the middle panel but for changes to the wind field (red line in Fig. 10). Both the right and middle panels are coloured by the effect the changes in flux vectors have on the magnitude

So, we have seen that in the RCP8.5 climate scenario, that moisture flux into the western Himalaya increases significantly, as a result of increasing moisture winning the tug-of-war over weakening dynamics. But what does this look like on the large scale? In Fig. 11, we vertically integrate the components from Fig. 10a over the whole outer domain, showing the mean vertically integrated moisture flux for the control experiment cases, as well as the first order perturbation terms ($q'\bar{{\mathbf {u}}}$ and ${\bar{q}}{\mathbf {u}}'$) for the RCP8.5 experiment cases. The control mean (left panel) indicates the expected features: cyclonic flow drawing in moisture from over the head of the Arabian Sea and pushing northwards into the Himalayan foothills. The orographic forcing there results in a channel of particularly high moisture flux that impinges on our region of interest (the inner domain is marked in red). Very little of that moisture makes it far into the mountains, causing the heavy precipitation that is the focus of this study.

The contribution to RCP8.5 moisture flux from changes to moisture alone is given in the middle panel of Fig. 11. Unsurprisingly, the additional moisture present in the warmer tropical atmosphere adds considerably to the existing moisture flux. While this effect is greatest over peninsular India, it has a marked northward extent that easily encompasses the inner domain. Finally in the right-hand panel is the contribution to moisture flux changes from dynamics. This pattern opposes the mean, acting to reduce the flux magnitude, and acts most strongly parallel to the barrier in our inner domain—but even there, it is smaller in magnitude than the contribution from added moisture, corroborating Fig. 10a. The anticyclonic structure of this perturbation lends further credence to its source being the projected weakening of western disturbance intensities in future climate.

4 Discussion and conclusions

The water cycle of the western Himalaya is strongly influenced by the passage of western disturbances (WDs), which provide almost all of the precipitation for the region during the winter months. The rainfall feeds the Indus River throughout the winter, and as the snow melts during the pre-monsoon months, it becomes the majority source of water for the river, and local irrigation in the region. Thus, quantifying how the water cycle here responds to climate change is vitally important, yet challenging due to the steep orography and complex flows surrounding it. In particular, previous studies have shown that CMIP5 models cannot agree on the sign of change for extreme precipitation, nor on whether increased atmospheric moisture in a future climate will conquer the projected weakening of the WDs bringing moisture into the region.

In this study, we have set up four experiments to assess the impacts of different climate change scenarios on WDs (shown in Fig. 12) in the western Himalaya. Each consisted of ten-day hindcasts based on forty previously studied WDs, run in a nested WRF suite with inner and outer domain resolutions of 4 and 20 km respectively. Initial and lateral boundary conditions were derived from ERA-Interim reanalysis data. For the pre-industrial, RCP4.5, and RCP8.5 experiments, the boundary conditions were perturbed using the difference between the CMIP5 historical multi-model-mean, and the corresponding CMIP5 experiment MMM. This allowed us to investigate how existing events would behave differently in a range of distinct climate scenarios.

This technique is called pseudo global warming (PGW), and is advantageous over the more commonly-used downscaling method because it allows us to mitigate against single-model biases by using multi-model mean data as an input, and allows us to inspect the weather associated with historical events (as they would behave in different climates) directly, so we do not need long climate-timescale simulations. The major disadvantage, then, is that those long-period, large-scale effects are not considered, so, for example, we do not see the effects of declining WD frequency in this study.

We found that increased temperatures between the pre-industrial and present-day have resulted in a direct conversion of snowfall to rainfall around the freezing level, which is confirmed by longer observational records; though there has been no significant change to the total precipitation. Furthermore, in the RCP8.5 cases, mean rainfall increases by one third, alongside a significant increase to total precipitation and further decrease to snowfall. In both experiments, however, snowfall actually increases over the glaciated areas further upslope.

The snowfall-to-rainfall conversion, or phase-change component of the difference is a direct result of the warmer atmosphere increasing the altitude of the freezing layer; between present-day and RCP8.5 experiments, the height of the mean freezing layer rose 660 m from 2.71 to 3.37 km, losing a significant amount of surface area. As such, we found that snowfall intensity generally does not seem to be sensitive to climate change at sufficiently high altitudes, notwithstanding the aforementioned glaciated areas.

In contrast, rainfall intensity on all timescales (ranging from hours to a week) increases significantly, even more so at the tail of the distribution. Two-hour rainfall maxima with an 8 WD return period were projected to increase tenfold between present-day and RCP8.5 scenarios, 1 week accumulations with an 8 WD return period, fiftyfold.

High frequency output for a single control case (January 1997) showed that moisture flux is predominantly cross-barrier (i.e. southerly, crossing the Himalaya), with a smaller along-barrier (i.e. easterly) component. Almost all flux entering the region is in the vapour phase, and a small fraction exits the region to the north/northeast, having been converted to snow after interacting with the orography. Event-centred composite analysis was used to show that this is true for the average case, and also highlighted the role of local humidity, which increases through large-scale convergence before precipitation occurs. This analysis also showed that both the contributing fluxes and atmospheric storage terms increase in RCP8.5, leading to the increase in net precipitation seen earlier, with increased net flux being the most important term.

We then decomposed the flux change into linear contributions from changes to moisture and dynamics terms, comparing the control experiment with RCP8.5 and pre-industrial regimes. We found that the increase in RCP8.5 was the sum of a large increase from the moisture contribution and a small decrease from dynamics. The pre-industrial, which had no significant change in net flux (corroborating the earlier finding that there was no change in total precipitation) had small, equal, and opposite contributions from moisture and dynamics. These differences highlight the importance of future work exploring the structural changes to WDs in future climate scenarios.

Finally, we explored the spatial structure of these components by considering their vertical integrals over the outer domain. We confirmed that control moisture flux was highest over the inner domain, at the foothills of the western Himalaya, and showed that the moisture and dynamics components added cyclonic (i.e. southerly) and anticyclonic (i.e. northerly) contributions to the overall moisture flux. The latter result is particularly interesting because it shows that the projected weakening of WDs—of which this is a clear signature—is an important factor in mitigating the load placed on the water cycle by additional atmospheric moisture in future clime projections.

Arguably the most important implication of the results presented here is the effect of increased rain (in both the climatology and extremes) on the surface—where it could result in greater flood and landslide risk. Further work, with dedicated land surface models, will be needed to assess these. Another potentially damaging impact is the loss of snow; in the present day, this snowmelt comprises a majority of the Indus water source, and is thus vital for irrigation; research using longer, seasonal-length simulations would be needed to quantify this robustly. Furthermore, while WDs provide a majority of the winter precipitation to this region, they do not provide the totality. Thus, future work investigating how non-WD precipitation (e.g. from orographically forced storms such as cloudbursts) changes under future climate scenarios is needed to make this assessment complete.

Fig. 12

Snapshots of each of the forty case studies: shown here over the outer domain: showing (modelled, control experiment) 350 hPa streamlines for the timestep in which the vortex crosses $70^\circ \,\hbox {E}$ and accumulated precipitation over the ten days used in the hindcast