Abstract

This study presents an improved constitutive model for concrete under uniaxial cyclic loading which considers the fatigue stiffness degradation, fatigue strength degradation, and fatigue residual strain increment of concrete fatigue damage. According to the constitutive model, the dynamic response and cumulative damage of the tunnel cross structure under various train operation years were analyzed. The results show that the vibration in the middle of the main tunnel is most violent. With the increase of train operation period, the acceleration in the middle of the transverse passage floor, both sides of the wall corner and the vault increase significantly, and the maximum principal stress increases significantly only in both sides of the wall corner. The compressive damage is mainly distributed at both sides of the wall corner, while tensile damage is distributed in both sides of the inner wall corner. The accumulative damage of the cross structure exhibits a two-stage profile. The size and range of accumulative tensile damage of the connecting transverse passage are greater than those of accumulative compressive damage.

1. Introduction

For long and large tunnels, a certain number of connecting transverse passages are normally set up to meet the needs of the operating ventilation, accident evacuation, fire rescue, and other functions. The use of these passages forms the cross tunnel structures. This kind of structure is complex, which lead to uneven distributed forces for the whole structure. The stress concentration most often appears at the intersection due to the train vibration loads [14]. In addition, during the service life of the tunnel, the cross tunnel structure is subjected to the vibration loads caused by the train for a long time. The concrete material deteriorates continuously, and the structure is damaged and cracked all the time, which ultimately results in structural damage and poses serious threat to the operational safety of the structure [510]. Therefore, it is very important to study the dynamic response and fatigue cumulative damage law of the special structure of the cross tunnel under long-term train which caused loads for the safety of long and large tunnels.

For the properties of the tunnel lining structure material, the tunnel lining structure is mainly made of reinforced concrete in China [1115]. The fatigue of concrete structures under cyclic loading causes the damage and cracking of concrete lining. At present, the fatigue performance of concrete structure is mainly studied using indoor fatigue tests. Aas-Jakobsen [16] proposed a general formula for logarithmic lifetime and cyclic stress; Tepfer and Kutti [17] identified the general formula of the basic parameters by the fatigue test; Holmen [18] studied compressive fatigue performance of concrete cylindrical specimens under fatigue load; Huang et al. [19] studied the propagation law of the main crack of reinforced concrete beams strengthened with prestressed CFRP sheet under the fatigue load; Cao et al. [20] studied the fractal characterization in the evolving damage of concrete structures based on physical model experiments and found the surface-crack distribution of the damaged concrete structures.

The fatigue test can accurately describe the fatigue performance of the material, but the test is very time-consuming, and it is difficult to conduct full-scale tests. Researchers started to use numerical methods to study the fatigue damage of the complex structure. Teng and Wang [21] proposed a two-dimensional damage constitutive model of a reinforced concrete structure; Petryna and Krätzig [22] proposed a calculation method for long-term performance evaluation of reinforced concrete structures considering the accumulated damage; Zhang and Shi [23] used the finite element method to study the interface peel stress and its influencing factors on reinforcement and concrete under the fatigue load; J.-S. Zhu and X.-C. Zhu [24] proposed a simplified method for numerical simulation of the fatigue failure process of reinforced concrete bridge structures under operating loads; Wang [25] established a stochastic damage constitutive model based on modified elastomeric Helmholtz free energy under tension and compressive conditions; Wang [26] proposed an equivalent static analysis method for the fatigue cumulative damage process of concrete components.

For the vibration effects caused by train, studies have been conducted on the dynamic response characteristics of tunnels under the train vibration. Gharehdash and Barzegar [27] used a complex elastoplastic 3D dynamic finite difference model by fully considering the joints to study the dynamic response of the shield tunnel buried in soft soil under the vibration loads; Gupta et al. [28] presented the experimental validation of a numerical model for the prediction of subway induced vibrations; Gupta et al. [29] used a coupled periodic finite element-boundary element model to study the vibration response from a Thalys high-speed train in the Groene Hart tunnel; Lin [30] studied the dynamic response of the tunnel under different conditions, such as the preconstruction of the train vibration load.

However, most of the current methods are complex in theory and cannot simulate the fatigue damage behavior of concrete structures under high cyclic loads. Most of the studies only focus on the fatigue analysis of concrete beams, aiming at the dynamic response of the tunnel under the train causing vibration loads. There is a lack of research on dynamic responses of the tunnel structure under the high cyclic loads and lack of the fatigue damage analysis of the cross tunnel structure, formed by the main tunnel and the transverse passage.

