Spatial variability of karst and effect on tunnel lining and water inflow. A probabilistic approach

Spatial variability of karst and effect on tunnel lining and water inflow. A probabilistic approach

Tunnelling and Underground Space Technology 97 (2020) 103248 Contents lists available at ScienceDirect Tunnelling and Underground Space Technology j...

5MB Sizes 1 Downloads 35 Views

Tunnelling and Underground Space Technology 97 (2020) 103248

Contents lists available at ScienceDirect

Tunnelling and Underground Space Technology journal homepage: www.elsevier.com/locate/tust

Spatial variability of karst and effect on tunnel lining and water inflow. A probabilistic approach

T



Karl Yaua, Chrysothemis Paraskevopouloua, , Spyridon Konstantisb a b

School of Earth and Environment, University of Leeds, Leeds, UK Ruler Consult Ltd, London, UK

A B S T R A C T

Probabilistic methods can provide more insight and facilitate rationalized risk allocation when used in the framework of geotechnical baseline reports (GBR) by assessing the influence of a spatially variable distribution of cavities to tunnelling projects. A novel, probabilistic approach of the spatial distribution of karstic cavities is presented to assess both the risk of liner instabilities and water inflow when tunneling or mining through karstified rock masses. The presented research work is based on the following geotechnical baseline statement: “There is a 10% probability of encountering a sediment-filled cavity, with a maximum size of 2 m3”, which was assumed to be correct for the analysis. A Matlab script generates a random distribution of karst cavities, modelled as a circular shape, which is transferred into FEM 2D and 3D numerical analysis. The measured variables are: (a) the total discharge velocity (TDV) and (b) the index for capacity utilization of linings in tunnels (CULT-I), which act as a measure for liner stability assessment and water inflow rates, respectively. The approach can provide added value when used as a probability density function, in assessing the risk of a certain threshold to be reached, such as liner failure or a 50% increase in TDV values from baseline conditions. The model is most useful when analysed in 3D because the entire length of the tunnel alignment can be simulated. However, the approach is limited by the disproportional (to the anticipated results) computational effort required of analysing a 3D mesh and is restricted in applicability to the baseline statement assumption. The methodology presented provides a framework for further investigations that can be applied to the contractual requirements and specifications of a project and to varying baseline statements in order to to assess the risk from karst cavities to tunnelling works.

1. Introduction Tunnel construction faces many technical challenges, most of which stem from the heterogeneity in the ground surrounding the tunnel. The spatial variability and the associated risk of uncertainty needs to be carefully considered and assessed prior to project design and construction, as encountering unforeseen conditions can lead to, sometimes significant and even catastrophic project delays and cost overruns (Paraskevopoulou & Benardos, 2013). An example of such uncertainty comes from the spatial variability of solution cavities, when tunnelling through karstic limestone, which comprises of approximately 10–15% of Earth’s ice-free continental surface (Ford & Williams, 2013). Methods for karst risk evaluation have largely been qualitative or semi-quantitative and are limited to risk classification. Quantitative approaches may assess cavity parameters but not in a probabilistic fashion. Probabilistic approaches could be implemented for tunnels excavated through karstic geology as these can provide insight and facilitate rationalized risk allocation when used in the framework of geotechnical baseline reports (GBRs). GBRs are increasingly being used in tunnelling projects worldwide, with the baseline statements introduced and included in the contract documentation being primarily based on qualitative or semi-quantitative assessments of the



investigation results, past experience and engineering expert judgment. However, the baseline statements are currently seen solely as a risk allocation and demarcation borderline and contractual means for dispute resolution. They may have direct impact on the contract price and risk contingencies but are currently considered only in relation to the particular risk described in the baseline statement without assessing its impact on and interaction with other baseline statements and failure mechanisms that may also have a significant impact on the project’s finances. With particular reference to karstic cavities, a Contractor currently prices accordingly to effectively manage with appropriate mitigation and/or reinstatement measures the impact of encountering the prescribed karsts up to the baseline statement without considering the impact for example on the tunnel lining and any potential groundwater ingress, which will most probably require additional pre and post mitigation measures with direct impact on the contract price. When a closed face Tunnel Boring Machines (TBMs) is used to excavate a tunnel, it is now becoming increasingly common as it is also required by the project owners, to use geophysical exploration methods in an effort to identify in advance and ahead of the tunnel face ground anomalies that may affect the boring sequence, such as karstic cavities. If any identified cavities are not inside the excavation contour, the following questions arise: Should the identified risk be mitigated or

Corresponding author. E-mail address: [email protected] (C. Paraskevopoulou).

https://doi.org/10.1016/j.tust.2019.103248 Received 1 April 2019; Received in revised form 17 December 2019; Accepted 17 December 2019 0886-7798/ © 2019 Elsevier Ltd. All rights reserved.

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

hydraulic pressures can affect the roof and face stability leading to collapse, whilst settlement of fine soil or weathered limestone can occur from consolidation of material as a result of groundwater lowering. Tunnel induced settlements are not usually a major concern from karst tunnels unless tunnelling is taking place in an urban environment, where structures may be sensitive to ground deformations and/or changes to the ground water level (e.g. timber piles). Limestone is generally a strong material that does not consolidate under reduced pore pressure. Lastly, continual flow can erode pathways and wash out fines, leading to material degradation (Marinos, 2001).

not? Are these cavities going to affect the lining stability and/or serviceability and result in water ingress? Hence, apart from the baseline statement which is triggered when a dispute arises, it is necessary to investigate the influence of the element prescribed in a baseline statement on other baseline statements (e.g. ingress of groundwater) or on failure mechanisms (e.g. ULS and SLS, such as gasket opening, of the segmental lining). The present contribution is attempting to provide a methodology to evaluate such interdependencies, currently not taken into account, that may have to be considered by the project stakeholders during the tender, design and construction stages of a tunnelling project. The primary goal of this paper is to assess and evaluate the applicability of probabilistic methods in quantifying the risk of encountering karstic cavities in a tunnel alignment and the potential effect on liner instability and water inflow, since having a more elaborated risk evaluation and rationalized risk demarcation lines can be proven vital for such projects. The proposed framework provides guidance on assessing potentially anticipated conditions when tunnelling through a rock mass with uncertain karst presence and spatial distribution. More specifically, this research work focuses on producing a generalized model to quantify the probability of encountering solution cavities during tunnelling through karst prone areas, using simplified assumptions and the baseline statement that: “There is a 10% probability of encountering a sedimentfilled cavity, with a maximum size of 2 m3”.

