Abstract
In order to explore the influence of internal water on the seismic response of hydraulic tunnel, the combined mechanical analysis models of multimaterial including surrounding rock, lining structure, and internal water are built. Based on the explicit central difference method, the dynamic finite element analysis methods for rock, lining, and water are discussed, respectively. The dynamic contact force method is used to simulate the rock-lining contact interaction, and the arbitrary Lagrange-Euler (ALE) method is used to simulate the lining-water coupling interaction. Then a numerical simulation analysis method for combined seismic response of rock-lining-water system in hydraulic tunnel is proposed, and the detailed solving steps are given. This method is used to study the seismic stability characteristics of the water diversion tunnel in a hydropower station, and the displacement, stress, and damage failure characteristics of the lining structure under the conditions of no water, static water, and dynamic water are comparatively analyzed. The results show that the hydrostatic pressure restricts the seismic response of the lining, while the hydrodynamic pressure exacerbates its seismic response and leads to damage, separation, and slip failure appearing on the haunch, which can provide a scientific reference for the seismic design of hydraulic tunnel with high water head and large diameter.
1. Introduction
In order to ensure the national energy supply and optimize the water resources allocation, a large number of long-distance and cross-basin water diversion projects have been built, under construction and planned in Southwest China, which is rich in water resources. Most of the water transmission lines of these projects have to pass through high mountains, and they mainly rely on hydraulic tunnel for water diversion. Hydraulic tunnel generally inevitably passes through the high seismic intensity area in the construction process due to its long extension, so its seismic safety needs detailed analysis and demonstration [1–4]. Different from the ordinary underground tunnel, the lining structure of hydraulic tunnel not only interacts with the surrounding rock, but also interacts with the internal water. In particular, the interactions between surrounding rock, lining structure, and internal water are more complex under seismic load. Especially for the hydraulic tunnel with high water head and large diameter, the hydrodynamic pressure formed by the violently shaking internal water has a great influence on the lining stress. Therefore, the research on the dynamic response characteristics of hydraulic tunnel structure under the combined actions of seismic load and hydrodynamic pressure is conducive to ensure the long-term safe and stable operation of hydraulic tunnel, which has great theoretical value and social economic benefit.
There are many research results on the seismic response of ordinary tunnel [5–8], but the seismic response analysis of hydraulic tunnel is still relatively less. Although the hydraulic tunnel characterized by the internal water is different from the ordinary tunnel, the seismic analysis methods of ordinary tunnels are still suitable for hydraulic tunnels, including analytical method, experimental method, and numerical method. Considering the relatively complex conditions such as the nonlinearity of material parameters, the irregularity of tunnel geometry, and the dynamic interaction of fluid-structure-rock system in hydraulic tunnel, the numerical method is still an efficient and economical way to simulate the seismic response of hydraulic tunnel, which has broad application prospects. In particular, the fluid-structure interaction should be taken into account in the seismic analysis of hydraulic tunnel. Based on the dynamic finite element analysis methods for single-phase fluid and solid mediums, Deng et al. [9] proposed a dynamic explicit finite element integration format for solving the fluid-lining coupling interaction in diversion tunnel. However, this study ignores the rock-lining interaction, which cannot accurately reflect the mechanical characteristics of the lining. Chen et al. [10] built a full three-dimensional (3D) finite element model of a double-line shield tunnel for water supply in Shanghai and studied the seismic performance of the tunnel by defining the nonlinear contact interface of soil-structure interaction. Wang et al. [11] conducted a full 3D seismic analysis of the portal structure of a hydraulic tunnel and investigated the effects of hydrodynamic pressure and rock-structure interaction on the seismic response and damage characteristics of the tunnel portal. Sun et al. [12] built a numerical model considering the fluid-structure-rock interaction in the commercial software ABAQUS and evaluated the nonlinear response effects of hydraulic arched tunnel induced by the mainshock-aftershock excitation. The above studies mainly adopt the simplified added mass method [13] to consider the influence of fluid, which do not reflect the dynamic coupling characteristics of fluid-lining system. Therefore, in order to scientifically evaluate the dynamic stability characteristics of hydraulic tunnel under seismic load, it is necessary to propose an effective method for simulating the combined dynamic interaction of rock-lining-water system.
On the other hand, the violent shaking of fluid under seismic load belongs to a large deformation problem of mesh. The methods adopted to describe the large deformation of mesh mainly include Lagrange method and Euler method in the finite element calculation of continuous medium. The pure Lagrange method and the pure Euler method have their own advantages, but also have serious defects. For the Lagrange method, the solution of motion equation is relatively simple, but the large deformation may lead to mesh deformity. For the Euler method, the calculation accuracy will not decrease due to the material distortion, but it needs to introduce a very complex mathematical mapping when dealing with the boundary problem, which may make a large difference between the calculation results and the real solutions. With the continuous development of computer technology, the fluid-structure coupling analysis method based on the ALE description [14, 15] has been applied in many fields. The ALE method combines the advantages of the Lagrange method and the Euler method and overcomes the defects of the two methods. For the ALE method, the mesh does not depend on the fluid motion, so it can not only track the free surface of fluid but also ensure the accuracy and stability of solving the large deformation problem. At present, some scholars have tried to apply the ALE method to the seismic analysis of hydraulic tunnel. Lou et al. [16] built a full 3D fluid-tunnel-soil coupling numerical model of a large-diameter parallel diversion tunnel and used the multimaterial ALE method to calculate the seismic response of the tunnel under the action of fluid. Wang et al. [17] used the ALE method to simulate the movement and sloshing of internal water and proposed a parallel calculation method based on coupling load balance to analyze the dynamic response of a water conveyance tunnel under nonuniform seismic excitation. Although the ALE method has the unique advantages in describing the fluid motion and the fluid-structure coupling problems, the research results of its application in the seismic analysis of large-scale hydraulic tunnel are still few, and there is no widely accepted dynamic coupling analysis theory for water-lining system.
In this paper, the Lagrange method is used to describe the motions of surrounding rock and lining structure, and the ALE method is used to describe the motion of internal water. Based on the mechanical analysis models of multimaterial in hydraulic tunnel, the dynamic explicit finite element analysis methods for rock, lining, and water are discussed, respectively. The dynamic contact force method [18] is used to simulate the dynamic interaction process of rock-lining system, and the ALE method is used to simulate the dynamic coupling process of lining-water system. Then a combined seismic response analysis method for rock-lining-water system is proposed by combining the above two analysis methods, which is very suitable for simulating the seismic response characteristics of hydraulic tunnel.
2. Combined Seismic Response Analysis Method for Rock-Lining-Water in Hydraulic Tunnel
Hydraulic tunnel is a complex underground structure composed of multimaterial including surrounding rock, lining structure, and internal water, so the seismic response of hydraulic tunnel can be regarded as the conditional combined seismic response of the three materials [19]. The mechanical analysis models of multimaterial in hydraulic tunnel under seismic load are shown in Figure 1. Figure 1(a) is the mechanical model of whole hydraulic tunnel. Figures 1(b)–1(d) are the mechanical models of rock, lining, and water, respectively. is the input load of seismic wave. is the interaction force exerted on rock by lining, and is the interaction force exerted on lining by rock. is the interaction force exerted on lining by water, and is the interaction force exerted on water by lining.
2.1. Dynamic Finite Element Analysis Methods for Multimaterial
2.1.1. Dynamic Finite Element Analysis Method for Rock and Lining
For rock or lining, its momentum conservation equation based on the Lagrange description can be expressed aswhere is the density of rock or lining. is the displacement. is the time. is the spatial coordinate. and are the coordinate directions. is the body force of unit mass. is the stress.
After the rock or lining is discretized by finite element in space, the dynamic balanced equation of finite element based on Lagrange can be obtained:where and are the mass and damping matrixes, respectively. and are the velocity and acceleration vectors of nodes, respectively. is the internal force vector of nodes, which is obtained by the internal force integration of elements. is the external force vector of nodes. For rock, , in which is the load vector of seismic wave and is the contact force vector exerted on rock by lining. For lining, , in which is the contact force vector exerted on lining by rock and is the coupling force vector exerted on lining by water.
The explicit center difference method is used here to solve (2), which does not require iteration and has no convergence problem. Therefore, this method is very suitable for dealing with the highly nonlinear problem with complex boundary and contact conditions. When the motion state of rock or lining at time is known, the time-domain integration formats of displacement and velocity at time can be obtained:where is the time step. The initial conditions are given as
2.1.2. Dynamic Finite Element Analysis Method for Water
For internal water, its mass and momentum conservation equations based on the ALE description can be obtained when it is assumed to be incompressible viscous fluid:where is the water density. is the water velocity. is the convective velocity, and , in which is the mesh velocity. is the water stress, which obeys the Stokes constitutive equation:where is the dynamic viscosity of water. is the water pressure. is the Kronecker function.
Similarly, the dynamic discrete equation of finite element based on ALE can be obtained by discretizing the water in space:where is the transfer matrix. is the external force vector of nodes, and , in which is the coupling force vector exerted on water by lining. The explicit center difference method is still used to solve (8).
From the above analysis, it can be seen that, for each material in hydraulic tunnel, the motion state of nodes at time can be obtained according to that at time . The top priority is to solve the external load of each material. is known here, while , , , and are unknown. Therefore, the key to the seismic response analysis of hydraulic tunnel lies in solving the interaction forces between rock, lining, and water at each time.
2.2. Combined Dynamic Analysis Method for Rock-Lining-Water
2.2.1. Dynamic Contact Analysis Method for Rock-Lining
Many researches [20, 21] have shown that the dynamic contact force method is very suitable for simulating the contact behaviors such as bond, separation, and slip between rock and lining, and this method is easy to be combined with the central difference method.
Node 1 and node 2 are assumed to be a contact node pair of rock-lining system. and are recorded as their displacements, and and are recorded as their contact forces at time . Then the rock-lining contact interface should satisfy the force balance condition:
When is known, the predicted displacements and of the contact nodes at time without considering the dynamic contact force can be calculated by using (3). Then, the normal contact gap and tangential contact gap of the contact nodes can be determined:where and are the unit normal and tangential vectors of the contact nodes, respectively. The direction of is from node 2 to node 1.
In order to solve the dynamic contact force , it is necessary to assume first that the contact nodes are in bond contact state at time . After is calculated, the contact state must be checked according to and , and then must be corrected accordingly. The calculation formula of can be expressed in the form of normal component and tangential component :where and are the normal and tangential stiffness of contact interface, respectively.(1)If and , it can be judged that the contact nodes enter into separation state, and then and are corrected as where is the cohesive force of contact interface. is the control area of the contact nodes.(2)If and , it can be judged that the contact nodes enter into slip state, and then is corrected aswhere and are the coefficients of static and dynamic friction for contact interface, respectively.
When , , and are known, the modified displacements , and modified velocities , of the contact nodes considering the dynamic contact force can be calculated by using (3) and (4).
2.2.2. Dynamic Coupling Analysis Method for Lining-Water
Node 3 and node 4 are assumed to be a coupling node pair of lining-water system. and are recorded as their velocities, and and are recorded as their coupling forces at time . Then the lining-water coupling interface should satisfy the motion continuity condition and the force balance condition:where is the Lagrange coordinate. is the ALE coordinate. The calculation formula of can be expressed as [22]where is the interpolation matrix of water shape function. is the unit external normal vector of lining-water coupling interface. is the stress of water element in direct contact with the lining.
The lining-water coupling process based on ALE can be briefly described as follows [23]: the water transfers force to the lining on coupling interface, while the lining transfers velocity to the water on coupling interface. Therefore, the coupling force is reflected as the velocity boundary condition exerted on the water by the lining.
When the motion state of the coupling nodes at time is known, the prediction-correction thought based on ALE proposed by Wei et al. [22] is used for a reference to solve the seismic response problem of lining-water coupling system at time . This thought is revised reasonably in order to match the dynamic contact analysis of rock-lining system, and its detailed calculation steps are as follows:(1)When is known, calculate the predicted displacement and predicted velocity of the lining node on lining-water coupling interface without considering the coupling force by using (3) and (4).(2)According to , calculate the mesh displacement of water, and then preupdate the mesh coordinate of water.(3)According to , calculate the predicted velocity of water by using (8) based on the preupdated water mesh, and then calculate the coupling force by using (15).(4)When and are known, calculate the corrected displacement and corrected velocity of the lining node, the updated water mesh, and the corrected velocity of water by repeating the previous steps.
2.2.3. Seismic Response Analysis Steps of Hydraulic Tunnel
It must be noted that the lining-water coupling force is ignored in the above rock-lining contact analysis, and the rock-lining contact force is ignored in the above lining-water coupling analysis. Therefore, in order to correctly simulate the seismic response characteristics of each material in hydraulic tunnel, it is necessary to combine the rock-lining contact analysis and the lining-water coupling analysis. The combined seismic response analysis steps of rock-lining-water system at time are given as follows:(1)When is known, calculate , , , , , and .(2)According to and , calculate and , and then check the rock-lining contact state and correct and .(3)According to , preupdate the water mesh. According to , calculate the predicted velocity of water and calculate .(4)When , , , and are known, calculate , , , , , and .(5)According to , update the water mesh. According to , calculate the corrected velocity of water.
3. Engineering Case Analysis
3.1. Calculation Model and Parameters
The installed capacity of the expansion project of a hydropower station is 80 MW, and the water diversion and power generation system is arranged according to a water diversion tunnel, a surge shaft, and two units. The water intake is about 2 km away from the upstream left bank of the dam. The total length of the diversion tunnel from the intake to the surge shaft is 7126 m, whose center elevation is reduced from 80.7 m to 32.0 m, and the total length of the diversion tunnel from the surge shaft to the power house is 148.1 m. The water diversion tunnel is located in the eroded middle and low mountain landform with undulating ground along the route, and the geological conditions are relatively poor.
In this paper, a tunnel section within 80 m range downstream of the surge shaft is selected for the seismic response analysis, whose center elevation is 32.0 m. The excavation section is circular, and the tunnel diameter is 10.5 m. The lining adopts C25 concrete structure, whose thickness is 0.75 m. The maximum piezometric head of internal water is 123.0 m. A 3D finite element model of the tunnel section is built, as shown in Figure 2. The model is discretized completely by the 8-node hexahedrons, which consists of 58080 rock elements, 6336 lining elements, and 14256 water elements. The x-axis of the model is vertical to the tunnel axis along the horizontal direction, and the calculation range is from −60.0 m to 60.0 m. The y-axis is consistent with the tunnel axis, and the calculation range is from 0.0 m to 80.0 m. The z-axis is consistent with the water level elevation, and the calculation range is from −28.0 m to 154.7 m. The contact node pairs are arranged between rock and lining, and the coupling node pairs are arranged between lining and water.
(a)
(b)
The lateral pressure coefficients of initial in-situ stress field are , , and . The values of the physical and mechanical parameters of rock, lining, and water are shown in Table 1. Both the coefficients of static and dynamic friction for rock-lining contact interface are 0.6, and the cohesive force of the contact interface is 0.5 MPa.
3.2. Calculation Conditions
The self-developed 3D dynamic finite element program for underground cavern [24] is used for calculation, in which the proposed combined seismic response analysis method for rock-lining-water system is embedded. Both rock and concrete belong to quasibrittle materials, which show obvious nonlinear deformation characteristics under dynamic cyclic load. They will produce certain plastic permanent deformation and show obvious fatigue damage characteristics before destruction. Therefore, the elastoplastic damage constitutive relation based on the Mohr-Coulomb yield criterion is used for both surrounding rock and lining structure. The yield criterion can be expressed aswhere represents the shear yield criterion, and represents the tensile yield criterion. and are the minimum and maximum principal stresses, respectively. is the tensile strength. is the cohesive force. is a parameter related to the internal friction angle. is the damage coefficient, which is expressed aswhere is the damage coefficient in the direction of i-th principal stress, . is the cumulative plastic deviator strain in the direction of i-th principal strain. is the damage constant.
Viscoelastic artificial boundary condition is applied at the bottom of the calculation model, and free field artificial boundary condition is applied at the four vertical boundaries. The top of the model is built to the free surface.
The representative Imperial-Valley wave is truncated for 20 s as the input load, and its peak acceleration is adjusted to 1.575 m/s2 to match the seismic intensity of the project area, as shown in Figure 3. The input wave is incident vertically from the bottom of the calculation model. Both x-direction and z-direction seismic excitation are considered in the calculation. The x-direction input load adopts the Imperial-Valley wave, and the z-direction input load adopts 2/3 that of the x-direction.
In order to ensure the stability of dynamic calculation, the mesh size in the direction of wave propagation should be smaller than 1/12-1/8 of the shortest wavelength [25]. Most energy of the Imperial-Valley wave is concentrated within 20 Hz (the high frequency components over 20 Hz are filtered out), so the maximum mesh size of the calculation model is determined to be 6.75 m, while the maximum mesh size of the finite element model built in this paper is within 6.75 m, which can satisfy the requirements of wave propagation in surrounding rock.
The calculation is divided into three working conditions: (1) no water condition; (2) static water condition (the hydrostatic pressure is considered, while the hydrodynamic pressure is not considered); and (3) dynamic water condition (both the hydrostatic and hydrodynamic pressures are considered). Three monitoring points are arranged on the top arch, haunch, and bottom arch at the middle section of the lining to monitor the hydrodynamic pressure exerted on the lining, the displacement, stress and damage characteristics of the lining, as shown in Figure 4. It must be noted that the tunnel excavation, lining support, and hydrostatic pressure application should be carried out before the seismic load is applied.
3.3. Calculation Results and Analysis
3.3.1. Hydrodynamic Pressure Analysis
When the dynamic coupling interaction between the internal water and the lining is considered under seismic load, the hydrodynamic pressure produced by the internal water exerted on the lining can be calculated. The hydrodynamic pressure time-history of the three monitoring points is shown in Figure 5.
It can be seen from Figure 5 that the variation laws of the hydrodynamic pressure time-history of the three points are basically consistent. At time 0–7 s, the hydrodynamic pressure gradually increases with time and increases to the maximum value due to the violent fluctuation of seismic wave. At time 7–11 s, it gradually decreases due to the decreasement of seismic strength. After time 11 s, it tends to be stable, while the seismic strength is already very small. The maximum hydrodynamic pressures of the top arch, haunch, and bottom arch are 0.256 MPa, 0.466 MPa, and 0.326 MPa, respectively, which all appear at time 5.2 s. On the whole, the hydrodynamic pressure of the haunch is larger than that of the other parts.
3.3.2. Lining Displacement Analysis
The peak displacements of the three monitoring points all appear at time 14.9 s both in x-direction and z-direction. Taking the peak displacements at time 14.9 s as the maximum displacements, the maximum displacements of the monitoring points under the three working conditions in x-direction and z-direction can be obtained, as shown in Figure 6.
(a)
(b)
It can be seen from Figure 6 that in x-direction, the maximum displacements of the top arch, haunch, and bottom arch under condition (1) are 6.00 cm, 6.19 cm, and 5.88 cm, respectively. The maximum displacements of the three parts under condition (2) are 0.33 cm, 0.31 cm, and 0.27 cm smaller than that under condition (1), respectively. While the maximum displacements of the three parts under condition (3) are 0.31 cm, 0.60 cm, and 0.39 cm larger than that under condition (1), respectively. In z-direction, the maximum displacements of the three parts under condition (1) are 3.55 cm, 3.74 cm, and 3.86 cm, respectively. The maximum displacements of the three parts under condition (2) are 0.15 cm, 0.19 cm, and 0.17 cm smaller than that under condition (1), respectively. While the maximum displacements of the three parts under condition (3) are 0.19 cm, 0.24 cm, and 0.38 cm larger than that under condition (1), respectively. These show that the lining deformation is restricted by the internal water when only the hydrostatic pressure is considered, while the hydrodynamic pressure exacerbates the displacement response of the lining. In addition, the transverse deformation of the haunch and the vertical deformation of the bottom arch are the mostly affected by the hydrodynamic pressure.
3.3.3. Lining Stress Analysis
In order to explore the influence of hydrodynamic pressure on the mechanical characteristics of the lining, the maximum principal stress is taken as an example to illustrate. The maximum principal stress time-history of the three monitoring points under the three working conditions are shown in Figure 7.
(a)
(b)
(c)
It can be seen from Figure 7 that under condition (1), the maximum principal stress of the lining is almost 0 MPa at time 0 s; at time 0∼7 s, it increases rapidly; after time 7 s, it tends to be stable. Under condition (2), the maximum principal stress of the lining is about 0.76 MPa at time 0 s; at time 0∼7 s, it decreases gradually; after time 7 s, it tends to be stable gradually. Under condition (3), the maximum principal stress of the lining is about 0.76 MPa at time 0 s; at time 0∼7 s, it fluctuates violently and increases rapidly to the maximum value; after time 7 s, it tends to be stable gradually. On the whole, the maximum principal stress of the haunch is larger than that of the other parts.
The maximum principal stress of the lining is mainly reflected as tensile stress. Under condition (1), the maximum tensile stresses of the top arch, haunch, and bottom arch are 0.65 MPa, 1.09 MPa, and 0.92 MPa, respectively. Under condition (2), the maximum tensile stresses of the three parts are 0.80 MPa, 0.88 MPa, and 0.79 MPa, respectively, most of which decrease due to the restriction of the hydrostatic pressure compared with condition (1). Under condition (3), the maximum tensile stresses of the three parts all reach the concrete tensile strength (1.27 MPa), which increase obviously compared with condition (1). These show that the hydrodynamic pressure exacerbates the stress response of the lining and may lead to the lining damage.
3.3.4. Lining Failure Analysis
According to the above analysis, when there is no water or only the hydrostatic pressure is considered, no damage appears on the lining under seismic load. While when the hydrodynamic pressure is considered, the tensile damage appears on the lining. Based on the calculation formulas of damage coefficient (equations (17) and (18)), the damage coefficient time-history of the three monitoring points can be obtained, as shown in Figure 8.
It can be seen from Figure 8 that at time 0∼3 s, the damage coefficient of each part is 0, which shows that the lining is in elastic stress state. At time 3∼7 s, it gradually increases with time due to the violent fluctuation of seismic wave. After time 7 s, it basically remains unchanged. Due to the plastic strain accumulation of concrete, the damage developing process of the lining shows irreversibility. On the whole, the damage of the haunch is more serious than that of the other parts.
The damage coefficient distribution of the lining (taking the middle section of the lining) under the hydrodynamic pressure after the earthquake is shown in Figure 9.
It can be seen from Figure 9 that all parts of the lining are damaged, and the damage degree of each part is different. From the view point of radial direction, the lining part in direct contact with the internal water is seriously damaged, and the closer the lining is to the rock, the smaller its damage degree is. From the view point of tangential direction, the haunch is the most seriously damaged, and the damage tends to expand to its both sides. These show that the hydrodynamic pressure has a great influence on the mechanical behaviors of the lining and can directly affect its seismic safety.
The rock-lining contact interface is basically in bond contact state before the earthquake. As the lining damage develops gradually under seismic load, the rock-lining contact interface is also gradually destructive. The separation and slip zone distributions of contact interface after the earthquake are shown in Figure 10.
It can be seen from Figure 10 that the separation and slip zones of rock-lining contact interface are mainly distributed in the haunch and its sides, while the top arch and the bottom arch are well bonded. This is mainly because under the combined actions of seismic load and hydrodynamic pressure, the contact interface of the haunch bears more load, and the damage of the haunch is the most serious, which lead to the separation and slip failure appearing on this part. Therefore, the influence of hydrodynamic pressure cannot be ignored on the seismic design of hydraulic tunnel, and the haunch is the seismic weak part of lining structure.
4. Conclusions
(1)Combining the dynamic contact analysis of rock-lining interaction and the dynamic coupling analysis of lining-water interaction, a combined seismic response analysis method for rock-lining-water system in hydraulic tunnel is proposed. Taking the water diversion tunnel project in a hydropower station as an example, the dynamic response characteristics and failure mechanisms of the lining structure under seismic load are simulated and analyzed, and the research results can provide a scientific reference for the seismic design of hydraulic tunnel with high water head and large diameter.(2)The hydrodynamic pressure produced by the internal water increases with the increasement of seismic strength, and the hydrodynamic pressure of the haunch is larger than that of the other parts. Compared with the no water condition, both the maximum displacement and maximum tensile stress of the lining decrease under the static water condition, while they both increase obviously under the dynamic water condition. These show that the hydrostatic pressure restricts the deformation and stress response of the lining, while the hydrodynamic pressure exacerbates the deformation and stress response of the lining.(3)When the hydrodynamic pressure is considered, the tensile damage appears on the lining, and the damage coefficient increases obviously with time during the violent fluctuation of seismic wave. The lining part in direct contact with the internal water is seriously damaged, and the haunch is the mostly seriously damaged. The separation and slip zones of rock-lining contact interface are mainly distributed in the haunch and its sides under the combined actions of seismic load and hydrodynamic pressure. Therefore, the haunch is the seismic weak part of lining structure.Data Availability
The data used to support the findings of this study are included within the article.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This study was supported by the IWHR Research and Development Support Program (Grant no. EB0145B082020). This support is greatly acknowledged and appreciated.