The improved uniaxial cyclic loading constitutive model for concrete is proposed based on the latest concrete uniaxial monotone load constitutive model given by “Code for Design of Concrete Structures” (GB50010-2010) [31], together with the concrete fatigue constitutive relation proposed by J. S. Zhu and X. C. Zhu [24]. This model is able to more accurately simulate the mechanical behavior of the commonly used concrete. The formulas for calculating the concrete fatigue stiffness variables, fatigue residual strength variables, and fatigue residual strain variables are included in the cyclic loading constitutive model. Based on the actual situation of the Shiziyang tunnel project of Guangzhou-Shenzhen-Hongkong Railway Passenger Dedicated Line, numerical analysis models were established, and the dynamic response and cumulative damage characteristics of the tunnel cross structures under train vibration load were analyzed.

2. Constitutive Model of Concrete Uniaxial Monotone Loading

The stress-strain curves of concrete under monotonic compression were obtained according to the test data fitting in the “Code for Design of Concrete Structures” (GB50010-2010) and are as follows [31]: wherewhere is concrete strain; is nondestructive elastic modulus of concrete; is compressive stress of concrete; is peak compressive stress; is peak compressive strain corresponding to the peak compressive stress and can be taken as .

When the concrete is monotonically tensile, the stress-strain curve is as follows [32]:wherewhere is concrete tensile stress; is peak tensile stress; is the peak tensile strain corresponding to the peak tensile stress and can be taken as .

3. Constitutive Model of Concrete under Uniaxial Cyclic Loading

The related research [33] shows that the fatigue damage of concrete structures under cyclic loading is mainly demonstrated in three aspects: stiffness decrease, strength degradation, and residual strain increase. Therefore, according to concrete uniaxial constitutive model of the above specification and the concrete fatigue constitutive relation proposed by J. S. Zhu and X. C. Zhu [24], the fatigue constitutive model of concrete under uniaxial compression can be proposed: that is, the stress-strain curve is as follows:where the residual strain of concrete after times fatigue loads, the peak compressive strain , and the modulus of elasticity after times fatigue loads are considered, and the relevant revised parameters for the constitutive model are as follows:

Considering the effect of peak compressive stress after times fatigue loads, the peak compressive strain can be obtained after times fatigue loads:

Similarly, when concrete is subjected to tensile loads, the formula can be proposed as follows:where the residual strain of concrete , peak tensile strain , and elastic modulus are considered, and the parameters are revised as follows:where is the structural concrete peak tensile stress after loading the times fatigue load and is the peak tensile strain after loading times fatigue load.

3.1. Concrete Fatigue Stiffness Related Variable

According to the relevant fatigue test, Holmen [18] proposed the degradation formula for the concrete elastic modulus:where is the concrete fatigue life.

3.2. Concrete Fatigue Residual Strength Variables and

The residual fatigue strength of concrete is related to the number of fatigue load cycles and the maximum and minimum stresses of the load [21].

The study [32] shows that the maximum total strain when concrete is broken under tensile and compressive fatigue loads is equivalent to the strain corresponding to the maximum stress of fatigue load in monotonic loading softening section, as shown in Figure 1, the point B in the stress-strain curve of concrete under uniaxial static load and fatigue process. It is assumed that the concrete fatigue residual strength envelope [34] is represented by the shape of softening section of monotonic loading stress-strain curve of the concrete. Therefore, the concrete residual strength envelope can be obtained by the softening section shape of monotonic loading stress-strain curve of the concrete.

According to the softening section shape of the monotonic stress-strain relationship curves of the concrete, the envelope equation of the residual fatigue strength of the concrete can be obtained as follows: [24]where is the function of fatigue load times : that is [24],The tensile residual strength envelope of concrete is shown by the following equation [24]:

Considering the initial conditions of concrete residual strength and failure criterion, the boundary conditions were taken into account in the concrete compressive and tensile fatigue residual envelope equation. The concrete compressive and tensile dependent variables A and B were available. Substituting those variables to formulas (11) and (13), the fatigue residual strength of concrete can be obtained.

3.3. Residual Strain of Concrete Fatigue