2.2. Uncertainty, hazard and risk Uncertainty can be defined as either aleatoric or epistemic. Aleatoric uncertainty arises from temporal or spatial, natural variation, such as discontinuity structure, material proper-ties and voids in karst terrain. Epistemic uncertainty stems from a lack of fundamental understanding and limitations in data collection, such as inappropriate models and laboratory testing procedures. Tools to deal with parameter uncertainty in geotechnical engineering include partial factors, quantification using statistical procedures and the observational approach (van der Pouw Kraan, 2014). Reliability analyses, as Langford and Diederichs (2014) performed, can also be carried out to assess the uncertainty in ground response and what support is required using a probabilistic approach. Acceptability criteria must be de-fined to establish the failure conditions for the project. A hazard can be defined as a situation that has the potential to cause harm to life, property, the environment or finances. The main potential sources for tunnelling hazards are ground and groundwater interaction with the construction, contaminants, existing structures and human errors in construction technique. These can lead to tunnel collapse or distortion, water inflow and external effects to people, the environment, property and finances. Hazards during tunnelling usually result from insufficient information about the ground or a failure to comprehend the geology in engineering terms (Eskesen et al., 2004). Risk can be defined as the probability of a hazard occurring in a given timeframe and the magnitude of the resultant consequences. The main consequences from encountering karsts are delays and cost overruns (Paraskevopoulou and Bernardos, 2012; Paraskevopoulou and Benardos, 2013; Benardos et al., 2013), as any unforeseen conditions can lead to water inrush and potential for collapse, whilst considering that unanticipated cavities may be encountered would need preventative mitigation measures employed, such as grouting of cavities and the surrounding rock, impervious liners, drainage and dewatering. The cost of such techniques can be considerable and typically their application and employment in a contractual framework requires a cost/benefit analysis to be carried out, i.e. balance between accepted risk and reducing costs (Eskesen et al., 2004). The risk of catastrophic water inrush is mitigated by segmental liners around the tunnel, that are typically watertight for water pressures of up to 600 kPa, with acceptable limits for water ingress of up to 2.0–2.5 m3/s (Kolymbas, 2008). However, absolute water tightness is rarely achieved in bored tunnels below the water table, as concrete liners are not impervious and not build perfectly which will allow water under pressure to gradually permeate the lining (ITA, 1991). Drainage to lower the water level around the tunnel and grouting of foreseen cavities are commonly used techniques to reduce risk of water inflow. Drainage is usually more effective and cheaper but has the risk of settlement of overlying soil layers. Grouting can be targeted to selected cavities or ground curtains can be projected ahead of the face and around the tunnel. Grouting lowers permeability and reduces risk of catastrophic inflow. Tunnels strongly interact with the surrounding rock mass, acting as a self-supporting system. The nature of the rock mass sources aleatoric uncertainty and large variation in parameters. Therefore, risk assessment in project management is vital for such projects. Minor delays in

2. Background 2.1. Karst Limestone is a good medium for tunnelling but when karstified it can introduce hazards and the associated risks in the construction process. Karst can be described as a type of process result that occurs from the dis-solution of highly soluble rocks, such as limestone and gypsum by water movement. This process forms extensive water networks made up of cavities and discontinuities that induce secondary porosity (fractures) and tertiary porosity (conduits). The latter contributes to spatial uncertainty as the formation of pathways is both a function of rock dissolution rate and the geological structure. In addition, karst systems continually evolve over time because of the feedback loop between the dissolution and material property. Continual flow can further erode pathways and wash out fines, leading to material degradation (Marinos, 2001). Most dissolution happens at or near the bedrock surface, resulting in surface karst features. However, underground karst can be present despite the absence of any surface indications and karst should always be assumed to be present in a carbonate terrain unless proved otherwise (Ford and Williams, 2013). Ground investigations in potentially karstified limestome for karst identification present some of the largest challenges in fore-seeing the complex characteristics of the system. Geophysical surveys can be used to detect voids and caverns in order to inform further investigation but can be limited by survey depth, resolution and extent. Hydraulic and tracer testing can lack proper assessment of the groundwater network. Also, probes usually need to be closely spaced for a sufficient model of the subsurface cavities. For a 90% probability of encountering one cavity of 2.5 m diameter, there would need to be 2500 probes per 10,000 m3 (Waltham and Fookes, 2003). Karst cavities can be infilled by air, water, sediment (Fig. 1) or a mixture, which can lead to potential issues, such as water inflow, mud flow and tunnelling through weak fill material or void. The main problems when tunnelling through karst arise from groundwater and any voids that intersect the tunnel alignment. Groundwater may seep into underground works, posing threats of contamination and disturbance of surrounding structures (Song et al., 2012). Even minor leakage into the tunnel can be proven to be problematic, especially if groundwater is contaminated. Large changes in 2

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

Fig. 1. Example of karstic void partially filled with clay and silt in Dodonis tunnel, Greece (modified after Marinos, 2001).

construction and maintenance of the works can have major financial implications. The risk in tunnelling is always significant due to high costs, limited design flexibility and ensuring safety during construction. The biggest issue surrounding risk evaluation of tunnelling through karst is simulating the conditions of conduits and cavities. This aleatoric uncertainty is difficult to replicate but can be quantified using simplifying models and assumptions. This will give a generalised estimate of the risk that can be applied across projects but lack accuracy. More accurate probabilistic models can be made by adapting site specific inputs but will be project dependent and less applicable to other situations. When modelling more complex karst interaction, difficulties may arise from the complexities of turbulent flow and 3D fracture patterns (Ford and Williams, 2013).

