Mapping surface motion parameters and liquefaction susceptibility in Tribhuvan International Airport, Nepal

ABSTRACT Local site effects and liquefaction are major causes of building and lifeline damage in Kathmandu Valley as depicted by several historical earthquakes including recent 2015 Gorkha earthquake (Mw 7.8). Critical facilities like airports play important role in post-earthquake relief, response and recovery efforts. Tribhuvan International Airport (TIA) is the only international airport in Nepal that is located in alluvial soil deposit regarded as Gokarna formation. To assess the site effects and liquefaction occurrence in TIA, this study performs one-dimensional ground response analysis as well as liquefaction susceptibility analysis using shallow borehole logs. Surface motion parameters estimated from one-dimensional ground response analysis are mapped for 0.3g and 0.5g peak ground acceleration inputs. In addition to this, qualitative liquefaction susceptibility map is prepared considering the liquefaction potential index estimated for each borehole log.


Introduction
Nepal Himalaya lies in one of the most active seismic regions in the world due to continuous convergence of Indian plate beneath the Eurasian plate. Current slip rate of 20 mm/year along with the slip deficit of »10 m in Nepal Himalaya (Ader et al. 2012) highlights the possibility of strong to major earthquakes that may cause massive destruction in buildings and lifelines. Historical earthquakes since 1255 caused enormous losses in terms of casualties, injuries, economy, infrastructures and environment. Several strong to major earthquakes have highlighted damage concentration in Kathmandu Valley even in case of far-field events like the 1833, 1934 and 2015 events. Strong to major earthquakes like the 1934 Bihar-Nepal (M w 8.1) and 1988 Udaypur (M w 6.8) and 2015 Gorkha (M w 7.8) depicted that damage would occur particularly in alluvial soft soil deposits recurrently in case of Kathmandu Valley (Rana 1935;Fujiwara et al. 1989;Gautam et al. 2016aGautam et al. , 2016bGautam & Chaulagain 2016;Gautam 2017aGautam , 2017b. In addition to this, major urban centres, critical facilities and residential areas are situated in these formations, so future earthquakes may also affect the same areas severely. Geologically, Kalimati and Gokarna formations comprise the largest fraction of the alluvial valley. Both of these formations were the most affected ones during several earthquakes when compared to other geological formations due to dominantly occurring soft unconsolidated soils.
Notably, 2001 Bhuj earthquake (M w 7.7) and 1985 Mexico City earthquake (M w 8.3) depicted intense damage, respectively, in Ahmedabad (»330 km epicentral distance) and Mexico City (»300 km epicentral distance) than epicentral regions due to local site effects. Soil liquefaction can be detrimental in terms of damage aggravation in building infrastructures and lifelines as observed during 2015 Gorkha earthquake (e.g. Angster et al. 2015;Gautam & Chaulagain 2016;Gautam 2017b). During 1964 Alaska earthquake (M w 9.2), 1964 Niigata earthquake (M w 7.6) and 1989 Loma Prieta earthquake (M w 6.9), soil liquefaction was one of the major causes of destruction. Rana (1935) reported effects of soil liquefaction especially in southern lowlands of Nepal as well as Kathmandu valley that affected both buildings and lifelines during 1934 Bihar-Nepal earthquake. Liquefaction occurrence is also supported by peak ground acceleration (PGA) occurrence in Kathmandu Valley even in the case of far-field earthquakes like the Bihar-Nepal and Gorkha earthquakes. Recorded PGA in Kathmandu Valley was greater than the global threshold PGA estimated by Santucci de Magistris et al. (2013Magistris et al. ( , 2014 during Gorkha earthquake. As TIA is the only centre for rapid response, relief and recovery efforts for both national as well as international agencies to reach up to the victims within and outside Kathmandu Valley, uninterrupted operation of TIA during and after seismic events is crucial. To the best of author's knowledge, there is not any work that addresses site effects and liquefaction susceptibility for TIA. To fulfil the gap of quantification of site effects and liquefaction potential of shallow soil columns in TIA quantification in terms of surface motion parameters and liquefaction potential are required. Thus ground response analysis and liquefaction susceptibility analysis are carried out in this study.

General geology and geotechnical characterization
Kathmandu valley lies in the Lesser Himalaya region in Nepal, an intermountain basin of Pliocene to Quaternary age with varying sediment depth up to 600 m (Yoshida & Igarashi 1984). Gravity measurements in Kathmandu Valley suggested that the maximum thickness of soft soil deposit is »650 m (Moribayashi & Maruo 1980). Other field investigations depicted more than 300-m-thick sequence of mud and sand. Sediment deposit in central valley is estimated to be nearly 550 m by Katel et al. (1996). Valley outskirt is characterized by outcropping bedrock whereas central valley comprises thick fluvio-lacustrine unconsolidated sediments ( Figure 1). Basement rock in Kathmandu Valley is composed of Phulchoki Group and Bhimphedi Group of the Kathmandu Complex (Stocklin & Bhattarai 1977). The southern part of the valley consists of hill terraces formed during late Pliocene to middle Pleistocene age (Yoshida & Igarashi 1984) with sediment exposure. The southern part of the valley is formed from the Tarebhir formation, Lukundol formation and Itaiti formation (Sakai 2001). The central part of the valley consists of Bagmati formation, Kalimati formation and Patan formation. Bagmati formation was active before the lake formation in valley and is considered to be responsible for sediment deposition in the most parts of the valley. The black clayey central portion is called the Kalimati formation with dark grey carbonaceous and diatomaceous beds of the open lacustrine facies. Patan formation is distributed in and around Kathmandu and Patan cities that consists of fine to medium sand and interlayered silt with clay and fine gravels in some areas. The northern and northeastern parts of Kathmandu valley consist of fluvio-deltaic or fluvio-lacustrine sediments that are mostly sandy facies regarded as Gokarna formation and Thimi formation (Yoshida & Igarashi 1984;Sakai 2001). Kathmandu Valley rocks are intersected by a number of fault systems. The Chandragiri fault and the Chovar fault that cut the colluvial slopes and terraces of the late Pleistocene age are considered to be active faults situated towards the southern part of Kathmandu Valley (Sakai 2001).
Geotechnical characteristics of entire Kathmandu Valley suggest dominant occurrence of fluviolacustrine soft sediments; however, variations in small distance are highly discrepant in many areas (Gautam & Chamlagain 2016). Borehole logs considered for this study depict the occurrence of well-graded loose silty sand in soil columns as shown in Figure 2. As shown by Figure 2, the loose silty sand is followed by medium to dense sand in TIA.

Ground response analysis
One-dimensional seismic site response analysis was carried out using an open-source program EERA (Equivalent Linear Earthquake Site Response Analysis) developed by University of Southern California (available at: http://www.ce.memphis.edu/7137/eera.htm) (Bardet et al. 2000). EERA is a user-friendly excel add-in program which provides platforms to input user-defined strong ground motions, soil damping curves, shear wave velocity profiles, depth of ground water table and maximum PGA scaling. For input motion, acceleration time history of 1991 Uttarkashi earthquake (M w 6.9) was used. PGA recorded at Uttarkashi station was 0.32g at 5.86 s as shown in Figure 3. Similarly, for shear wave velocity profiling, Imai and Tonouchi's (1982) correlation was used to correlate the standard penetration resistance with shear wave velocity. Surface motion parameters for 21 borehole logs obtained from Civil Aviation Authority of Nepal ( Figure 4) were calculated for 5% damping using backbone curves suggested by Seed and Sun (1989) and Idriss (1990) for clay, and Seed and Idriss (1990) for sand. The calculated PGA, peak values of spectral acceleration, soil amplification factor and soil fundamental period were mapped in ArcGIS environment.

Liquefaction susceptibility analysis
Safety factor against liquefaction for each layer was calculated using the simplified method of soil liquefaction analysis as suggested by Seed and Idriss (1971). After the calculation of safety factor at each layer (i.e. 1.5-m depth interval), liquefaction potential index was estimated using the formula suggested by Iwasaki et al. (1982) as in expression (1). All of the 21 borehole logs used for ground response analyses were also subjected for liquefaction susceptibility analysis: where F = 1 ¡ F L for F L 1.0 and F = 0 for F L > 1.0, W(z) = 10-0.5z, wherein z is the depth in metres. F L is the factor of safety against liquefaction calculated as follows: where CRR is the cyclic resistance ratio and CSR is the cyclic stress ratio. CSR is calculated for each depth according to the following expression: where g s is the unit weight of soil, a max is the peak horizontal acceleration at the ground surface that is taken as 0.3g for this study, r d is the stress reduction coefficient and its value is limited to less than 1 (Seed & Idriss 1971), z is the depth below ground in metres and g is the acceleration due to gravity. Estimated liquefaction potential indexes were subjected for further analysis to depict liquefaction susceptibility of q particular site using the system of classification developed by Iwasaki et al. (1982) as outlined in Table 1. Liquefaction susceptibility map was prepared for TIA by extrapolation of results obtained for 21 borehole logs in ArcGIS environment.

Results
Surface motion parameters were calculated for desired maximum acceleration of 0.3g and 0.5g separately for 21 borehole logs. Soil amplification factor, PGA, peak spectral acceleration and soil fundamental period were estimated for each borehole log along with liquefaction potential index. Each of the parameters is described below in detail along with ArcGIS mapping.

Soil amplification factor
The maximum soil amplification factor was obtained in the range of 4.5-6.7 for the desired maximum acceleration of 0.3g as shown in Figure 5(a). Higher values of soil amplification factor were obtained   Iwasaki et al. 1982). Liquefaction potential index (I L ) Susceptibility I L = 0 Very low 0 < I L 5 Low 5 < I L 15 High 15 < I L Very high particularly for runways and taxiways of TIA. Similarly, for the maximum desired acceleration of 0.5g, the maximum amplification was obtained in a close range of 4-5.2 ( Figure 5(b)). Both 0.3g and 0.5g maximum PGA inputs depicted higher values of soil amplification factor in runways and taxiways of TIA. As shown in Figure 5, the northern and central parts of TIA are more likely to amplify under strong to major earthquakes, whereas the southern part shows lower amplification than any other areas.

Peak spectral acceleration
The maximum spectral acceleration was obtained between 0.97g and 1.67g for 0.3g desired PGA input m (Figure 6(a)); similarly in the case of 0.5g scaling, variation in maximum values of spectral acceleration was estimated between 1.30g and 2.48g (Figure 6(b)). Higher values of maximum spectral acceleration were obtained particularly for the borehole logs close to the gas storage location and airport buildings. As depicted by Figure 6, the central-southern areas would be exposed to greater shaking than any other areas in the case of 0.3g PGA input, whereas the north-eastern part is found to be exposed to greater shaking in the case of 0.5g PGA input.

Peak ground acceleration (PGA)
PGA was estimated in the range of 0.30g-0.51g for maximum desired acceleration of 0.3g as shown in Figure 7(a). Apart from this, PGA variation in the case of 0.5g input was obtained between 0.39g and 0.72g (Figure 7(b)). Figure 7 depicts that the greatest shaking would occur exactly on runway course for both 0.3g and 0.5g PGA inputs. In addition to this, the northern fringe of airport shows high PGA occurrence for both PGA inputs.

Soil fundamental period
Analysis depicted the variation of soil fundamental between 0.38 and 0.58 s for 0.3g input acceleration (Figure 8(a)). In the case of 0.5g PGA input, the soil fundamental period was obtained between 0.32 and 0.68 s (Figure 8(b)). As shown in Figure 8, both PGA inputs depict higher soil fundamental period in central parts of airport.

Liquefaction potential index
Liquefaction potential index in TIA was estimated between 0 and 19.9 for the analysed 21 borehole logs. Using Iwasaki et al. (1982) criteria, the liquefaction potential indexes were converted into qualitative susceptibility grades as very low, low, high and very high. Out of the 21 borehole logs, 6 (29%) were estimated to be having very low liquefaction susceptibility. In addition to this, nine (43%) of borehole logs depicted low liquefaction susceptibility, whereas four (19%) of borehole logs were found to be having high liquefaction susceptibility and two (9%) were designated to be under very high liquefaction susceptibility. A qualitative liquefaction susceptibility map as shown in Figure 9 was also prepared by extrapolation of the liquefaction potential indexes (I L ) estimated for 21 borehole logs.

Discussion
Ground response analysis highlights the possibility of severe shaking in TIA in the case of strong to major earthquakes including high soil amplification. Areas like taxiways, runways, gas storage and airport buildings were found to be exposed to greater shaking than the other areas. Soil fundamental period obtained in this study highlights the higher seismic demand for three-to seven-storied buildings in TIA. As most of the buildings in TIA fall under this category, additional considerations in   terms of seismic safety are needed for the existing structures too. Ground response analysis presented in this study has used soil backbone curves, ground motion and shear wave velocity-standard penetration resistance correlation developed for some other parts of the world thus in situ tests in terms of shear profiling are backbone curves are required to better estimate ground motion parameters. Recently developed correlation between standard penetration resistance and shear wave velocity (Gautam 2016) can be important to incorporate in future studies. General trend of liquefaction shows that areas like runway, taxiways and gas storage and building sites are not exposed under high susceptibility, thus liquefaction susceptibility in general can be taken as low for TIA. However, this study considers the water table depth measured during dry season, and seasonal variation in water table depth could be important to include in future studies too. On 25 April 2015, strong earthquake of M w 7.8 struck Central Nepal and adjoining areas. The main shock and several strong aftershocks damaged nearly 750,000 houses and caused 8790 casualties and 22,300 injuries (NPC 2015). None of the building was damaged in TIA due to Gorkha seismic sequence. In addition to this, no case of soil liquefaction was reported; this may be due to moderate shaking of »0.2g. Details regarding ground motion and field observations during Gorkha seismic sequence are reported by several researches (e.g. Angster et al. 2015;Parameswaran et al. 2015;Ahmad & Singh 2016;Gautam & Chaulagain 2016;Gautam et al. 2016aGautam et al. , 2016bRai et al. 2016;Gautam 2017aGautam , 2017bRupakhety et al. 2017). Recent studies highlight the possibility of major earthquakes in Nepal Himalaya (e.g. Bollinger et al. 2014;Wesnousky et al. 2017) that may lead to severe shaking in Kathmandu Valley. In this case, the PGA values can reach 0.3g or even 0.5g leading to widespread damage including the TIA due to ground shaking, soil amplification and liquefaction.

Conclusion
Ground response analysis carried out in this study shows that strong shaking in terms of PGA and soil amplification can occur in TIA in case of strong to major earthquakes. Analysis performed for maximum PGA inputs of 0.3g and 0.5g depicts that runways, taxiways, gas storage locations and airport building may be exposed to greater shaking or higher amplification than other areas. Variation in soil fundamental period suggests that the seismic demand for three-to seven-storied buildings is higher than others, thus additional seismic safety considerations are required for design and construction of such buildings. Liquefaction susceptibility for 70% borehole logs considered in this study is estimated to be low in case of strong shaking; however, some runway areas possess high liquefaction susceptibility too. Mapping presented in this study that is based on ground motion parameters and liquefaction susceptibility can provide insights for condition assessment and preparedness initiatives.