Holmen obtained the formula of the fatigue residual strain of concrete by the curve fitting the experimental data, without considering the stress ratio. The formula is as follows: [19]:where is the residual strain of the structure after 1 time fatigue load; , are the maximum stress and minimum stress of the structure when the fatigue load reaches the upper or lower limit of stress, respectively. is the longitudinal total strain when the residual strain of concrete enters the final stage. The experimental result shows that it is close to the peak strain in the uniaxial stress-strain curve of concrete [33]. Wang et al. [35] obtained fatigue residual strain formula with the material constants through the fatigue test data fitting analysis and is as follows:

4. Life Analysis Method of Concrete Structure

Life estimation was based on FE-SAFE fatigue analysis software. Firstly, the ABAQUS calculation was used to obtain the dynamic response of the cross structure in the process of the train operation in the tunnel. Secondly, the concrete stress-life curve (i.e., S-N curve) of concrete was determined by formulas (16) and (17). Based on FE-SAFE’s nominal stress prediction method, the dynamic response and N-S curve of the cross structure were analyzed, so that the fatigue life of the structure can be obtained.

The fitting formula for the tensile fatigue life and fatigue load curves of concrete [17] was obtained by the uniaxial compression test of concrete, proposed by Tepfers and Kutti: where is the material constant of concrete, the value range is 0.064~0.080, Teng and Wang [21] suggested an average of 0.072; , are the maximum compressive stress value and the minimum compressive stress value of cyclic load, respectively. is the compressive strength of concrete.

The tensile fatigue life and fatigue load curves of concrete were obtained by the fitting formula based on concrete fatigue splitting test [17]:where and are the maximum tensile stress and the minimum tensile stress of the cyclic load, respectively; is the tensile strength of concrete.

5. Damage Analysis Theory

Plastic flow, microcracks, and microvoids are the fundamental reasons of nonlinearity of concrete. From the macroscopic performance, it shows the obvious difference of concrete tensile strength and compressive strength and the residual deformation of concrete [36]. The plastic damage constitutive model of concrete (CDP model) based on continuous medium was adopted to better simulate the nonlinear properties of concrete subjected to external load [37].

According to the energy equivalent principle proposed by Sidoroff [38], the structural damage factor is as follows:where is the elastic modulus when the concrete is damaged.

According to the plastic damage theory of concrete, when the concrete is tensile, the cracking strain is [37]where is the tensile stress of concrete.

When the concrete is compressive, the inelastic strain is [37]where is the tensile stress of concrete.

6. Dynamic Fatigue Damage Analysis Progress and Model of Tunnel Structure

The complex cross structures of the tunnel with the designed service life up to 100 years are subjected to high cycle fatigue problems under the train caused vibration loads. It is uneconomical to calculate the dynamic impact of train on the tunnel every time. Petryna and Krätzig proposed the idea of high cycle structural fatigue [32, 39]. The dynamic response of the train operation in the cross tunnel structure under the vibration load for the first time can be simulated. Then, based on the response law of the structure under the vibration load of train for the first time, the fatigue life of the cross structure can be calculated to obtain the concrete uniaxial cyclic loading constitutive model of the cross tunnel after th train vibration load. Finally, based on the concrete constitutive model after times of train operation, the dynamic response and fatigue cumulative damage of tunnel cross structure after the specific operation years can be simulated. Specific analysis process is shown in Figure 2.

7. Project Overview

Guangzhou-Shenzhen-Hong Kong Railway Passenger Dedicated Line is a fast railway channel connecting Guangzhou, Shenzhen, and Hong Kong, which is an important part of the intercity railway network in the Pearl River Delta. The full length of the Shiziyang subsea tunnel is 10.8 km, and it is the longest and highest standard subsea railway tunnel in China. Guangzhou-Shenzhen-Hong Kong Passenger Dedicated Line Shiziyang Tunnel is located at Dongyong Station, Humen Station interval. The ground layers where the tunnel run through are mainly soil, mucky soil, silty clay and fine sand, coarse sand, weathered and weak weathered argillaceous siltstone, siltstone, and fine sandstone. The planar graph and the vertical sectional profile of the Shiziyang subsea tunnel are shown in Figures 3 and 4.