factual and what is interpretative (van der Pouw Kraan, 2014). Probabilistic analyses of geotechnical uncertainty can be implemented into the GBR, following van der Pouw Kraan’s (2014) approach. 2.4. Risk evaluation Qualitative and semi-quantitative methods include risk registers, fault tree analysis, event tree analysis and bowtie diagrams. Risk registers are produced for most geotechnical engineering projects, which identify the likelihood and consequences of potential risk events and what controls and mitigation measures are to be put in place. The likelihoods are specified over a given time period and each risk is given a qualitative or semi-quantitative risk rating. Quantitative risk assessment in geotechnical engineering usually involves probabilistic analysis to quantify the uncertainty in input parameters. These methods may include Monte Carlo simulations or discrete sampling in the form of numerical modelling (Eskesen et al., 2004).

2.3. Geotechnical baseline reports (GBRs) GBRs are a widely used framework for tunnelling projects that contractually define the anticipated ground conditions (Hatem, 1998). They can be used for risk allocation, bid preparation, technical and commercial evaluation and solving financial and technical disputes during construction (The International Tunnelling Insurance Group (ITIG), 2006). Any anticipated ground conditions stated in the GBR (i.e. in the baseline statements) are the financial responsibility of the contractor, whilst any encountered conditions that exceed the baseline statements shift responsibility onto the owner. This framework can be used during the tendering process to reduce project bid costs, provide a reference baseline for bid evaluation and enable rational financial contingency allocation and during the construction phase as a dispute resolution mechanism. A GBR must be adapted for each project but needs to be objective with no or very little space for interpretation and may typically contain conditions on man-made structures, contamination, geological features, groundwater, etc. Developing accurate or representative numerical baselines can be challenging because of the extensive variability in most geologic formations. This usually requires expertise in interpreting the site investigation, as constraints in project budget and access, limit the extent of data available. Site investigation for karst tunnels will usually include geophysical surveys, groundwater monitoring and borings and the contractors need to be aware of what is

2.5. Previous studies The present contribution proposes an alternative probabilistic approach applied in the framework of GBRs, to quantify the risk of encountering karst cavities which may result in liner instability and excessive water inflow, contrasting most other literature on karst tunnel risk evaluation, which use more qualitative and semi-quantitative methods. Shi et al. (2018) proposed a hazard evaluation model for tunnels in karstified rockmasses by defining the main influencing factors of water inrush, such as aquifer thickness, average temperature, groundwater level, fault and fold geometry. Each index was then given a weightage to produce a comprehensive attribute measure to categorise the risk grade. Another example is the use of a two-tier fuzzy comprehensive evaluation to assess the risk of water inrush in such tunnels (Chu et al., 2017), where, 4 first-tier, influencing factors, including geological structure, hydrogeology, surrounding rock and tunnel characteristics, were established as weighted indices using analytical hierarchy process. Both methods give a good practical evaluation before construction but are prone to subjectivity as engineering judgement is required to 3

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

Fig. 2. FLAC models demonstrating the variation of cave position around the tunnel (after Shan et al., 2017).

encountered. Regardless of the actual geological conditions to be encountered, this baseline statement defines what the contractor bases its bid on and what the owner has to pay on one hand and sets the commercial and contractual demarcation line when disputes may arise during construction due to differing site conditions. This baseline statement defines for contractual purposes that the cavity size range has an upper bound of 2 m3. During construction, cavities of appreciable size within or near tunnel alignments may be dewatered, grouted or filled, so that the rock mass integrity is not compromised. Fig. 3 summarizes the procedure in the form of a flowchart which was used in the analysis and presented probabilistic approach. A 2D model is first created to assess the viability of the modelling procedure and to define the influence zone of the baseline statement and is later adapted to the 3D model, upon which, a more accurate representation of the entire tunnel alignment is assessed. The methodology and approach presented is also applicable to nonmechanized (non TBM) tunnel excavation, i.e. NATM, SCL, SEM methods. In the present contribution however, the full face excavation of a closed face TBM was considered due to the simplicity in modelling the excavation sequence (one single excavation phase, wished in place lining assuming that sufficient face pressure is applied in order to result in relaxation of the advance core in the elastic domain).

determine the weightage factors of indices and have been designed for use on certain case studies. These methods only give a general risk classification and do not allude to any mechanism of water inrush or potential for liner instability or collapse. Numerical modelling approaches have been implemented to karst tunnels but they have assessed different aspects of risk. Lu et al. (n.d.) presented a parametric, stability study of a generalised limestone cave using FLAC. The depth below ground level, exposed span of cave and rock mass rating were varied to determine the failure load of the roof, showing that the mechanism and rates of cave roof failure were dependent on rock structure. Shan et al. (2017) determined the safe thickness between the tunnel and concealed cave, which was simulated on FLAC 3D and verified by monitoring data from the Yuanliangshan tunnel. Locations without a safe thickness caused high water filled caverns to burst, due to squeezing damage to the rock mass by a narrow space of water pressure that results in crack formation. The analysis considered 3 positions for the cave: above, below and the sides of the tunnel (Fig. 2). Numerical models can also be used to simulate the actual movement characteristics of water in rush into a tunnel constructed in karstic limestone. Boahua (2015) presented such a model on FLOW-3D software based on the Baofeng tunnel, demonstrating the key driving forces from gravity and external hydrostatic pressure that must not exceed the internal pressure. Addenbrooke and Potts (2001) considered the influence of tunnel position, spacing, rest period and sequence of excavation on interaction between two tunnels, predicting ground movement and tunnel distortion. Similar interaction analysis has been conducted by Zhang & Anthony (2015). These examples demonstrate the use of parametric studies to evaluate the expected behaviour of other void structures interacting with the tunnel. The approach and framework proposed in this contribution attempts to build on such approaches by assessing cavity parameters in a probabilistic fashion. Further research may include probabilistic models on the range of hazards from karst tunnels, including water inrush of conduits, mudflows, cavity collapse and settlement of surrounding ground.

