Abstract
In order to study the nonlinear characteristics of self-excited aerodynamic forces of bluff body bridge section with the change of motion parameters, a numerical wind tunnel is established by the dynamic mesh technique of computational fluid dynamics (CFD). A state-by-state forced vibration method is used to identify the self-excited aerodynamic forces of single degree-of-freedom (DOF) heaving and pitching motion. Fast Fourier transform (FFT) is adopted to obtain frequency-domain data for analysis. The reliability of the obtained aerodynamic results is verified by wind tunnel tests. The results show that the high-order harmonic components are found in the self-excited aerodynamic forces of semiclosed box deck section, which are more significant in aerodynamic lift than in aerodynamic moment. The proportion of aerodynamic nonlinear components increases with amplitude. The effect of amplitude on the nonlinear components of heaving motion is generally higher than that of pitching motion, and aerodynamic moment is highly sensitive to the increase of vertical amplitude. The variation of the nonlinear components of the deck section with frequency is not a simple monotonic relationship, and there is a stationary point at 10 Hz frequency. The existence of wind attack angle makes the proportion of nonlinear components reach more than 30% and greatly increases the proportion of second harmonic. In addition, the high-order harmonic components, which are not integer multiples, are found at large amplitude and positive angle of attack.
1. Introduction
In recent years, with the progress of engineering technology, the bridge has developed in the direction of large span, light weight, and flexibility [1]. As a result, the sensitivity of bridges to wind increases, and the problem of aerodynamic stability caused by the nonlinear characteristics of bridges becomes more and more significant. Bleich [2] is a critical dynamic instability phenomenon and hazardous to structural reliability and safety. In 1935, for the need of studying wing flutter, Theodorson and Mutchler [3] derived an analytical expression of flutter self-excited force for an ideal plate. The current bridge flutter analysis is based on the linear unsteady self-excited force model proposed by Scanlan and Tomko [4]. The self-excited force model is based on the assumption of linear and small amplitude motion and ignores the influence of bridge aerodynamic shape of the bluff body and other factors such as amplitude, frequency, and wind attack angle. In fact, bridge deck sections are bluff bodies, and their self-excited aerodynamic forces inevitably have nonlinear characteristics [5–10]. The nonlinear effect of the self-excited force causes an oscillation phenomenon called Limit Cycle Oscillation (LCO) [5, 11, 12] with lower starting wind speed and constant amplitude oscillations range. The LCO phenomenon cannot be explained by the above classical flutter self-excited force model. The nonlinear components of self-excited aerodynamic forces of deck sections can affect the applicability of Scanlan’s linear model and possibly cause an unsafe deviation of flutter critical wind speed based on the linear theory.
At present, the research on nonlinear self-excited aerodynamic forces of bridges mainly focused on the nonlinear dependence of aerodynamic forces on amplitude [13, 14] and the establishment of nonlinear aerodynamic models [15–18]. Many scholars such as Falco et al. [8] and Lee and Su [19] found high-order harmonics in wind tunnel tests on particular bridge decks, which confirmed the existence of nonlinear aerodynamic forces. Through the forced vibration method of wind tunnel test [20, 21], the obvious high-order harmonic components of a thin flat plate, a streamlined box deck section, and a twin-side-girder section were found by Chen and Yu [22], and the high-order harmonic components of bluff body sections were more significant. The amplitude range of obvious high-order harmonics was defined by Long et al. [23] by studying the high-order harmonics components in aerodynamic forces of a streamline box deck section. When studying the nonlinear aerodynamic forces of a rectangular section with an aspect ratio of 5 : 1, Lin et al. [24] proposed the amplitude limit of torsional motion where the proportion of higher harmonic components exceeds 5%. The wind tunnel test is the main method to study the aeroelastic stability of long-span bridges, but the forced vibration wind tunnel test involves many difficulties, such as the design and arrangement of driving device and the control of driving frequency. Thanks to the recent fast increase of computer power and improvement of numerical codes, CFD is becoming an increasingly effective tool for studying unsteady aerodynamics. It can overcome the above difficulties and also facilitate the setting of motion parameters and boundary conditions, although much work is still necessary before it can be considered fully alternative to the wind tunnel test. Larsen [25] discusses the aerodynamic and structural response of four generic cross-section shapes using two-dimensional numerical simulation. Huang and Liao [26] analysed the nonlinear self-excited forces of thin plate by using the CFD method integrated with the continuous wavelet transform method. Mannini et al. [27] investigated the influence on the self-excited forces of the amplitude of oscillation, mean angle of attack, and Reynolds number of a trapezoidal box girder bridge section by two-dimensional numerical simulations. Zhu et al. [28] found that the nonlinear proportion of the streamline box deck section increased with the increase of amplitude by using the step-by-step forced vibration method of numerical wind tunnel, and the synchronous relationship was stronger in torsional state. Using the numerical simulation method, Wang et al. [29] verified the superposition of self-excited aerodynamic forces of a streamlined box deck and the feasibility of step-by-step forced vibration with multifrequency.
Hence, for the nonlinear self-excited aerodynamic forces of bridges, most of the current researches focus on the streamlined deck section and less on other bridge sections with more significant characteristics of the bluff body. In this paper, a semiclosed box deck with poor streamline property is selected as the research object. It is more and more used in the design of modern long-span cable-stayed bridge because of its advantages of lightweight, convenient construction lifting, and economy. The obvious LCO phenomenon appears in the deck section, which is more conducive to capturing nonlinear characteristics and research. The influence of amplitude, frequency, and angle of attack on nonlinear aerodynamic forces of the bridge section is analysed by a combination of wind tunnel test and numerical simulation. By revealing the characteristics of aeroelastic response and associated influencing factors, this study is helpful to understand the dependence of nonlinear flutter on different motion parameters and vibration forms, provide a reference for the establishment of nonlinear flutter calculation method, provide direction for scientific evaluation standards for LCO, and promote more reasonable and perfect aerodynamic stability evaluation method of bridges.
2. Materials and Methods
The classic flutter theory of Scanlan et al. believes that whether it is a streamlined section or a blunt section, if only the torsion and vertical motion of the structure are considered, self-excited forces resulting from structural motion can be expressed linearly by 8 flutter derivatives:where L and M are the aerodynamic lift and moment, respectively; h and α are the heaving motion and pitching motion, respectively; ρ is the air density; U is mean wind speed of oscillation; B is the width of decks; K = ωB/U is the reduced frequency; ω is the circular frequency of oscillation; and () are the flutter derivatives. The superscript “·” is the symbol for derivation.
Equations (1) and (2) imply the following four assumptions:(1)There is a linear relationship between the self-excited aerodynamic forces and the motion parameters of bridge decks. If a deck makes a simple harmonic motion at a fixed frequency, its self-excited aerodynamic forces should be a simple harmonic of the same frequency, and there should be no high-order harmonic.(2)Flutter derivatives are functions of reduced velocity, which is uniquely determined by the aerodynamic shape of bridge decks, independent of other motion factors such as amplitude and frequency.(3)The self-excited aerodynamic force can be superimposed.(4)There is only small vibration in bridge decks, so there is no influence of attack Angle of incoming wind by default.
Hypothesis (3) has been studied by Matsumoto et al. [21]. This paper focuses on the validity of hypotheses (1), (2), and (4)
Generally, if a bridge deck makes simple harmonic vibration for a fixed period in an incompressible low-speed uniform flow field, the fluid will be disturbed by the boundary of the same period, then the macroflow field will also change according to this period, and the flow field forces on the bridge deck can be regarded as a function of this period.
According to the idea of the Fourier series, a function F (t) with a period of T can be represented by trigonometric series composed of sine function, denoted aswith and , , and (n = 1, 2, 3, …) are constants.
The periodic self-excited aerodynamic forces can be divided into the superposition of self-excited aerodynamic forces of each integer multiple frequencies by Fourier transform. Therefore, the time courses of aerodynamic forces are converted to the frequency domain using a fast Fourier transform, and then the high harmonic components of the self-excited aerodynamic forces of bridge decks can be obtained and analysed.
2.1. Numerical Simulation
2.1.1. Deck Model
Based on the assumptions, we studied a semiclosed steel box composite beam section. The deck shape is taken from an actual cable-stayed bridge named Jiaojiang Second Bridge in Taizhou City, Zhejiang Province. The actual bridge has a width and height of 42.58 m and 3.5 m, respectively. The bridge width to height ratio is 12.17. The section model is shown in Figure 1, and 1/80 is selected as the size of the model for numerical calculation. The deck section model has width and height of 0.532 m and 0.044 m, respectively. The shear center is located at 0.393 times the height of the upper surface. The influence of railings, crash barriers, and other ancillary facilities is excluded in the numerical simulation. Generally, the aerodynamic characteristics of the deck section with sharp edges and corners and obvious passivation are not sensitive to Reynolds number [30].
2.1.2. Grid Division and Solution Setting
In the present paper, the forced vibration of vertical and torsional two DOF of the deck section is conducted. Therefore, considering the calculation efficiency and quality comprehensively and referring to the relevant studies of other scholars [25–29], a two-dimensional wind tunnel model is established. The grid is drawn by professional preprocessing software ICEM, and the computing domain of the grid is shown in Figure 2. Blockage ratio refers to the ratio of the projected area of the section model to the cross-sectional area of the wind channel. The blockage ratio at zero angles of attack is approximately 1.03%, which meets the requirement of 3% blockage ratio of the model [31].
In order to improve the accuracy and efficiency of calculation, the method of “rigid moving region + dynamic grid region + static grid region” is adopted to divide the grid. The rigid moving region moves with the deck section model to ensure the good quality of the near-wall grid. The dynamic grid region is updated and remeshed with the movement of the deck section model. The static grid region is always kept unchanged, and appropriate encryption is carried out at the end of the flow field to ensure calculation accuracy. Appropriate magnification is adopted for the whole grid from inside to outside to control the total amount of the grid and improve the efficiency of simulation calculation. The turbulence model of numerical calculation is the SST k-ω model, which needs to calculate the near-wall region. Therefore, it is necessary to reasonably control the height of the first-layer grid and encrypt the near-wall grid, making y+ ≤6 in general [32]. In this paper, the y+ value is set as 1, the height of the first-layer grid near the wall is calculated by formula (4) as = 6 × 10−5 m, and the normal growth rate is set as 1.01.where is the wall friction velocity; is the wall shear stress; is the coefficient of kinematic viscosity; is the distance to wall.
The commercial software ANSYS FLUENT is used for numerical calculation. The left boundary of the calculation domain is the velocity inlet, the turbulence intensity is 0.5%, and the turbulence viscosity ratio is 2. The right boundary is the pressure-outlet with the reference pressure of zero, and the upper and lower boundaries are symmetrical boundary conditions. The SIMPLE algorithm of velocity and pressure decoupling is selected for calculation. The second-order scheme is used to discretize the pressure equation, and the second-order upwind scheme is selected to discretize the momentum equation, turbulent kinetic energy equation, and turbulent dissipation rate equation.
The forced vibration is realized based on the dynamic grid technology in FLUENT software. The grid update of forced vibration was realized by specifying the vertical and torsional vibration values to DEFINE_CG_MOTION macro. A combination of two dynamic grid updating methods of “Smoothing + Remeshing” is adopted to ensure the efficiency and quality of grid updating.
2.1.3. Independence Verification
The drag coefficient, CD, and lift coefficient, CL, of the semiclosed box girder section are defined aswhere FD and FL are the resistance and lift acting on the bridge girder section, respectively; H is the height of deck section; B is the width of deck section; U is mean wind speed.
Four sets of computational grids are divided according to the simulation object to verify the independence of grid quantity. It can be seen from Table 1 that when the number of grids reaches 505334, the influence of encrypted meshes on the drag and lift coefficients is extremely small. It indicates that the grid resolution of 505334 meets the calculation requirements. Three groups of time steps are selected for independence verification. Table 2 shows that the result deviation among the three selected time steps can be ignored, but too small time step may cause the deviation of higher harmonic data under forced vibration, which is adverse to the analysis of nonlinear aerodynamic forces. In conclusion, the grid numbers of 505334 are selected by comprehensively considering the computational efficiency and accuracy, as shown in Figure 3. The time step of 0.0001 s is used for numerical calculations.
(a)
(b)
(c)
2.2. Wind Tunnel Test
In order to investigate the reliability of the grid and ensure the validity of the numerical calculation results, static forces wind tunnel test was carried out because of the lack of relevant data on the deck section of Second Jiaojiang Bridge for comparative verification. The wind tunnel test was conducted in a wind tunnel laboratory of a closed loop located in the Northeast Forestry University, China. The wind speed can be adjusted from 5 m/s to 70.5 m/s. The test section is approximately 1.5 m long, 1 m high, and 0.9 m wide, and the flow is nearly laminar with a uniform velocity profile at the inlet and with turbulence intensity of less than 0.5%. The flow field is stable and reliable.
2.2.1. Overview of Wind Tunnel Test
Comparability is taken into account for the results of the wind tunnel test and numerical simulation. Therefore, the physical parameters of the deck section model in the wind tunnel test are consistent with those in the numerical simulation; that is, the scale ratio between the deck section model and the real bridge is 1 : 80. The model has a streamwise width, B, of 0.532 m, a height, H, of 0.044 m, and a span, L, of 0.82 m. Two identical end plates are installed on both sides of the deck section model to enable a close to the bidimensional flow field. The force test is completed by a six components balance, which is placed on a rigid base with adjustable height, as shown in Figure 4. The bottom of the model is raised a sufficient distance from the wind tunnel wall to ensure uniform flow field, avoiding the influence of wind speed profile as far as possible. The section model placed in the wind tunnel is shown in Figure 5.
2.2.2. Wind Tunnel Test Results
According to the requirements of this paper, force tests are carried out at −3°, 0°, and +3° angles of attack, respectively. In order to avoid the accidental error of the experiment, three groups of wind tunnel tests were carried out. The sampling time of each group is 60 s, and the sampling frequency is 1000 Hz. As shown in Figure 6, the average of three experimental results was taken as a final result, which proves the reliability of the numerical simulation grid in this paper.
3. Results and Discussion
By loading User Defined Function (UDF) in fluent, single pitching and heaving motions of the deck section are carried out, respectively, and the torsional and heaving displacements are harmonic functions with the same frequency f:
In order to avoid the influence of the Reynolds number effect, different reduced wind velocities () are obtained by changing the vibration frequency f, and the mean wind speed is 24 B/s (B = 0.532 m). The specific test conditions are shown in Table 3. The aerodynamic time history curve obtained by calculation was taken 5 stable periods for FFT to obtain spectra, which were used for the study of the variation law and nonlinear characteristics of aerodynamic forces.
3.1. Effect of Amplitude on Nonlinear Aerodynamic Force
Firstly, the influence of amplitude on the aerodynamic force of the semiclosed box deck section at 0° angle of attack is analysed at the reduced velocity of 6 (according to formula , the driving frequency is 4 Hz; i.e., the fundamental frequency of aerodynamic forces is 4 Hz).
The time history curves of the aerodynamic forces of five stable periods are shown in Figures 7 and 8. From Figure 7, in the state of small amplitude heaving motion, the aerodynamic time history curves have good harmonic characteristics, and slight waveform distortion appears in the time history curve of aerodynamic lift. Figure 8 shows that, in the state of small amplitude pitching motion, varying degrees of clear deformation can be observed in the time history curves, which intuitively indicates that nonlinear responses exist in the vibration of bluff body section.
(a)
(b)
(a)
(b)
For single frequency and single degree-of-freedom heaving motion, as shown in Figures 9 and 10, there are integer multiple higher-order harmonic components other than fundamental frequencies in the aerodynamic spectrum of the deck section. Although the higher harmonic components are not obvious under small amplitude, the second harmonic is first exposed along with the increase of amplitude, and the proportion of third and higher harmonics is also increased. The third harmonic component is significant as the amplitude continues to increase, which shows the importance of the third harmonic under this degree of large vibration. It is notable that the fourth and higher harmonics also appear obviously under large amplitude, which is different from streamlined box deck sections. (In this paper, the term “proportion” represents the proportion of the nonlinear aerodynamic components in the total aerodynamic components, and the total aerodynamic components include both linear and nonlinear components.)
(a)
(b)
(a)
(b)
For the single frequency and single degree-of-freedom pitching motion, it can be seen in Figures 11 and 12 that the aerodynamic components of the section are also dominated by the fundamental frequency and have high-order harmonic components. The high-order harmonics in small amplitude motion are mainly the second harmonic, but its components are not clear. The high-order harmonic components increase with the amplitude. The third harmonic is the main part of the aerodynamic lift and the second harmonic is the main part of the aerodynamic moment as the torsional amplitude reaches 10°. As the amplitude continues to increase, the second harmonic (8.4%) is dominant in the aerodynamic lift, the third harmonic (8%) is dominant in the aerodynamic moment, and the fourth and higher harmonic components of both are also obvious. It is indicated that the changes of aerodynamic lift and aerodynamic moment are not completely synchronized. It is notable that when the torsional amplitude reaches 20°, some small components, which are not integer frequency multiplication, appear in the aerodynamic spectrum.
(a)
(b)
(a)
(b)
In Figures 13 and 14, the second to fifth harmonic components of aerodynamic lift and aerodynamic moment are shown. “Total” represents the total content of high-order harmonic components. The sum of the first to fifth harmonics is used as the cardinal number for normalized statistics. The proportion of high-order harmonic components of aerodynamic lift is higher than that of the aerodynamic moment in both motion cases. The difference of high-order harmonic components between them is larger under small amplitude vibration. It is revealed that the nonlinearity of the aerodynamic lift force is stronger.
(a)
(b)
(c)
(a)
(b)
(c)
Under the condition of vertical amplitude h = 0.2B and torsional amplitude α = 20°, the proportion of nonlinear components of the semiclosed box deck section reaches 9.81% and 26.26%, respectively, which is higher than the research results of Chen and Yu [22] on thin flat plate, streamlined box deck, and twin-side-girder. When the torsional amplitude α = 2°, the proportion of high-order harmonic components has exceeded 5%, which is lower than the limit of 10° torsional amplitude proposed by Lin et al. [24] when the proportion of high-order harmonics exceeds 5% in the study of streamlined box deck. It shows that the blunter profile will enhance the nonlinear characteristics of the deck section. The nonlinear components of the semiclosed box deck section are higher than those of other deck sections, which can also explain the result that the semiclosed box girder has more obvious soft flutter characteristics found by Zhu and Gao [33] in the analysis of typical bridge deck sections.
The variation of aerodynamic high-order harmonic when the amplitude increases at the same magnification is used as a reference to compare the high-order harmonic variation with amplitude between single frequency and single DOF pitching and heaving motion, as shown in Figure 15. The influence of amplitude on the aerodynamic nonlinearity of pitching motion is generally higher than that of heaving vibration, which is similar to the research results of a streamline box deck by Zhu et al. [28]. However, when the amplitude increases to a certain extent, the sensitivity of aerodynamic moment to vertical amplitude increases significantly.
(a)
(b)
3.2. Effect of Frequency on Nonlinear Aerodynamic Force
In order to explore the influence of vibration frequency on the nonlinear characteristics of self-excited aerodynamic forces of the semiclosed box deck section, the analysis is carried out under the two motion states of vertical amplitude 0.02B and torsional amplitude 2o. Since the high-order harmonic components are mainly concentrated in the first 5 times wave, the sum of the proportions of 2nd to 5th harmonics is taken as the proportion of high-order harmonics.
Figure 16(a) shows that, under the state of single heaving motion (f = 2 to 12, i.e., U = 12 to 2), the proportion of nonlinear components decreases with the increase of reduced velocity as the driving frequency is lower than 10 Hz, and the proportion of nonlinear components increases with the reduced velocity as the driving frequency is higher than 10 Hz. The maximum proportion of high-order harmonic components of aerodynamic lift and aerodynamic moment is reached at the driving frequency of 10 Hz; that is, there is a stationary point in the curve of the proportion of high-order harmonic components with the change of driving frequency.
(a)
(b)
Figure 16(b) shows that, under the state of single pitching motion (f = 2 to 12, i.e., U = 12 to 2), the high-order harmonic components of aerodynamic lift reach the maximum at 10 Hz, and the high-order harmonic components of aerodynamic moment reach the maximum at 8 Hz. Thus, there are stagnation points in the graph line of the proportion of higher-order harmonics, which is not a simple monotonic relationship. The result explains the differences among monotonic in the conclusions of different scholars when studying the nonlinear aerodynamic variation characteristics of bridges under different reduced velocity.
3.3. Effect of Wind Attack Angle on Nonlinear Aerodynamics
At the reduced velocity of 6 (i.e., the driving frequency is 4 Hz), the influence of wind attack angle on the nonlinear aerodynamic forces of the semiclosed box girder is studied under the large amplitude heaving motion (h = 0.2B) and large amplitude pitching motion (α = 20°) of single DOF. Because more than the 5th-order harmonics are found in the spectrum analysis under the condition of large amplitude motion, the 1st- to 10th-order harmonics are taken as the aerodynamic research object, and the 2nd- to 10th-order harmonics are taken for local amplification analysis.
Figures 17 and 18 show that the high-order harmonic components are very obvious under the single heaving motion of large amplitude. When the wind angle of attack exists, the second-order harmonic component increases significantly, especially when the angle of attack is −3°. The importance of the second-order harmonic in large vertical vibrations with angles of attack is revealed. For aerodynamic lift, the high-order harmonic components accounted for 11.82% at 0° angle of attack, 17.61% at +3° angle of attack, and even 32.9% at −3° angle of attack. For aerodynamic moment, the high-order harmonic components are 7.28% at 0° angle of attack, 13.54% at +3° angle of attack, and the proportion is the highest at −3° angle of attack, reaching 20.96%. It is evident that the nonlinear components of aerodynamic forces are enhanced significantly in large vertical vibration with a negative wind attack angle.
Figures 19 and 20 show that the obvious high-order harmonic components can be obtained under the single pitching motion of large amplitude. When the angle of attack exists, the second harmonic is significantly increased, which is in line with the law of heaving motion state. For aerodynamic lift, the +3° angle of attack has a more significant impact on the second-order harmonic. However, the proportion of high harmonic components is reduced by wind attack angle, which is against the law of heaving motion state. The high-order harmonic components are 37.46% at 0° angle of attack and 37.1% at the +3° angle of attack, and the −3° angle of attack is reduced to 32.62%. For aerodynamic moment, the −3° angle of attack significantly enhances the second harmonic. At the 0° angle of attack, the high-order harmonic components are 25.7%, the −3° angle of attack increases to 29%, and the +3° angle of attack increases to 34.24%. In addition, noninteger high harmonics of aerodynamic forces are found at 0° angle of attack. The positive angle of attack makes these noninteger high-order harmonics more significant, and the negative angle of attack reduces them to almost disappear.
4. Conclusions
The aerodynamic nonlinear characteristics were presented and discussed on a semienclosed box deck section model by numerical simulations and wind tunnel test. The effects of amplitude, frequency, and wind attack angle on the nonlinear aerodynamic forces components of the section are studied, and the main findings are shown as follows:(1)The different motion amplitudes of single frequency and single degree-of-freedom vibration are controlled by a step-by-step forced vibration method. It is clear that as amplitude increases, the proportion of section aerodynamic high-order harmonic components increases; that is, the proportion of bridge nonlinear components increases. Under vertical heaving and pitching motion, the maximum value reached 9.81% and 26.26%, respectively. The nonlinear characteristics of aerodynamic lift are more obvious than that of aerodynamic moment, and the difference is 5 times in small amplitude heaving motion. The variation degree of single frequency and single degree-of-freedom pitching motion with amplitude is generally more obvious than that of heaving motion. It should be noted that the high-order harmonic of aerodynamic forces under pitching motion has components that are not integer multiples when the amplitude is large enough.(2)The nonlinear aerodynamic forces of the deck section are affected by motion frequency, but the relationship between them is not simple monotone or linear, and its variation curve has stationary points at about 10 Hz. The proportion of nonlinear components peaks at this point.(3)The existence of wind attack angle greatly increases the second-order harmonic component, which makes the second-order harmonic particularly important to the nonlinear aerodynamic forces of the deck section. Under the pitching motion of −3° angle of attack, the component of the second harmonic can be increased to 6.6 times. In most cases, the existence of wind angle of attack will strengthen the nonlinear characteristics of the deck section. However, the proportion of nonlinear components of aerodynamic lift decreases under pitching motion when the angle of attack exists.(4)As a significant bluff body bridge section, the semiclosed box deck section has obvious high-order harmonic components. Under the condition of large amplitude with angle of attack, the proportion of aerodynamic nonlinear components can reach more than 30%. It is shown that the bluff body bridge has more significant nonlinear characteristics and is more prone to nonlinear flutter (LCO phenomenon) theoretically. When establishing a nonlinear self-excited aerodynamic model, the influence of motion parameters on self-excited force parameters should be properly considered according to the research results, and the nonlinear components should be considered in the analysis of bridge aeroelastic response.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
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 Natural Science Foundation of Heilongjiang Province (no. QC2014C043) and the Fundamental Research Funds for the Central Universities (no. 2572015CB31).