According to the actual situation of the project, the connecting transverse passage cross structure of Shiziyang railway shield tunnel was selected and studied here. The dynamic response of cross structure was simulated when the marshalling train is running in the A tunnel. It is assumed that the train is running in the main tunnel A. The clear distance between two tunnels is 5.0 m and the design speed is 300 km/h. The buried depth of the selected section of the tunnel is 19.0 m, located in the weak weathered muddy siltstone, topsoil layer covered with lighter silty clay layer, and fine sand layer. The outer and inner diameters of the shield tunnel are 10.8 m and 9.8 m. The lining is assembled in a 7 + 1 block way with a universal wedge ring reinforced concrete single segment. In order to consider the impact of the segment on the structure, the stiffness reduction ring is set at the main tunnel spacing, and the reduction factor is 0.8. The width and height of connecting transverse passage are 4.0 m and 5.0 m, respectively. The length, width, and height of the stratigraphic structure model are 800.0 m, 80.0 m, and 50.0 m, respectively. All the boundaries except for the upper boundary are simulated using a continuously distributed parallel spring-damper system. This boundary treatment can effectively solve the near field fluctuation problem at soil-structure dynamic interaction.

Marshalling train does not consider the connection between the carriages. The train contains 8 carriages with single carriage length of 25.0 m. Each carriage at the front and rear part has two pairs of axles, a total of 32 pairs of axles. The physical and mechanical parameters of surrounding rock, lining concrete, and track are shown in Table 1. The numerical analysis models of cross tunnel are shown in Figure 5.

The profile irregularity of a railway line is one of the essential vibration sources for vehicles and track [40]. Assuming that the high-speed trains run 40 times per day in the tunnel, the vibration load of the train was selected by the measured load curve of train vibration at 300 km/h, as shown in Figure 6.

In order to simulate the variation of the spatial position and the vibration load time in the upper tunnel where the high-speed train is running, the travel speed of 300 km/h was applied to the marshalling train, to simulate the space driving effect of the train, as shown in Figure 7.

Track and train were simulated using linear elastic materials, the surrounding rock was simulated by the elastic-plastic model with damping, and the vibration system damping used Rayleigh damping. According to the above analysis process and the dynamic response analysis model of the cross tunnel structure under the train vibration, the train vibration response and fatigue cumulative damage analysis of the shield tunnel cross structure were carried out.

8. Result Analysis

Firstly, the fatigue life of the tunnel cross structure under the vibration load of the train was analyzed. The stress time-history of the tunnel cross structure obtained by the calculation during the train operation in the tunnel for the first time was introduced into the FE-SAFE software. The logarithmic life distribution nephogram of the cross tunnel structure was calculated as shown in Figure 8.

From Figure 8, the shorter part of the cross structural life is mainly concentrated in the main tunnel A near the middle of the vault and connecting transverse passage on both sides of the wall corner. The minimum logarithmic life of the cross structure is 6.491. The minimum life of the structure is about 3.10 × 106 times.

8.1. Dynamic Response and Cumulative Damage of Main Tunnel
8.1.1. Structural Acceleration

Based on the response law and fatigue life of the structure under the train vibration load for the first time, the fatigue constitutive model of the structure can be obtained. Then, the dynamic response and cumulative damage effect of the cross tunnel structure after operating a certain period of time can be calculated.

Four positions of the tunnel arch bottom were selected as the analysis points, and the analysis points layouts are shown in Figure 9. The analysis points are located at the arch bottom of the tunnels. The longitudinal spacing between the points is 20 m.

The acceleration amplitude of the main tunnel analysis points A1, A2, A3, and B were extracted, after high-speed trains ran in the tunnel 1, 1 × 103, 1 × 104, 1 × 105, 2 × 105, 5× 105, 1 × 106, 1.5× 106 times, as shown in Table 2.

From Table 2, the vibration acceleration amplitudes of tunnel analysis points have a certain degree of increase with the increasing train operation time. The acceleration amplitude of the main tunnel A with the train operation is obviously larger than that of the adjacent tunnel B. The maximum growth amplitude of the tunnel A is 0.21 m/s2.

8.1.2. The Maximum Principal Stress of Structure

The maximum principal stress time-history curve of point A2 of the main tunnel A is shown in Figure 10, when the train is running in a cross shield tunnel for the first time.

Figure 10 shows that when the train is running in the main tunnel for the first time, the closer the train is to the point of analysis, the more intense it vibrates. At 3.7 s, the vibration at the arch bottom A2 is the most violent, and the maximum principal stress reaches 0.94 MPa.