3.1. Numerical modelling A base model was created, that aimed to define the initial conditions, test the viability of such an analysis and determine the sensitivity of the ground response to the characteristics of a sediment-filled void. As the potentially karstified rock mass was assumed to be homogenous, isotropic and a continuous medium, FEM was selected for this analysis. Elements in the interaction area between tunnel and voids were meshed with smaller elements (more discretised mesh) in order to more accurately represent the in-situ conditions. The movement of underground water through karst cavities and network can be conceptualised by either understanding flow through individual conduits (discontinuum) or by modelling the flow in the fracture network as an equivalent porous medium (continuum). For a discontinuum model to be accurate, information, such as fracture density, orientation, roughness and interconnectivity needs to be known. Modelling water filled conduits presents certain difficulties as the fractures can have a dendritic nature that produces 3D spatial variability. Additionally, the complexity in fluid dynamics between the fracture patterns and turbulent flow poses great challenges in conceptualisation. For the continuum model, laminar flow is considered across the entire cross-section, which fails to capture the restricted flow through voids but is considerably less data and time intensive (Domenico and Schwartz, 1998). For simplicity, the flow simulated in this contribution follows the principles of a continuum model, that uses

3. Methodology The primary aim of the proposed methodology is: (a) to create an approach to probabilistically assess the spatial distribution of karst cavities and (b) to assess what implications this may have when tunnelling (full-face excavation) through such geological medium based on a typical geotechnical baseline statement that can be set prior to the tendering stage of the tunnelling project, such as “There is a 10% probability of encountering a sediment-filled cavity, with a maximum size of 2 m3”. This baseline statement means that in maximum 10 m of a typical 100 m section of the tunnel alignment, karst cavities with a maximum size of 2 m3 filled with sediments may or may not be 4

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

Fig. 3. Flowchart of the methodology used in the analysis presented herein.

typical Limestone formation located in the Arabic peninsula and in particular the so called Simsima limestone located in Qatar the initial model assumed the rock mass was elastic and isotropic. Rock mass parameters are based on Mohr-Coulomb criterion evaluated from the following types of limestone: Slightly Weathered Limestone (SWL) is assumed for the rock mass and Highly Weathered Limestone (HWL) parameters are used on the sediment fill of karst voids. HWL represents the common derivation of infill for karst cavities, which is usually a weakened, more permeable version of the source rock or soil above. Table 1 details some of the parameters selected for this analysis. The liners used are unreinforced concrete segmental linings and modelled as elastic, with default proper-ties of concrete (Table 2). As the analysis investigates also the effects of water inflow into the tunnel from karst cavities, the relative permeabilities between the sediment-filled cavity and rock mass were important for comprehension. Permeabilities of 1 * 10−4 m/s and 1 * 10−3 m/s are used for the limestone mass and cavity infill respectively. These values were selected so that the relative permeabilities could be easily related and had a noticeable impact on the results, whilst being realistic. These values seemed appropriate, when compared to the variation of permeability values with decreasing rock-quality designation (RQD), with in-situ permeability values for limestone generally ranging from 1 * 10−3 m/s to 1 * 10−9 m/s (Qureshi et al., 2014). The assigned permeabilities are assumed to not vary with depth. As with most geotechnical analyses, a range for input values is usually applied so that any geotechnical uncertainty in data acquisition or aleatoric variation in materials can be accounted for. Only individual values for each parameter are implemented because a generalised rock mass model was being considered and this allows focus of the study to be directed on the effect of the random spatial variation of cavities.

an equivalent permeability across the limestone body that is representative of the flow through the fracture network. Karstic caves typically have smooth boundaries defined by the natural processes of dissolution, meaning that the vertices of rectangular cavities would be uncommon and may emphasise instability conditions at the corners. The cavities here were modelled as circular voids, as this is a relatively accurate representation of the structure that is easy to generate. This simplification only requires a centre and radius for input, reducing the number of parameters needed. Other alternative geometry simplifications may include rectangles or ellipsoids, these would bring further uncertainty from determining the length to height ratio and orientations used. In reality, the shape of cavities is irregular as the dissolution process is arbitrary, 3D and interconnecting, so the simplified shape adopted fails to accurately represent that. The interconnection of voids has been considered, as the randomly generated circles can coalesce to form extended globular shapes. Circular voids present the most stable geometry, as stress concentrations at the corners are distributed. This may produce conservative stress estimates around voids and around any voids near the tunnel. Although, the effects on the voids are not directly considered in this study, the effect from the voids to the tunnel will affect the results. The model has initially been simplified to 2D, plane strain conditions in RS2 for the tunnel cross section and voids. Karst cavities, modelled as circular in 2D and idealised as spherical in 3D, would assume a cylindrical shape trending into the z-axis in 2D. This means that only the largest cross section of the void would be considered in the 2D plane.

3.1.1. Material properties The geotechnical rock mass parameters used are derived from a 5

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

Table 1 A summary of the geotechnical parameters of HWL and SWL, *the UCS of HWL is undetermined from the samples used as they were too weak, (simplified from Karagkounis et al., 2016). Rock Unit

Highly Weathered Limestone

Slightly Weathered Limestone

Unit Weight (kN/m3) UCS (MPa) Young’s Modulus (GPa) Poisson’s Ratio Permeability (m/s) Cohesion (kPa) Friction Angle (°) RMR Range

19 – 0.15 0.33 1.6 * 10−6–1.7 * 10−4 20 21 10–30

22 6–20 1 0.33 3.1&10−6–3.0 * 10−4 200 26 30–50

introduced by (Kostantis and Spyridis, 2013). CULT-I provides a measure for comparison of varying conditions and different tunnel configurations when considering capacity limit curves of tunnel liners, whereas the TDV measures the rate of inflow of water into the excavation and could also be used as an indicator of the water pressures that may be expected on the tunnel liner. Fig. 5 shows the principle of the CULT-I derived from a capacity limit curve, indicating the envelope of values of axial force and bending moment inside which a liner can be considered as structurally ‘sound/stable’. Values within the curve indicate a ‘stable’ lining (i.e. there is reserve in the material utilization), whilst failure occurs outside the envelope (i.e. the material has utilized its full capacity). The limit curve used to derive CULT-I in this analysis was specified by the Eurocode 2 EN 1992-1 (EN 1992-1-2, 2004). The CULT-I ranges from 0 to 1. A zero value indicates that the axial force/bending moment pair of values lies on or outside the envelope, hence the full material capacity is utilized or exceeded, respectively and the liner is considered to be on the limit of having been failed and be ‘unstable’. Positive values (from zero towards 1) indicate that the axial force/bending moment pair of values lies in-side the envelope, hence there is material capacity to be further utilized and the liner is considered ‘stable’. It is important to note that in a framework of optimization and probabilistic approach, the objective is to have a ‘stable’ liner for all the possible scenarios/combinations, including material and loading variability, with the maximum capacity utilization possible. Using the base model, the initial conditions of TDV and CULT-I are determined. A single karst cavity of 0.78 m radius (the equivalent of 2 m3 in 2D) was spatially varied around the excavation by changing the orientation and position of the cavity relative to the tunnel. The effect of the cavity to TDV and CULT-I values are recorded for each configuration. When the values returned to the conditions seen in the base model (allowing for 1% difference), the location was recorded. The boundary of the zone of influence was produced from connecting these locations. Fig. 6 shows the zone of influence in relation to the tunnel

Table 2 Liner properties used in the 2D base model. Thickness (m)

UCS (MPa)

Young’s Modulus (GPa)

Poisson’s Ratio

Tensile Strength (MPa)

0.3

35

30

0.15

3

3.1.2. Boundary conditions and geometry of numerical model The analysis considered the full-face excavation with a Tunnel Boring Machine (TBM) of a 6 m diameter tunnel, located at 25 m depth. The boundary conditions are set by a 70 × 50 box, with the tunnel located at the center and all edges restrained, except for the top boundary to represent the surface level (Fig. 4). This boundary was determined to be sufficient for analysis, without being computational inefficient. Hydraulic conditions are set using a steady state groundwater analysis. The groundwater level was assumed at ground surface (0 m hydraulic head applied to the top/surface boundary) and a zeropressure boundary around the excavation to simulate free water discharge into the tunnel. Given the permeability of the rock mass, ground water was allowed to enter the tunnel from the face.

3.2. Zone of influence The zone of influence determined in the base model can be defined as the maximum distance from the tunnel excavation contour where an individual 2 m3 cavity does not have any impact on the tunnel. The variables chosen in the analysis to determine the impact on the tunnel are the index for capacity utilization of linings in tunnels (CULT-I) and the total discharge velocity (TDV). The Capacity Utilisation of Lining in Tunnels Index (CULT-I) is a geometrical comparison of two points in a structural interaction diagram (actual point from the analysis and the max utilisation reserve), through the Pythagorean theorem and was

Fig. 4. RS2 Model geometry model of the 6 m tunnel excavation (white) and 1.56 m karst cavity (the equivalent of 2 m3 in 2D). 6

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

Fig. 5. A capacity limit curve, showing the components needed to calculate CULT-I (modified after Spyridis et al., 2016).

3.3. From to 2D to 3D

location for TDV as a smoothed kite shape, with a longer bottom and tapers at the top and sides. The zone for CULT-I is much smaller, with only the top and bottom showing any appreciable effect. Defining the zone of influence allows the analysis to focus within the boundary where there is any appreciable effect from the cavities, reducing computation time. Moreover, defining the zone allows for analysis of which cavities are affecting the alignment.

Matlab is used to generate the random distribution of karst voids, for subsequent transferal to a suitable format compatible with RocScience software for analysis in 2D analysis. The model aims to quantify the karst cavity distribution properties encountered within the zone of influence of the tunnel. An appropriate boundary is selected so that computation time would not be significant but sufficient for the

Fig. 6. Simplified diagram of the tunnel (black), zone of influence for CULT-I (red) and TDV (blue) in relation to the 6 m diameter tunnel location. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) 7

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

Fig. 7. Examples of different karst cavity configurations (black circles) and tunnel (blue circle) in 2D; b. RS2 model showing the implemented cavity configuration (orange circles) from the DXF file; c. a conceptual 100 m 3D tunnel alignment with a close-up of a 10 m section occupied by karst; d. different lengths of the tunnel alignment (cylinder shape) encountering karst cavities (spheres); e. the mesh of the external boundaries (box) with a 3D tunnel and a cluster of karst cavities in the center; f. an angled view of the tunnel (brown), with 11 query lines around the excavation, surrounded by karst cavities (yellow). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

In a typical 100 m tunnel alignment (Fig. 7c), there may be several combinations of encountering karst cavities that would fulfil the 10% chance of encountering karst. The combinations chosen in this analysis were 1 × 10 m, 2 × 5 m, 4 × 2.5 m and 10 × 1 m sections of the tunnel alignment that encounters karst cavities. Fig. 7c conveys a 10 m section occupied with karst, within a 100 m tunnel alignment. The different lengths of 10 m, 5 m and 2.5 m that represent a section of the tunnel alignment that encounters cavities are shown in Fig. 7d. In numerical analysis the model is meshed (Fig. 7d) and then is imported into RS3. The model used the same default liner properties as the 2D and applied to the surface of the cylindrical excavation. The steady state groundwater conditions are the same, with zero-pressure applied to the entire tunnel surface. The boundaries are restrained on all faces, edges and vertices, with the exception of the surface face, representing ground conditions. The same gravity loading as the 2D model is applied. The 3D model is then computed and 11 query lines (Fig. 7e) are evenly spaced 1 m apart, starting at each end of the 10 m section of karst. It should be noted that the same 2D mesh, finer around the excavation and coarser closer to the boundaries, is used in the 3D analysis. This procedure is repeated 10 times using different alignment lengths and with every model having a completely randomized cavity configuration. The results for the differences in TDV and CULT-I values of each alignment length are then plotted as statistical functions.