The maximum principal stress curves of the analysis points are shown in Figure 11.

From Figure 11, the maximum principal stress peaks of the tunnel analysis points increase with the increasing train operation time. During the tunnel design period, the maximum principal stress at the analysis points presents a two-stage change characteristic: in the first stage, from 0 to 13.70 years, the maximum principal stress amplitude of the analysis points A1, A2, and A3 increases rapidly. In the second stage, from 13.70 to 102.74 years, the growth rates of analysis point of the main tunnel arch bottom are relatively slow and close to the linear development. The maximum principal stress of the analysis point B in the tunnel operation period is basically linear trend, which is only 17.4% of the point A2. It is concluded that the adjacent tunnel B belongs to passive vibration; therefore, the vibration load of the train has little effect on the adjacent main tunnel B.

Figure 12 shows the development trend of the maximum principal stress of the main tunnel segment. It can be seen from Figure 10 that the maximum principal stress of the main tunnel arch bottom is the largest, and the right analysis point near the transverse passage is the second, and the left analysis point away from the transverse passage is the smallest.

8.1.3. Structural Cumulative Damage

The distribution nephogram of the cumulative tensile damage of the cross structure after tunnel operation 102.74 years is shown in Figure 13.

Figure 13 shows that the cumulative tensile damage is mainly concentrated near the invert of the main tunnel A. The size and range of the cumulative tensile damage value in the middle invert of the main tunnel A are obviously higher than those of other positions. In addition, the cumulative tensile damage develops from near the arch bottom to the side wall of the connecting transverse passage in the middle of the main tunnel A, due to the presence of the connecting transverse passage.

The development trend of the tensile damage of the analysis points C1, C2, and C3 is shown in Figure 14 to further analyze the cumulative damage development of the main tunnel structure.

Figure 14 shows that the cumulative tensile damage of the main tunnel segment increases with the increase of the train operation time. The cumulative tensile damage of the point C3 in the arch bottom is the largest. The cumulative tensile damage of the right analysis point C2 is obviously larger than the left analysis point C1 when the high-speed train operation is 102.27 years, because the right analysis point C2 is connected with the transverse passage, leading to the structural stress concentration.

8.2. Dynamic Response and Cumulative Damage of Connecting Transverse Passage Structure
8.2.1. Structural Acceleration

As can be seen from the above, the connection location is a weak part of the structure due to the stiffness singularity of main tunnel and connected transverse passage. Therefore, the interface of the connected transverse passage and the main tunnel was taken as the analysis section, and the maximum envelope of the acceleration in the time-history range was obtained by extracting the acceleration of the interface of the typical time point (the train first operation, operating for 0.07, 70.68, 49, and 102.44 years), as shown in Figure 15.

Figure 15 shows that the train is running for the first time in the tunnel, the acceleration of the connecting transverse passage mainly focuses on the middle of the transverse passage floor, and the acceleration maximum value is 4.26 m/s2. The minimum value of acceleration is about 0.76 m/s2 in the vault of connected transverse passage. With the increase of the operation time, the acceleration of the middle of the connecting transverse passage floor, both sides of the wall corner and vault increase significantly.

8.2.2. Maximum Principal Stress of the Structure

The maximum principal stress values of the analysis point of the left side wall corner of the transverse passage were extracted, when the train is running for the first time, as shown in Figure 14.

Figure 16 shows that the vibration of the left side wall corner of the transverse passage is more intense when the distance of the train from the connecting transverse passage becomes closer. When the train runs to the analysis point A2, at 3.7 s, the vibration of the analysis point is the most intense, and the maximum principal stress reaches 1.06 MPa. It is concluded that the vibration response of the left side wall corner of the transverse passage is more intense than that of arch bottom of the main tunnel, due to the singularity of the main tunnel and connecting transverse passage stiffness.

The maximum envelope of maximum principal stress in time-history range when the train is running for the first time was obtained by extracting the maximum principal stress at the typical time point interface, as shown in Figure 17.

Figure 17 shows that the principal stress of the connecting transverse passage is mainly concentrated near the side wall corner. The maximum principal stress at the right side wall corner reaches 1.06 MPa. As the operation time increases, the maximum of maximum principal stress moves upward from the skewback along the side wall to upside and its range gradually increases.