cavities to be randomly distributed within the zone. The higher bound is the area that is equivalent to the radius of a 2 m3 cavity (r = 0.78). A random area is first selected using a truncated Gaussian distribution, with −/+3 standard deviation set at 0.1 and 1.92. This is later adapted to a linear distribution, with a linear distribution that generates fewer circles with a higher average area than the Gaussian, with more circles and a higher concentration of average sized circles. For the 2D model, cavities are plotted until 10% of the total area is occupied by the circles (Fig. 7a). Although the set baseline statement (i.e. there is a 10% chance of encountering a sediment-filled cavity, with a maximum size of 2 m3) is to be interpreted in the 3D tunnel alignment, the same assumption was adopted in the 2D for the karst cavities distribution. The plots from the Matlab script are then saved as a scalable vector graphic (SVG), which is converted to a drawing exchange format (DXF) and re-scaled in CAD software. This CAD file is then imported to RS2, where analysis of the model is undertaken. The tunnel excavation is then placed on top of the karst distribution, which removes any voids within that area. This is representative of conditions of a tunnel that has been constructed without any ground relaxation (wished-in-place), which would likely be minimal (or its elastic portion would take place) in a competent limestone rock mass and face pressure balance. This procedure is repeated ten times by running the Matlab script again to produce a new cavity configuration, which was transferred to RS2 for analysis (Fig. 7b). It should be noted that the difference in CULT-I values consisting of all the analyses are plotted as statistical functions for analysis. The 3D model has the same basic outline as the 2D analysis, except in 3D the entire tunnel alignment can be simulated and the baseline statement can be realistically represented. In 3D analysis the tunnel now becomes a cylinder and karst cavities can be modelled as spheres.

4. Numerical results The analysis of the 2D results focuses more on reviewing the TDV contour diagrams investigating the mode of water discharge through the karst to the tunnel. CULT-I values cannot be directly seen in the 8

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

Fig. 8. RS2 analysis of TDV values of one arbitrary karst configuration, excavation (white).

a single cavity that multiple cavities will have a larger zone of effect. The difference in TDV and CULT-I values from the base conditions (%), from all 10 m alignment length models are plotted on a histogram. Fig. 10a shows how the frequency of TDV change (%) can be slightly equated to a normal distribution, with a peak close to zero. However, the negative change in values do not extend as far as the positive changes because the TDV values cannot be negative and can only trend towards zero, meaning there is a limit on negative change. The positive values extend with a sparse distribution beyond 30% and have a maximum of 875%. This is condensed on the histogram as a high frequency of values greater than 100%. Fig. 10b shows the frequency distribution of difference in CULT-I values. The positive values (blue) indicate an increase in difference in CULT-I and a ‘more’ stable lining (with reduced lining utilization) whereas the negative values (orange) represent a decrease and a ’less’ stable lining (but higher capacity utilization). It is reiterated here that in a framework of optimization and probabilistic approach, the objective is to have a ‘stable’ liner for all the possible scenarios/combinations, including material and loading variability, with the maximum capacity utilization possible. From Fig. 10b, it can be noted that the maximum positive difference value is less than the negative maximum difference value. The

analysis with contours. The results in this section are examined and presented in terms of cavity size, distance from the tunnel vicinity and related to the frequency of the difference in TDV and CULT-I values (%). 4.1. 2D analysis The TDV values are highest when cavities are in direct contact with the tunnel, as it can be observed in Fig. 8 (i.e. outline of the pink (semicircle) cavity that is in direct contact with the excavation profile. Cavities that are located within the zone of influence of the tunnel seem to have the contours extend through the cavities (from the tunnel wall and outwards), such as those in the bottom right of Fig. 8. The cavity area and the cavity’s distance from the tunnel, TDV and CULT-I influence zones are investigated herein. The cavity area is randomised by the model using a linear distribution, therefore there is an even spread of the cavity size. The distance of the cavities from the tunnel (Fig. 9) and zones of influence are plotted as histograms that show the variable distribution of distances in the model because it is randomised. It is also indicated that the model produced an excess of cavities, and it should be noted that as a majority were not within the zone of influence. It must be considered that the zone is determined by

Fig. 9. Histogram of the distances of cavities from the TDV zone of influence from the 2D model. Cavities within the zone (orange), cavities outside the zone (blue). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) 9

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

Fig. 10. Frequency distribution of the difference in a. TDV and b. CULT-I values derived from the 2D analysis.

lower values, as seen in the cavity in contact with the bottom of the tunnel in Fig. 11. In 3D, the model is more challenging to fully visualise but the cavities are observed to amplify the TDV values if they are close or agglomerated together. In general, the contour of axial force on the liner shows low values on the tunnel walls and higher values at the crown and bottom. The cavities in contact with the tunnel then seemed to locally affect the axial force on the liner. High values are seen on the direction parallel to the tunnel alignment on the cavity edges, whilst low values occur perpendicular to that (Fig. 12a). The bending moment contours show a slightly different pattern that have medium to low values at the sides of the cavity that are perpendicular and parallel to the tunnel alignment. Between the low values on each side, the bending moment is high, which can be seen as the four red areas around the cavity in Fig. 12b. For a lower CULT-I value to occur, the liner needs to have a higher bending moment and lower axial force, although, there is more influence from a change in axial force. Therefore, the cavity sides that are perpendicular in direction (showing a low axial force) and the cavity corners (red contours on Fig. 12b) produce a lower CULT-I value. The 3D cavity parameters show many of the same traits as the 2D model, apart from the fact that more cavities are produced in the 3D model, meaning that distributions are better-defined. This can be seen with the comparison of Fig. 9 and Fig. 13. A more certain increasing linear trend can be seen in the 3D plots; however, the increase may not necessarily be linear and could be exponential. Fig. 14a allows for a better fit of a normal distribution curve than the 2D model (Fig. 10a). However, there is a noticeable dip between the positive and negative values that does not fit the normal trend. The 3D analyses had more cavities to analyse and so had smoother transitions between the bins on the histogram, however the randomisation of the model still creates discrepancies in the distribution. The highest

distribution of negative values shows a high frequency of small negative changes, with a decreasing trend towards the larger negative values. The negative values do not reach −100% change (indicating a CULT-I equal to zero), meaning that liner failure does not occur in any of the 2D results. The distribution does not have a discernible trend that could be inferred. 4.2. 3D analysis The 3D contour plane intersecting the tunnel cross-section shows similar traits to the 2D contour of TDV values (Fig. 11), with high TDV values where cavities are in direct contact with the tunnel. The cavity sides perpendicular to where the high TDV values occur have much

Fig. 11. Front view of RS3 analysis, showing an orthogonal contour plane of TDV values, excavation (white). 10

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

Fig. 12. RS3 results: a. Axial force and b. Bending moment values of the tunnel liner.

range, although, the negative values have a similar range to the other alignments. The median of the 1 m values shifts the interquartile range up and indicates that the majority of the data points are positive, meaning that there are more values with an increase in CULT-I values. The 2D values have a smaller range, with more values that are positive, like the 1 m alignment length.

frequency bins are the low positive values, where the mean and median values lie. In the 3D analysis (for the 10 m tunnel alignment), there are several CULT-I values that reach a difference equal to 100%, meaning that the liner is outside the capacity curve limit and would fail under the combined effect of the axial force and bending moments (Fig. 14b). The majority of data points are negative, with a large peak in frequency in the bin containing −100% to −75% difference in CULT-I values. The values steadily decline from there, with no discernible change between negative and positive values. However, the positive changes indicate an increase in CULT-I values, meaning the area of that liner is more stable as it is closer to the centre of the capacity curve limit. All values greater than 500% were grouped into an upper bin for clarity. The analysis using different alignment lengths of the baseline statement (10 m, 5 m, 2.5 and 1 m) is simplified as box plots, including the 2D model for comparison. These plots shown in Fig. 15 are used to observe the differences in TDV and CULT-I distribution and the accuracy of the model to smaller alignment lengths. For the TDV values, the longer alignment lengths (10 m and 5 m) have a higher median that shifts up the interquartile range (Fig. 15a) and a slightly larger range, excluding extraneous data points. The 10 m, 5 m and 2.5 m lengths have large maximum values, whilst the 1 m alignment is restricted to a small range. The 2D values have values beyond the range of the whiskers but do not extend as far as the 3D models. However, the negative values of the 2D model extend slightly further than in 3D. The box plot for CULT-I values show a similar distribution for 10 m, 5 m and 2.5 m lengths, although, the 10 m values have a larger range (Fig. 15b). The 1 m length has a smaller overall extent and interquartile

5. Discussion and concluding remarks The model created using the Matlab script generates a true random distribution of cavity sizes for both 2D and 3D that is in line with the geotechnical baseline statement set. It is difficult to determine whether a linear distribution would be likely to occur in the karst geology, as the distribution and size of cavities is dependent to on the specific location. Small cavities are more easily formed; however, cavities can connect to form larger structures (Ford and Williams, 2013). The simplification of cavity shape to a circle/sphere is an appropriate assumption if the geology and hydrogeology are likely to create isolated, oblong-shaped cavities. If this is the case, it should be noted that the irregularities and roughness of cavities are disregarded in the simulation. Other karst morphologies may have structurally controlled cavities or cave and conduit systems, which have much more complex hydrogeological systems that this model would fail to capture. However, the allowance of cavities to agglomerate in the model did allow for larger cavities to form and create more interconnected flow. Such an assumption may or may not be representative of in situ conditions as the formation of karsts and karstic cavities configuration is arbitrary. However, this assumption is in line with the GBR framework and the accepted principle

Fig. 13. Histogram of the distances of cavities from the TDV zone of influence from the 23 model. Cavities within the zone (orange), cavities outside the zone (blue). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) 11

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

Fig. 14. Difference in a. TDV and b. CULT-I values in the 3D analysis. Negative changes (orange) indicate a decreased TDV or CULT-I, positive changes (blue) indicate an increased TDV or CULT- I and difference in TDV or CULT-I values of 100% (red) indicate lining failure. The frequency refers to the recorded value of equally spaced nodes across the 3D excavation liner surface for all the analyses. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

computational effort. The reliability of the results is dependent on the accuracy of the model and to the extent that it is representative of the actual conditions, as well as on the number of model runs. The cavities are modelled as spheres, which presents a shape without irregularities that may not accurately portray the effect of irregular void morphologies. The 10 m model presents the highest values of TDV and CULT-I and demonstrates the worst expected conditions. This means that this could be assumed to be the most applicable for projects with a similar baseline statement as the model presents the most conservative estimate for the baseline statement. Any other combination of alignment lengths would produce a lower likelihood of higher values occurring. The use of the results as probability density functions is limited by the amount of analyses completed because the model may produce random discrepancies, which become less pronounced when more model runs are carried out. However, the results could still be used for assessment of the expected risks that could occur. The main limitation of this approach in assessing risk is that it is very time consuming to set up each model and process it. This is because the 3D model can have complications with the meshing of elements, when karst cavities are in contact with the tunnel and the number of elements in 3D can take considerable time to compute. This could be improved by producing cavities within the zone of influence only or by using shorter alignment lengths that require less karst to be generated. Another model constraint is the use of queries to select the data points used for analysis. The queries fail to capture all values

that the baseline statements must contain measurable contractual descriptions of the geotechnical conditions to be anticipated or to be assumed to be anticipated during construction. The applicability of the proposed model to tunnelling projects is limited to the baseline statement set. The methodology however allows the Matlab script to be easily adjusted for varying cavity parameters and probabilities of karst occurring. The histograms produced using the 2D and 3D model for TDV and CULT-I values could be applied as probability density functions and subsequently transformed to cumulative distribution functions to assess the probability and hence the risk of certain conditions occurring around the excavation affecting the lining and the excavation process, considering the baseline statement. For TDV, this could be useful if there is clause or contract requirement that specified the tolerable probability or a certain amount of water discharge allowed to occur into the tunnel. If the project specified an allowable amount of water pressure, then the TDV values can be converted to water pressure values, if the area concerned is specified. For CULT-I, this process could be useful to determine how likely the considered liner is to fail under the assumed geological conditions and hence to define the probability of requiring special segments or special lining arrangements and accessories to be designed and constructed for certain scenarios of karst/tunnel interaction. The greater frequency of negative values in the 3D analysis completely alters the distribution seen in the 2D version, as expected. The 10 m, 5 m and 2.5 m alignment lengths in the 3D analysis yield relatively comparable results, however with significant differences in