8.2.3. Cumulative Damage of Structure

Cumulative compressive and tensile damage of connecting transverse passage was extracted to further analyze damage situation of connecting transverse passage, as shown in Figures 18 and 19.

Figures 18 and 19 show that the cumulative damage of the connecting transverse passage is mainly distributed in both sides of side wall corner. As the operation time increases, the damage develops towards the side wall and the floor; especially after the tunnel operation for 102.74 years, the damage is particularly evident. In addition, the damage size and range of the side wall corner are greater than those of the right side wall corner. The cumulative compressive damage of connecting transverse passage is mainly distributed in the outer wall corner position, while cumulative tensile damage is mainly distributed in the inner wall corner. The size and range of tensile damage of structure are larger than those of the compressive damage. The results show that the tensile damage of transverse passage is the main factor of structure damage.

Considering that the connecting transverse passage is mainly affected by the tensile damage, the cumulative maximum damage values of the transverse passage at the different operation times were extracted and the maximum development curve was shown in Figure 20.

Figure 20 shows that, in the first stage of the train operation of 0~6.85 years, the damage of connecting transverse passage develops rapidly with the duration of 6.8 years. In the second stage of 6.85~68.49 years, the cumulative damage development is relatively gentle and close to the linear development, which can be used to predict structure damage development. In the third stage of 68.49~102.74 years, cumulative damage development is intense. It shows nonlinear development trend, with the duration of 68.5 years. It is concluded that the initial microcracks in the interior of the material caused by the defects of the concrete itself are developed due to the reciprocating vibration load of the train. The longitudinal development is faster in the first stage of the train operation. In the second stage, because the concrete is constrained by other aggregates, the number of microcracks does not increase significantly. But the existing cracks continue to expand, the material strength decreases continuously, and the damage is close to linear accumulation. In the third stage, the microcracks of concrete are interconnected and expanded with each other, and they continue to damage the bond between the aggregate and the mortar, and cracks are rapidly expanding.

9. Conclusions

Considering the driving effects of high-speed train, the vibration fatigue life of tunnel cross structure was calculated using fatigue analysis software. The dynamic response and cumulative damage characteristics of cross tunnel structure of Shiziyang railway shield tunnel at various operation years were analyzed, which meet the requirements of the designed life of 100 years according to the Chinese standard. The following main conclusions are obtained:

(1) According to the latest concrete design code, the uniaxial cyclic loading constitutive model of concrete is proposed by taking into account the factors such as the stiffness degradation of concrete, the strength decrease of concrete, and the increase of fatigue residual strain. The proposed model is suitable for high cycle vibration fatigue analysis of train and can reflect the current commonly used concrete mechanical properties.

(2) The high-speed train ran in the cross tunnel structure and the middle area of the main tunnel in which the train runs were the most violent. The dynamic response of the arch bottom in the middle of the main tunnel opposite the transverse passage was the largest. The vibration response of the hance of the main tunnel near the transverse passage was larger than that of the opposite of the main tunnel hance.

(3) The stiffness singularity between the connecting transverse passage and the main tunnel caused large stress concentration phenomenon at the interface, and the stress and acceleration were relatively large at the interface. The maximum principal stress of the transverse passage mainly appeared near the side wall corner, while the maximum acceleration mainly appeared near the middle of the transverse passage floor.

(4) The maximum principal stress and acceleration of the connecting transverse passage increased with the increasing train operation years. The acceleration increased significantly in the middle of the transverse passage floor, both sides of side wall corner and vault, while the maximum principal stress developed obviously only in both sides of side wall corner.

(5) The cumulative damage of connecting transverse passage was mainly distributed in both sides of side wall corner. With the increasing operation years, the accumulative damage developed towards the side wall and floor of connecting transverse passage. The damage value and range of the left side wall corner were larger than the corresponding position of the right side wall.

(6) The cumulative compressive damage of connecting transverse passage was mainly distributed in the outer wall corner position, while cumulative tensile damage was mainly distributed in the inner wall corner, which was close to the linear development. The size and range of tensile damage of structure were larger than those of the compressive damage.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

Financial support for this work provided by the National Key R&D Program of China (Grant no. 2016YFC0802205) and the National Science Foundation of China (Grants nos. 51178400, 51278425, and U1361210) is gratefully acknowledged.