12

Tunnelling and Underground Space Technology 97 (2020) 103248

K. Yau, et al.

References Addenbrooke, T., Potts, D., 2001. Twin tunnel interaction: surface and subsurface effects. Int. J. Geomech. 1 (2), 249–271. Benardos, A., Paraskevopoulou C., Diederichs, M., 2013. Assessing and benchmarking the construction cost of tunnels. In: Proceedings of 66th Canadian Geotechnical Conference, GeoMontreal on Geoscience for Sustainability, September 29–October 3, 2013, Montreal, Canada. Boahua, H., 2015. Numerical Simulation Analysis of Karst Tunnel Water Bursting Movement. EJGE, vol. 20, Bund. 25. Chu, H., Xu, G., Yasufuku, N., Yu, Z., Liu, P., Wang, J., 2017. Risk assessment of water inrush in karst tunnels based on two-class fuzzy comprehensive evaluation method. Arab. J. Geosci. 10 (7). Domenico, P., Schwartz, F., 1998. Physical and Chemical Hydrogeology. Wiley, New York, pp. 506. EN 1992-1-2, 2004. Eurocode 2: Design of Concrete Structures – Part 1–2, 1st ed. BSi, Brussels. Eskesen, S.D., Tengborg, P., Kampmann, J., Veicherts, T.H., 2004. Guidelines for tunnelling risk man-agement: International Tunnelling Association, Working Group No. 2.Tunnelling and Underground Space Technology, vol. 19, pp. 217–237. Ford, D., Williams, P., 2013. Karst Hydrogeology and Geomorphology. Wiley, Hoboken, N.J., pp. 1–208. Hatem, D.J., 1998. Geotechnical baselines: professional liability implications. Tunn. Under-ground Space Technol. 13, 143–150. International Tunnelling Association, 1991. Report on the damaging effects of water on tunnels during their working life. Tunn. Underground Space Technol. 6 (1), 11–76. International Tunnelling Insurance Group, 2006. A code of practice for risk management of tunnel works. Karagkounis, N., Latapie, B., Sayers, K., Mulinti, S.D., 2016. Geology and geotechnical evaluation of Doha rock formations. Geotech. Res. 3 (3), 119–136. https://doi.org/ 10.1680/jgere.16.00010. Konstantis, S., Spyridis, P., 2013. Design optimization of tunnel linings based on risk and probabilistic finite element analyses. 11th International Probabilistic Workshop, Brno – Czech Republic. Langford, J.C., Diederichs, M.S., 2014. Support design for excavations in brittle rock using a Global Response Surface Method. Rock Mechanics and Rock Engineering. https:// doi.org/10.1007/s00603-014-0567-z. Lu, Z., Reddish, D., Waltham, T., n.d. A numerical modelling stability assessment for construction over a generalised limestone cave. Marinos, P., 2001. Tunnelling And Mining In Karstic Terrane; An Engineering Challenge. Geotechnical & Environmental Applications of Karst Geology & Hydrology, Beck & Herring (Eds.). Paraskevopoulou, C., Benardos, A., 2012. Construction cost estimation for Greek road tunnels in relation to the geotechnical conditions. In: Proc. Int. Symp. Practices and Trends for Financing and Contracting Tunnels and Underground Works (Tunnel Contracts 2012), March 2012, Athens. Paraskevopoulou, C., Benardos, A., 2013. Assessing the construction cost of Greek transportation tunnel projects. Tunn. Undergr. Space Technol. 38, 497–505. Qureshi, M., Khan, K., Bessaih, N., Al-Mawali, K., Al-Sadrani, K., 2014. An Empirical Relationship between In-situ Permeability and RQD of Discontinuous Sedimentary Rocks. EJGE, vol. 19 [2014], Bund. R. Shan, R., Zhang, X., Lu, M., 2017. Numerical application of safe thickness between a tunnel and surrounding concealed caves. Geotech. Geol. Eng. 36 (1), 95–104. Shi, S., Xie, X., Bu, L., Li, L., Zhou, Z., 2018. Hazard-based evaluation model of water inrush disaster sources in karst tunnels and its engineering application. Environ. Earth Sci. 77 (4). Song, K.I., Cho, G.C., Chang, S.B., 2012. Identification, remediation, and analysis of karst sinkholes in the longest railroad tunnel in South Korea. Eng. Geol. J. 135–136, 92–105. https://doi.org/10.1016/j.enggeo.2012.02.018. Spyridis, P., Konstantis, S., Gakis, A., 2016. Performance indicator of tunnel linings under geotechnical uncertainty. Geomech. Tunn. 9 (2), 158–164. Van der Pouw Kraan, M., 2014. Rockmass Behavioural Uncertainty: Implications For Hard Rock Tunnel Geotechnical Baseline Reports. MSc Thesis. Queen’s University Kingston, Ontario, Canada. Waltham, A., Fookes, P., 2003. Engineering classification of karst ground conditions. Quart. J. Eng. Geol. Hydrogeol. 36 (2), 101–118. Zhang, W.J., Anthony, T.C., 2015. Numerical study of pillar stresses and interaction effects for twin rock caverns. Int. J. Numer. Anal. Meth. Geomech. 39, 193–206.

Fig. 15. Box plot of the difference in a. TDV and b. CULT-I values (%) for different alignment lengths.

around the excavation, resulting in a less accurate representation of all the values that may occur. Notwithstanding the difficulties and the limitations presented, the authors believe that the proposed methodology can yield very useful conclusions and can be followed to determine the risk of certain conditions occurring in the framework of GBRs. A next step could concentrate in automating the process in reducing the time required to set up the model and produce a time effective karst configuration as per the baseline statement assumed. Declaration of Competing Interest The authors declared that there is no conflict of interest.

13