Loading...

#### Table of Content

18 March 2020, Volume 52 Issue 2
Theme Articles on ”Mechanical Behaviors of Disordered Solids”
DOES STRUCTURE DETERMINE PROPERTY IN AMORPHOUS SOLIDS? 1)
Wang Yunjiang, Wei Dan, Han Dong, Yang Jie, Jiang Mingqiang, Dai Lanhong
2020, 52(2):  303-317.  DOI: 10.6052/0459-1879-19-368
Asbtract ( 599 )   HTML ( 38)   PDF (4835KB) ( 428 )   PDF（mobile） (4835KB) ( 68 )
Figures and Tables | References | Related Articles | Metrics

The mechanical properties and plastic deformation mechanisms of crystalline solids are mainly determined by their structural defects, e.g., the motion of the versatile dislocations. However, how structures determine properties in non-crystalline solids remains as a major unsolved issue in both solid mechanics, materials sciences, as well as condensed matter physics. Structure determines property is the traditional paradigm of materials science. Following this rule, there are vast experimental characterizations, theoretical studies, and computer simulations appeared in the literature, trying to establish a one-to-one correspondence between a specific structural feature with a unique dynamic property in the general amorphous solids. However, up to date, people gain very little understanding of the structure-property relationships in amorphous solids, not to mention whether there exists any hidden rule behind the structure-property relationships. For this purpose, we focus on the unique features of deformations mechanisms in amorphous solids as well as their microstructure characteristics. Thorough proper samplings of the activation energies of the excitation of these structural parameters by an advanced molecular dynamics technique, we are trying to quantitatively assess the validity of simple short-range structures and medium- to long-range structures in determination of their properties. This is done by examination of the possible correlation between parameters of structures with their activation energies, which implies the level of difficulty in activation of the events. By this we find that the hidden governing rule of structure-property relationship in amorphous solids involves a critical role of spatial autocorrelation length of the specific structural parameter. Constraint is more relevant than geometry itself. If only one structural descriptor presents spatial autocorrelation length up to sub nanometer level, it might effectively predict the mechanical property of amorphous solids; otherwise, the short-range local structures lacking such correlation length fails to predict property. Furthermore, we present a general metric to assess the utilities of structures in determining functions of the amorphous solids, which can be served as a screening rule to seeking for effective structures in amorphous solids.

REVIEW ON THE DEFORMATION BEHAVIOR AND CONSTITUTIVE EQUATIONS OF METALLIC GLASS MATRIX COMPOSITES 1)
Zhang Juan, Kang Guozheng, Rao Wei
2020, 52(2):  318-332.  DOI: 10.6052/0459-1879-20-038
Asbtract ( 350 )   HTML ( 15)   PDF (25179KB) ( 179 )   PDF（mobile） (4835KB) ( 19 )
Figures and Tables | References | Related Articles | Metrics

Metallic glass and metallic glass matrix composites have good application prospects because of their excellent mechanical properties, and now more and more researches have been conducting on them. The deformation behavior, toughening mechanism and constitutive relationship of metallic glass matrix composites are summarized and reviewed in this paper, based on the existing research results in literature by other groups and the latest work done by the authors. Firstly, the research progress in the deformation behavior, failure mechanism and constitutive relation of metallic glass in recent decades is briefly reviewed. Then, the state-of-the-arts in the deformation behavior and failure mechanism of metallic glass matrix composites are introduced from the aspects of experiments and numerical simulation, and the plastic deformation, toughening mechanism and their correspondent influencing factors of metallic glass matrix composites are summarized. Furthermore, the existing studies on the constitutive equations of metallic glass matrix composites are reviewed, with emphasis on the application of homogenization method in this field. In addition, a two-stepped homogenization method proposed by the authors is introduced in more details as a representative approach, and then the constitutive model established on the two-stepped homogenization method and with a help of a failure criterion obtained by introducing a concentration of nano-voids as an internal variable is addressed. The deformation and failure behavior of metallic glass matrix composites are predicted reasonably by the proposed constitutive model. Finally, the research progress of this field is briefly summarized, and some future topics are suggested.

MECHNICAL PROPERTIES AND BEHAVIORS OF HIGH ENTROPY ALLOYS 1)
Li Jianguo,Huang Ruirui,Zhang Qian,Li Xiaoyan
2020, 52(2):  333-359.  DOI: 10.6052/0459-1879-20-009
Asbtract ( 740 )   HTML ( 45)   PDF (85411KB) ( 742 )
Figures and Tables | References | Related Articles | Metrics

High-entropy alloys (HEAs) are a class of new metallic materials that have revolutionized alloy design over the past ten years. Unlike conventional alloys with one and rarely two base elements, HEAs contain multiple principal elements (at least four principal elements) with equal or nearly equal atomic concertation to promote the formation of simple solid solution phases. Due to the presence of multiple principal elements, multiple deformation mechanisms (including dislocation activities, deformation twinning, and phase transformation) activate during deformation of HEAS. Therefore, HEAs usually exhibited many excellent mechanical properties, such as ultrahigh hardness, high tensile strength, good ductility, high thermal softening resistance, remarkable irradiation resistance, and good wear resistance. HEAs are thought to be the most promising structure materials and have attracted tremendous attention over worldwide in the fields of solid mechanics and material sciences. In this review paper, we first briefly introduce the unique and complicated microstructural features of HEAs, i.e. HEAs have both chemically short-range orderings and severe lattice distortion. Then, we review the recent experimental studies on mechanical properties, behaviors and deformation mechanisms of HEAs with face-centered cubic, body-centered cubic, hexagonal close-packed, dual or meta-stable phases. We also mainly emphasize some effective strengthening and toughening strategies, including solid solution, grain refinement, second phase or precipitation. We further summarize some advanced atomistic simulations/modelling on microstructures, mechanical properties and deformation of various HEAs. Finally, we address a list of open problems and challenges for the future studies about design, fabrication and mechanics of HEAs, and provide some important mechanistic insights into design and fabrication of HEAs with excellent mechanical properties and performances.

DYNAMIC RELAXATION CHARACTERISTICS AND HIGH TEMPERATURE FLOW BEHAVIOR OF ZR-BASED BULK METALLIC GLASS 1)
Hao Qi, Qiao Jichao, Jean-Marc Pelletier
2020, 52(2):  360-368.  DOI: 10.6052/0459-1879-20-004
Asbtract ( 417 )   HTML ( 17)   PDF (699KB) ( 209 )
Figures and Tables | References | Related Articles | Metrics

Dynamic mechanical relaxation processes of amorphous alloys are very important to understand plastic deformation, glass transition phenomenon, diffusion behavior and crystallization. How to establish the correlation between mechanical properties and mechanical relaxation modes is one of key issues. In the current research, with the help of dynamic mechanical analysis (DMA), dynamic mechanical behavior of Zr$_{50}$Cu$_{40}$Al$_{10}$ bulk metallic glass from room temperature to supercooled liquid region was probed. In parallel, based on the uniaxial tensile tests, high-temperature flow behavior of Zr$_{50}$Cu$_{40}$Al$_{10}$ metallic glass around glass transition temperature were investigated. Dynamic mechanical behavior and high temperature deformation behavior were discussed in the framework of quasi-point defects theory. The results demonstrated that main $\alpha$ relaxation process of metallic glass can be well described by the quasi-point defects theory. Based on internal friction of Zr$_{50}$Cu$_{40}$Al$_{10}$ metallic glass, activation energy of elementary movement of atoms $U_\beta$ is 0.63 eV. In addition, correlation factor $\chi$ corresponding to concentration of the quasi-point defects in solid glass remains almost constant below the glass transition temperature. When the temperature above the glass transition temperature, the correlation factor $\chi$ increases by increasing the temperature (below the crystallization temperature). Finally, high temperature flow behavior in tensile mode near the glass transition temperature of Zr$_{50}$Cu$_{40}$Al$_{10}$ metallic glass was studied. The normalized viscosity decreases with increasing strain rate at low temperatures or high strain rates, indicating a non-Newtonian flow behavior. Whereas Newtonian flow behavior is observed at higher temperatures and lower strain rates. The apparent viscosity is affected by temperature and strain rate. High-temperature flow behavior of Zr$_{50}$Cu$_{40}$Al$_{10}$ metallic glass was described by stretched exponential function and free volume theory. Specifically, experimental master curve of the high temperature flow behavior of metallic glass is in good agreement with the prediction of the quasi-point defects theory, which provides a new insight on understanding of viscous effects during high temperature deformation of solid glasses.

PREDICTION OF SHEAR TRANSFORMATION ZONES IN METALLIC GLASSES BASED ON LAPLACIAN OF ATOMIC VOLUME 1)
Shi Ronghao, Xiao Pan, Yang Rong
2020, 52(2):  369-378.  DOI: 10.6052/0459-1879-19-369
Asbtract ( 269 )   HTML ( 11)   PDF (20609KB) ( 81 )
Figures and Tables | References | Related Articles | Metrics

Shear transformation zone (STZ), as a basic characteristic unit of plastic events in metallic glasses (MGs), has been widely accepted by researchers, but the source of its origin and activation mechanism are still controversial. Deformation behaviours of Cu$_{64}$Zr$_{36}$ MGs under simple shear loadings are investigated using molecular simulation method in this paper. The results indicate that the activation locations of STZ are related to the initial configuration of MGs. Though the field of atomic volume and its gradient are a direct representation of the local atomic structural heterogeneity of MGs, they lack an obvious correlation to the regions of STZ activation. A new local structural parameter $\xi$ is proposed in this paper based on the initial configuration of MG to predict the potential regions of STZ. $\xi$ is the product of two factors: the Laplacian of atomic volume field (AVF) and the absolute difference between components of the gradient of AVF. Vectors of the AVF gradient present a distribution pattern of pointing inside if the Laplacian of AVF is negatively large, representing the localized soft regions in MGs. The absolute difference of AVF gradient components is used to select different patterns of the AVF gradient distribution. Furthermore, the relationship among structural parameter $\xi$, nonaffine displacement and shear localization is established, revealing that only certain patterns of AVF gradient distribution would lead to nonaffine displacements field strengthening shear localization, which is more likely to result in activation of STZs. The correlation analysis shows that the averaged spatial correlation index of $\xi$ and STZ is larger than 78%, so $\xi$ can be used as an effective parameter for predicting the activation regions of STZs in MGs. Moreover, the ideology of using Laplacian of local AVF in predicting potential STZ regions in MGs would bridge the analysis between atomic simulations of MGs, the mechanism of STZ activations and the traditional mechanical theory.

SHEAR-BAND DYNAMICS IN METALLIC GLASSES 1)
Dong Jie, Wang Yutian, Hu Jing, Sun Baoan, Wang Weihua, Bai Haiyang
2020, 52(2):  379-391.  DOI: 10.6052/0459-1879-19-378
Asbtract ( 294 )   HTML ( 11)   PDF (13388KB) ( 139 )
Figures and Tables | References | Related Articles | Metrics

Shear banding is a localized plastic deformation form that wildly exist in amorphous systems, and plays the decisive role in controlling the failure and plasticity. Due to the localization in space and instantaneous in motions, clarifying the map of shear bands is still a challenge mission. For traditional amorphous systems as rocks, colloid, oxide glasses and polymers, the studies on shear bands are impeded due to the poor mechanical properties or complex structures of these systems. In recent years, the developments of metallic glasses (MGs) overcome this barrier with the abundant mechanical experiments on shear bands in laboratories, which greatly promotes the insights to shear bands. In addition, the dynamics of shear bands have a significant influence to the mechanical properties and behaviors of MGs, which is also of significant importance for understanding the deformation mechanisms of MGs. The insights for the structural origin, morphology, affect area, properties and motions of shear bands has made great progress. With the investigations on MGs, shear mands are found inhomogenous in space distribution and non-steady in time. The behaviors of shear bands show a similarity with that of the complex systems in nature as well as physic areas. This article mainly reviews our recent studies on the complex behaviors of shear bands in MGs, including the serration follow behaviors and the self-organizing criticality (SOC) behaviors of shear bands. The models that quantify these behaviors of shear bands, such as stick-slip model, machine-sample coupling model, are also introduced. The review also identified the key questions remaining to be answered, and presents an outlook for the field.

A GENERAL METHOD TO INTRODUCE PRE-CRACK IN BULK METALLIC GLASSES FOR PLANE STRAIN FRACTURE TOUGHNESS TEST 1)
Xie Yiling, Liu Ze
2020, 52(2):  392-399.  DOI: 10.6052/0459-1879-19-377
Asbtract ( 310 )   HTML ( 4)   PDF (19766KB) ( 74 )
Figures and Tables | References | Related Articles | Metrics

We report a simple, low-cost but robust method to introduce pre-crack in bulk metallic glasses (BMGS) for plane strain fracture toughness test. In recent years, bulk metallic glasses have shown potential and promising engineering application due to their excellent properties such as high elasticity, high strength, wear resistance and soft magnetism. As an important material parameter for engineering application, fracture toughness has also attracted great attention in the BMG community. However, there is still challenge on the fracture toughness test of BMGs because of the metastable nature and limited maximum castable size of BMGs. On the one hand, the casting induced thermal history difference and the defects such as internal micropores and impurity inclusions in BMGs, and the way of crack prefabrication will significantly affect the reliability of the fracture toughness test. On the other hand, limited by the sample size, most of the measured fracture toughness of BMGs are not the plane strain fracture toughness, resulting in the reported values showing large deviation, even for the same amorphous alloy. In this work, pre-notched BMG samples were thermoplastically compressed at the notched region to form ideal crack by creep flow. As an example, the fracture toughness of Zr-based amorphous alloy is tested by this method. The experimental results show that as the sample thickness increasing, the measured toughness decreases quickly and tends to a constant value. It should be pointed out that in the experiments, the local thermoplastic compression is designed to obtain a neck shaped pre-crack, which makes the minimum thickness of the specimen being far less than the required thickness by the plane strain fracture toughness testing standard.

STRUCTURE AND TOUGHNESS MODULATION OF A Zr$_{52.5}$Cu$_{17.9}$Ni$_{14.6}$Al$_{10}$Ti$_{5}$ METALLIC GLASS BY SURFACE MECHANICAL ATTRITION TREATMENT 1)
Chen Ken, Huang Bo, Wang Qing, Wang Gang
2020, 52(2):  400-407.  DOI: 10.6052/0459-1879-20-030
Asbtract ( 290 )   HTML ( 4)   PDF (17792KB) ( 97 )
Figures and Tables | References | Related Articles | Metrics

As a new type of structural material, the toughness of metallic glasses (MGs) needs to be further improved. The methods of improving the toughness of MGs include introducing dendrite phase, tuning their compositions to change the Poisson's ratios in order to affect the formation and spread of shear bands and cracks, {etc}. In this paper, we use the method of surface mechanical treatment to alter the microstructure and toughness of MGs. A Zr$_{52.5}$Cu$_{17.9}$Ni$_{14.6}$Al$_{10}$Ti$_{5}$ (at. %) MG (Vit105) plate was prepared by arc melting in vacuum and centrifugal casting system for thin plates in the metastable state. Surface mechanical attrition treatment (SMAT) is introduced to form nanoscale local crystal-like ordered structure in Vit105. Through differential scanning calorimetry and nano-indentation experiments, we find that the relaxation enthalpy near the surface of the SMAT-treated Vit105 plate is reduced, and its microstructure is more homogenous and stable The analysis by Vickers hardness tester shows that the hardness of the regions near the surface is increased and the hardness values distribute more narrowly after the SMAT treatment. Three-point bending fracture experiment reveals that notch toughness of the plate is also improved by SMAT. By SMAT treatment, the notch toughness increases from $70.7\pm 4.7$ MPa$\cdot$m$^{1/2}$ to $112.8\pm 3.7$ MPa$\cdot$m$^{1/2}$. Meanwhile, the density of shear bands becomes larger near the fracture surface as compared to the untreated sample. The enhancement of the toughness of Vit105 plate treated by SMAT originates from the promotion of the formation of shear bands. Our studies show that surface mechanical treatment leads to the formation of local crystal-like ordered structure in MGs with the enhancement of structural homogeneity. The hardness and toughness of MGs are improved, being associated with the formation of profusive shear bands. As a novel approach of improving the properties of materials, surface mechanical treatment has a broad application prospect in future.

Fluid Mechanics
NUMERICAL SIMULATION OF THE HYDRODYNAMIC CHARACTERISTICS OF VIOLENT AERATED FLOWS NEAR VERTICAL STRUCTURE UNDER WAVE ACTION 1)
Deng Bin,Wang Mengfei,Huang Zongwei,Wu Zhiyuan,Jiang Changbo
2020, 52(2):  408-419.  DOI: 10.6052/0459-1879-20-029
Asbtract ( 316 )   HTML ( 14)   PDF (27585KB) ( 268 )
Figures and Tables | References | Related Articles | Metrics

The air entrainment of water body under wave breaking is easy to produce pressure oscillation on the coastal structures. Therefore, in order to calculate the forces on the structures accurately, it is necessary to understand the motion characteristics of the aerated flows near coastal structure under wave action. In this study, a three-dimensional numerical wave flume is developed based on the OpenFOAM open source package together with the modified velocity inlet wave method. The model uses the S-A IDDES turbulence model to close the equations and the modified VOF method to capture the free surface. The interaction process of regular wave with a vertical structure on the smooth slope of 1:10 is simulated. Based on the comparison with the experimental results of the laboratory model, the hydrodynamic characteristics of violent aerated flows near the vertical structure are analyzed. The results show that the numerical model established in this study can accurately capture the space-time distribution of the free surface near the vertical structure and the transport process of bubbles under the action of regular waves. At the same time, the model can not only better describe the shape of the air cavity formed by the air entrainment, but also capture the process of merging and splitting among multiple cavities. The interaction between breaking waves and the vertical structure produces violent aerated flows, and its movement process is very complicated. The air-water interface is always accompanied by vortices in the wave propagation process. The bubble splitting in the air-water mixture is closely related to the positive and negative vortices. The transport trajectory of bubbles is mainly affected by the flow field around the bubbles. It is found that the distribution of turbulent kinetic energy in the vicinity of a vertical structure under the action of waves has a linear relationship with the characteristics of bubbles (bubble number, void fraction).

SELF-ASSEMBLY OF BUBBLE SWARM IN LARGE CAVITIES IN STEP-TYPE PARALLELIZED MICROCHANNELS AND ITS FEEDBACK ON BUBBLE FORMATION 1)
Zhang Zhiwei, Yin Xiangyu, Zhu Chunying, Ma Youguang, Fu Taotao
2020, 52(2):  420-430.  DOI: 10.6052/0459-1879-19-251
Asbtract ( 336 )   HTML ( 8)   PDF (26827KB) ( 59 )
Figures and Tables | References | Related Articles | Metrics

Step-emulsification microfluidic devices attract attentions due to the ability for the high-throughput production of bubbles and droplets. In this study, the bubble behavior of the bubble swarm in the cavity in a step-type parallelized microchannel and its feedback effect on bubble formation were studied by using a high-speed camera. The operating variables were the position of gas/liquid phase inlet, gas flow rate and liquid flow rate. Two bubble generation modes were observed: single-tube generation mode and multi-tube generation mode. The variation trend of the complex behavior of bubble swarm in the cavity with operating conditions was studied. It is found that in the confined space, the bubble can be self-assembled into a two-dimensional lattice with geometric features in the horizontal plane: including an ordered row triangular lattice, an ordered vertical triangular lattice, and a disordered triangular lattice. As the gas pressure changes, the crystal lattice changes, the bubble surface area becomes larger, and the interfacial energy becomes larger. Then, the effects of complex behavior of bubble swarm on bubble generation were analyzed by using the concepts of mesoscale, energy and activation. These fully explain the mesoscale characteristics of bubble behavior in confined spaces. The determinants of the bubble self-assembly path are revealed. It is indicated that the self-assembly path of bubbles is determined by the size and distribution of bubbles. The coefficient of variation of CV is used to represent the crystal structure of bubbles, and the coefficient of variation for the ordered triangular lattice is less than 5%, while that for the disordered triangular lattice is greater than 5%.

A THEORETICAL AND EXPERIMENTAL STUDY ON THE STATIC AERO-ELASTIC INSTABILITY OF AN INVERTED CANTILEVERED PLATE IN A CONFINED SUBSONIC FLOW 1)
Zhang Dechun, Li Peng, Liang Sen, Yang Yiren
2020, 52(2):  431-441.  DOI: 10.6052/0459-1879-19-255
Asbtract ( 285 )   HTML ( 5)   PDF (7383KB) ( 81 )
Figures and Tables | References | Related Articles | Metrics

The plate and shell structures have been widely used in many engineering areas such as aerospace, high-speed trains and energy harvesting. An inverted cantilever plate in an axial airflow confined by a rigid wall has the adjustable critical airflow velocity and is one of the most important ways for the optimization of energy harvesters. However, the mechanism of instability of inverted cantilever plates in rigid wall effect still needs further studies. This paper aims at the static aero-elastic instability of an inverted cantilevered plate in a confined subsonic flow, and such problem is studied by both the theoretical and the experiment analysis. In theory, first the effect of rigid wall is modelled by using the mirror image function and the fluid force on the plate is presented as a Possio integral equation within the frame work of differential operators, and the wall effect is featured by a composite operator involving the shifted Tricomi operator; then the solution of the instability equation is changed into the problem of approximation of function at a given interval; and finally the optimal solution is obtained by using the least square method with the help of the Wererstrass theorem. In experimental analysis, a test method for the static instability of the plate developed from the theory of column is applied for the wind tunnel test. The theoretical results show that the plate experiences the static aero-elastic instability; the critical dynamic pressure increases with the the increasing spacing between the wall and the plate but finally becomes steady (the case without rigid wall). A comparison of the present theoretical results with other existing theory and the experiment shows that the present theory is reliable and accurate. In this paper, the instability problem, which is not solved with the discretization of partial differential equation and by the eigenvalue analysis, has been transformed into the problem of approximation of function, and this can serve as another new thinking way to the problem of static aero-elastic instability of plates.

DEM ANALYSIS OF FACTORS INFLUENCING THE GRANULAR CAPILLARITY 1)
Zhang Huateng, Fan Fengxian, Wang Zhiqiang
2020, 52(2):  442-450.  DOI: 10.6052/0459-1879-19-301
Asbtract ( 271 )   HTML ( 5)   PDF (7664KB) ( 82 )
Figures and Tables | References | Related Articles | Metrics

Granular capillarity refers to the phenomenon that when a narrow tube is vertically inserted into a container filled with particles and then set into vertical vibration, the particles rise up along the tube and eventually reach a certain height. This provides a potential technical method for the transportation of granular materials against gravity. To explore the factors influencing the granular capillarity, the processes of granular capillarity were numerically investigated using the discrete element method (DEM). On this basis, the evolutions of vertical velocities of particles at different tube diameters were shown, and the dependences of the final capillary height of the particles on the tube diameter at different container widths and vibrational parameters were examined. The results obtained under the conditions with the container-width-to-particle-diameter ratio of 40, vibration-amplitude-to-particle-diameter ratio of 14.33 and vibration frequency of 12 Hz show that at the tube-to-particle diameter ratio $D/d = 3.33$ severe jamming occurs for the particles in the tube, which makes the particles rise slowly and leads to a discontinuous granular column within the tube. At $D/d = 8.33$, the granular capillary height rises rapidly at the beginning, and then the increasing rate of the capillary height decreases gradually. In this case, there is almost no particle velocity gradient along the tube radius. However, at $D/d=15$, as the granular column height within the tube increases, the particles in the tube separated into two layers. In the upper layer, almost no particle velocity gradient along the tube radius is observed, whereas in the lower layer obvious velocity gradient can be found. It is also found that for the tube diameter range in which the granular capillarity can occur, there exists a critical tube diameter that corresponds to the maximal final capillary height. When the tube diameter is less than the critical tube diameter, the final capillary height increases with the tube diameter, whereas when the tube size is greater than the critical tube diameter, the final capillary height tends to decrease with the tube diameter. Moreover, increasing the container width leads to an increase in the critical tube diameter, while increasing the vibration amplitude of the tube and appropriately increasing the vibration frequency can effectively promote the increase of the critical tube diameter.

Solid Mechanics
STABILITY AND RECOVERABILITY OF LIQUID-GAS INTERFACES ON SUBMERGED HIERARCHICALLY STRUCTURED SURFACES 1)
Yang Songmo, Wang Gang, Cao Yanlin, Huang Zhongyi, Duan Huiling, Lü Pengyu
2020, 52(2):  451-461.  DOI: 10.6052/0459-1879-20-025
Asbtract ( 341 )   HTML ( 8)   PDF (25479KB) ( 105 )
Figures and Tables | References | Related Articles | Metrics

The liquid-gas interface formed by immersing the surfaces with microstructures underwater is of great importance for applications such as drag reduction and so on. The stable existence of the liquid-gas interface is a prerequisite for the function of microstructure function surfaces. Therefore, how to enhance the stability of the liquid-gas interface to resist the wetting transition process, and how to implement the de-wetting process to improve the recoverability of the liquid-gas interface after collapse, both have scientific research significance and practical application value, and have triggered extensive investigations at home and abroad. This work is dedicated to investigating the stability and recoverability of liquid-gas interfaces formed on hierarchically structured surfaces after immersion in water. Different kinds of hierarchically structured surfaces were firstly designed and fabricated in order to investigate the influence of the sublevel structures respectively on the sidewalls and the bottom on the stability and recoverability of the liquid-gas interface. Experiments using laser scanning confocal microscopy to in-situ investigate the collapse and recovery process of liquid-gas interfaces were then performed. Theoretical analysis based on minimizing the total free energy of the system was further completed with the aim to better understand the inner mechanism. The experimental results agreed well with the theoretical analysis. This work reveals the mechanism of hierarchically structured surfaces resisting wetting transition and improving liquid-gas interfaces recoverability: sublevel structures (nanoparticles, fins) on the sidewalls enhance the stability of the liquid-gas interface by increasing the apparent advancing contact angle; sublevel structures (nanoparticles and "closed'' sublevel structures) on the bottom surface are able to maintain the existence of nanoscale gas pockets, which is conducive to the diffusion of dissolved gas in the bulk water into the microstructure, and eventually helps the liquid-gas interface to recover. The research in this paper provides ideas for designing hierarchically structured surfaces to obtain liquid-gas interfaces with good stability and recoverability.

TANGENTIAL RESPONSE MODELING OF JOINT SURFACE BASED ON IWAN MODEL 1)
Zhan Wanglong, Li Wei, Huang Ping
2020, 52(2):  462-471.  DOI: 10.6052/0459-1879-19-343
Asbtract ( 280 )   HTML ( 4)   PDF (4609KB) ( 68 )
Figures and Tables | References | Related Articles | Metrics

Aiming at the problems of tangential displacement response of common used pre-tightening lap joints under tangential displacement excitation in engineering, a new slip force density distribution function based on actual morphological parameters and material performance parameters was proposed. The tangential response constitutive model of lap joint was obtained by using the distribution function, and the hysteresis curve and energy dissipation value per loading cycle were obtained. Compared with the published experimental results, it is found that the simulation value was consistent with the experimental, which proves the rationality of the model. Based on this, the relationship between the tangential displacement of the joint surface and the tangential force, tangential contact stiffness, and energy dissipation was studied. The results showed that the model could well describe the relationship between tangential force and tangential displacement. The trend of critical slip force function rise rapidly and then converges to zero after reaching its maximum value. The tangential force and tangential displacement show a non-linear characteristic, and with the increase of tangential displacement, the tangential contact stiffness shows softening phenomenon; the initial tangential stiffness is related to normal load, roughness parameters and plastic index. For a given contact surface, the greater normal force is given, the greater initial tangential stiffness is observed; the initial tangential stiffness also increases with the increase of plastic index.

A NEW METHOD FOR SOLVING THE GRADIENT BOUNDARY INTEGRAL EQUATION FOR THREE DIMENSIONAL POTENTIAL PROBLEMS 1)
Dong Rongrong, Zhang Chao, Zhang Yaoming
2020, 52(2):  472-479.  DOI: 10.6052/0459-1879-19-308
Asbtract ( 337 )   HTML ( 7)   PDF (3402KB) ( 57 )
Figures and Tables | References | Related Articles | Metrics

In the boundary element analysis of three-dimensional potential problems, it is a very difficult task to calculate the boundary potential gradients with respect to the space coordinates instead of normal one. Several techniques have been proposed to address this problem so far. They, however, usually need complex and lengthy theoretical deduction as well as a large number of numerical manipulation. In this study, a new method, named auxiliary boundary value problem method (ABVPM), is presented for solving the gradient boundary integral equation (GBIE) for three dimensional potential problems. An ABVPM with the same solution domain as the original boundary value problem is constructed, which is an over-determined boundary value problem with known solution. Consequently, the system matrix of the GBIE, which is the most important problem for boundary analysis, will be obtained by solving this ABVPM. It can be used to solve original boundary value problem. The solution procedure is very simple, because only a linear system need to be solved to obtain the solution of the original boundary value problem. It is worth noting that when solving the original boundary value problem, it is not necessary to recalculate the system matrix, so the efficiency of the auxiliary boundary value method is not very poor. The proposed ABVPM circumvents the troublesome issue of computing the strongly singular integrals, with some advantages, such as simple mathematical deduction, easy programming and high accuracy. More importantly, the ABVPM provides a new idea and way for solving the GBIE. Three benchmark examples are tested to verify the effectiveness of the proposed scheme.

AN ABSORBING BOUNDARY CONDITION FOR SH WAVE PROPAGATION IN VISCOELASTIC MULTILAYERED MEDIA 1)
Wu Lihua, Zhao Mi, Du Xiuli
2020, 52(2):  480-490.  DOI: 10.6052/0459-1879-19-315
Asbtract ( 302 )   HTML ( 2)   PDF (1580KB) ( 91 )
Figures and Tables | References | Related Articles | Metrics

A high-accuracy absorbing boundary condition in the time domain is proposed, which can be coupled with the finite element method seamlessly to simulate the propagation of the transient scalar wave in the D'Alembert viscoelastic multilayered media. First, a semi-discrete displacement equation of the semi-infinite domain and the force-displacement relationship in the artificial boundary are obtained by semi-discretizing the semi-infinite domain along the vertical depth. Modal decomposition is utilized to convert the field of the semi-infinite domain in the physical space into the modal space. Then the dynamic stiffness of the semi-infinite domain in the frequency domain in the modal space is obtained according to both the displacement equation and the force-displacement relationship in the modal space. Second, a scalar continued fraction, which is convergent over the whole frequency domain, is proposed to describe the scalar dynamic stiffness in the modal space of the D'Alembert viscoelastic single-layered medium. The scalar continued fraction is extended to the matrix form to represent the dynamic stiffness in the modal space of the D'Alembert viscoelastic multilayered media. Last, by introducing auxiliary variables, a time-domain absorbing boundary condition in the modal space is constructed based on the proposed continued fraction. Subsequently, considering the relationship of the field in the modal space and in the physical space, a time-domain absorbing boundary condition in the physical space is obtained by converting the absorbing boundary condition in the modal space into the physical space. Two numerical examples of a single-layered medium and a five-layered media verify that the proposed method is accurate and stable for the D'Alembert viscoelastic single-layered medium, and for the D'Alembert viscoelastic multilayered media, in order to ensure the proposed method's property of high-accuracy, the distance from artificial boundary to the region of interest needs to be about 0.5 times of the total layer height of the infinite domain.

INVESTIGATIONS ON THE THERMAL BEHAVIOR AND ASSOCIATED THERMAL STRESSES OF THE FRACTIONAL HEAT CONDUCTION FOR SHORT PULSE LASER HEATING 1)
Xu Guangyin, Wang Jinbao, Xue Dawen
2020, 52(2):  491-502.  DOI: 10.6052/0459-1879-19-331
Asbtract ( 262 )   HTML ( 3)   PDF (4481KB) ( 53 )
Figures and Tables | References | Related Articles | Metrics

For complex heat transfer process and related thermal stress in materials subjected to short pulse laser heating, the existing thermal stress theory based on Fourier law or Cattaneo-Vernotte relaxation equation combined with elastic theory has serious defects in describing its thermophysical process.In this paper, based on the fractional calculus theory, the fractional Cattaneo type heat conduction equation and the corresponding thermal stress equation with appropriate initial and boundary conditions are established for a semi-infinite space irradiated by non-Gaussian lase. The analytical solutions of the temperature field and the thermal stress field are obtained via Laplace transform method, and the thermophysical behaviors are illustrated. Firstly, the theoretical solution is verified, then the variations of temperature field and thermal stress field are studied under the fractional order $p=0.5$, and the influence of laser parameters on temperature and thermal stress field are also researched. Finally, the effects of fractional order parameters on temperature and thermal stress field are calculated. The calculation results show that the temperature and thermal stress fields described by the fractional Cattaneo type heat transfer equation and thermal stress equation have wave diffusion the characteristics. Compared with the classical Fourier heat transfer model and the standard Cattaneo type heat transfer model, the larger the fractional order is, the smaller the thermal wave velocity is, the more significant the thermal wave dynamics is. On the contrary, the larger the thermal wave velocity is, the stronger the thermal diffusivity is. The faster the laser heating and cooling rate is, the faster the temperature rises and falls, the faster the alternating change of compressive stress and tensile stress is, the smaller the temperature change amplitude is, and the variations of thermal stress amplitude is not obvious

Dynamics, Vibration and Control
RESEARCH ON SEPARATION STRATEGY AND DEPLOYMENT DYNAMICS OF A SPACE MULTI-RIGID-BODY SYSTEM 1)
Luo Caoqun, Sun Jialiang, Wen Hao, Hu Haiyan, Jin Dongping
2020, 52(2):  503-513.  DOI: 10.6052/0459-1879-19-307
Asbtract ( 357 )   HTML ( 7)   PDF (4862KB) ( 149 )
Figures and Tables | References | Related Articles | Metrics

The paper focuses on the separation and deployment dynamics of an on-orbit compactly connected multi-rigid-body (MRB) system, which could separate autonomously from a carrier spacecraft. Based on the focused MRB system, it is not necessary to repeatedly use the launcher of the carrier spacecraft or install multiple launchers in the spacecraft to separate the MRB system. This is advantageous because it can effectively improve the space utilization rate of the spacecraft, simplify the separation deployment operations and reduce the risk of collision between rigid bodies. To realize the separation of such a MRB system, the paper presents an investigation on its on-orbit dynamics and the design of collision-free separation deployment schemes. Firstly, a dynamic model of a single rigid body is established based on the principle of virtual work and the Natural Coordinate Formulation (NCF) method accounting for the relative motion between rigid bodies and attitude changes of each rigid body. Considering the orbital motion, the variations of connecting constraints of the MRB system and the interactions between rigid bodies during the separation, the governing nonlinear dynamic equations including constraints of the system are obtained with a method of Lagrange multipliers. With practical engineering applications taken into consideration, the separation deployment of MRB system is realized through ejection mechanisms mounted on the four corners of each contact surface between rigid bodies. Secondly, the timing sequences of separation maneuvers are specially programmed and two separation schemes are developed by adjusting different ejection directions and ejection sequences to guarantee the non-collision between rigid bodies in the separation deployment. Finally, numerical case studies are presented for investigating the nonlinear dynamic behaviors of rigid bodies and demonstrating the effectiveness of separation schemes.

PRIMARY AND SUBHARMONIC SIMULTANEOUS RESONANCE OF DUFFING OSCILLATOR 1)
Li Hang, Shen Yongjun, Li Xianghong, Han Yanjun, Peng Mengfei
2020, 52(2):  514-521.  DOI: 10.6052/0459-1879-19-349
Asbtract ( 370 )   HTML ( 5)   PDF (7927KB) ( 135 )
Figures and Tables | References | Related Articles | Metrics

In this paper, the dynamics and stability of the Duffing oscillator subjected to the primary resonance together with the 1/3 subharmonic resonance are studied. At first, the approximate analytical solution and amplitude-frequency equation are obtained through the method of multiple scales, and the correctness and satisfactory precision of the approximate solution are verified by simulation. Then, the amplitude-frequency equation and phase-frequency equation of steady-state response are derived from the approximate analytical solution, and it can be found there are at most seven different periodic solutions, which are called multi-value characteristics and can be used to switch the state of the system. Moreover, the stability condition of steady-state response is derived based on Lyapunov theory, and the amplitude-frequency curves of steady-state response are compared with the cases where the primary or 1/3 subharmonic resonance exists alone, and it is found that the system contains both resonance characteristics. At last, the effects of nonlinear factor and excitations on the system response are analyzed by simulation. The particular phenomena in this system are revealed, i.e., the nonlinear factor affects the response amplitude, multi-value characteristics and stability of the system with stiffness softening. However, for the stiffness hardening system, the nonlinear factor only affects the response amplitude, which is similar to the cases of single-frequency excitation. These results are important for the study on the Duffing system or other similar systems.

COMPARISON OF VIBRATION CHARACTERISTICS OF THREE TYPICAL AXIALLY MOVING STRUCTURES 1)
Liu Xingguang, Tang Youqi, Zhou Yuan
2020, 52(2):  522-532.  DOI: 10.6052/0459-1879-19-304
Asbtract ( 361 )   HTML ( 7)   PDF (18727KB) ( 86 )
Figures and Tables | References | Related Articles | Metrics

The transverse vibration of the axially moving structure is always one of the hot topics in the field of dynamics. At present, most literatures are considering the study of one model. There are few researches considering the comparative analysis of several models. The vibration characteristics of three typical axially moving structures, such as the Euler beam, the panel, and the plate with two opposite sides simply supported and other two free, are compared and analyzed in this paper. In view of different structural parameters in engineering, this paper provides a reference for choosing a more reasonable model in the study of vibration theory. The governing equations of the three models are solved by the complex mode method. The corresponding natural frequencies and mode functions can be obtained. For the plate model, two rigid displacements and the coupled flexural and torsional vibration are both considered. The variations of the first four order natural frequencies of the three models with the axial velocity and the aspect ratio are given by numerical examples. At the same time, the analytical solution obtained by the complex modal method is verified by the differential quadrature method. The influence of different axial speed, damping, stiffness, and aspect ratio on the first order natural frequency of three models is analyzed by adopting the form of three-dimensional diagram for the first time. The effects of different aspect ratios and axial speed mixing on the relative errors of the first natural frequencies the first order natural frequency of the plates and beams are emphatically studied. The results show that the natural frequency of the three structures decreases gradually with the increasing axially speed. The panel is a simplified model of the plate. The damping has the least effect on the first order natural frequency when other parameters change. The complex model can be simplified into a simple model when the moving structure has a large aspect ratio and a small axial speed.

STUDY ON COUPLED VIBRATIONS OF ROTOR WITH AN ARBITRARY SLANT CRACK 1)
Jiao Weidong, Jiang Yonghua, Li Gang, Cai Jiancheng
2020, 52(2):  533-543.  DOI: 10.6052/0459-1879-19-292
Asbtract ( 267 )   HTML ( 5)   PDF (5953KB) ( 66 )
Figures and Tables | References | Related Articles | Metrics

The coupled vibrations of rotor with different type of cracks were studied, including transverse crack, transverse-slant crack and arbitrary slant crack. It aims to reveal the variation rule of stiffness parameters in different directions, their cross-coupling mechanism and especially the feature of resulted vibrations of the cracked rotor. The rotor segment with different type of cracks, including transverse crack, transverse-slant crack and arbitrary slant crack, was respectively modeled using the Timoshenko beam element, considering all six degrees of freedom including longitudinal, bending and torsional directions. Flexibility coefficients and then stiffness matrix was derived, based on strain energy theory. On this basis, motion equation of the cracked rotor was solved by numerical algorithm `Newmark-beta' to obtain the time-domain response of coupled vibrations of the cracked rotor under different excitations by single or multiple faults such as unbalanced excitation, torsional excitation and combined excitation of unbalanced with torsional. The spectral characteristics of coupled vibrations were then analyzed. Compared to either transverse crack or transverse-slant crack, the stiffness matrix governed by arbitrary slant crack is more populated with additional cross coupled coefficients, and its cross coupling effect is more significant, which contributes to stronger coupled vibrations in bending-torsional direction or even longitudinal-bending-torsional direction. The amplitudes of both bending vibration and torsional vibration are larger than that at the situation of either transverse crack or transverse-slant crack, under the action of either unbalanced excitation or torsional excitation. Moreover, the characteristic frequencies of coupled vibrations of rotor with different type of cracks including transverse crack, transverse-slant crack and arbitrary slant crack, are of different sensitivities to oriented angle of crack surface, for example the rotating fundamental frequency and its second harmonic frequency, the torsional excitation frequency and its sideband components. These research results can provide theoretical basis for both characteristic parameters identification and diagnosis of rotor cracks.

Biomechanics, Engineering and Interdiscipliary Mechanics
MECHANICAL STABILITY ANALYSIS OF STRATA AND WELLBORE ASSOCIATED WITH GAS PRODUCTION FROM OCEANIC HYDRATE-BEARING SEDIMENTS BY DEPRESSURIZATION 1)
Yuan Yilong,Xu Tianfu,Xin Xin,Xia Yingli,Li Bing
2020, 52(2):  544-555.  DOI: 10.6052/0459-1879-19-178
Asbtract ( 305 )   HTML ( 8)   PDF (16269KB) ( 92 )
Figures and Tables | References | Related Articles | Metrics

Gas production from hydrate-bearing sediment by depressurization will result in the changes of its porosity, permeability, pore pressure, effective stress and cementing strength, which will reduce the shear strength and carrying capacity of the hydrate-bearing sediments. As a consequence, the possible geohazards could trigger, such as wellbore instability, submarine landslide, and seafloor subsidence. Given this, based on the extended 3D Biot consolidation theory, the coupled thermo-hydro- mechanical (THM) model was built in the framework of TOUGH+Hydrate. The new THM model considers the coupling processes, including hydrate dissociation and formation, heat conduction, convection, and mechanical behavior of hydrate-bearing sediments. Taking the geologic data obtained from the Shenhu area SH2 site of South China Sea as an example, the THM coupling model for analyzing the mechanical stability of strata and wellbore using vertical well by depressurization was constructed. The evolution rules of reservoir temperature, pore pressure, stress, and hydrate dissociation zone during depressurization were predicted. In addition, the dominant sand producing region and seabed subsidence trend were revealed. The results show that the drawdown of reservoir pressure controls the increase in effective stress and strata subsidence around the production well. The subsidence mainly occurs in the early stage of depressurization, and the maximum subsidence is located around the production interval. The dissociation of hydrate results in significant decrease in the mechanical strength of sediments, which aggravates the subsidence. Wellbore depressurization results in the stress concentration in the perforation section significantly, which causes the potential wellbore damage. What's more, the stress concentration region is the key area for the prevention and control of sand production associated with hydrate production.

AN ELASTOPLASTIC CONSTITUTIVE MODEL FOR GAS HYDRATE-BEARING SEDIMENTS 1)
Liu Lin, Yao Yangping, Zhang Xuhui, Lu Xiaobing, Wang Shuyun
2020, 52(2):  556-566.  DOI: 10.6052/0459-1879-19-184
Asbtract ( 215 )   HTML ( 7)   PDF (10509KB) ( 86 )
Figures and Tables | References | Related Articles | Metrics

The density is one of the most important factors for the mechanical behavior of soil. The content and occurrence modes of hydrates obviously affect the density of gas hydrate-bearing sediments (GHBS) because hydrates exist in the pore of sediment as solid phase. Therefore, it is necessary to consider the effect of the hydrate content and occurrences to describe the mechanical properties of hydrate sediments well. In this paper,based on the unified hardening model for clays and sands (CSUH model), the relation between the volume fraction of hydrate and the compressive hardness parameter is firstly established to reflect the influence of hydrate on the compressibility of sediments. Secondly, in order to consider the influence of hydrate content and occurrence modes on sediment density, we propose a formula to describe the effective initial void ratio, and it is then introduced into the state parameter to describe the influence of hydrate on dilatancy of sediment. Finally, combining the drop-shaped yield surface of the CSUH model, an elastoplastic constitutive model for GHBS is developed. Compared with the laboratory test results, it is verified that the model can reasonably describe the mechanical behaviors of GHBS containing hydrates with different occurrence modes and contents. For the same occurrence mode but different contents, the set of parameters is the same.

SUBMARINE SLOPE STABILITY DURING DEPRESSURIZATION AND THERMAL STIMULATION HYDRATE PRODUCTION WITH HORIZONTAL WELLS 1)
Tan Lin,Liu Fang
2020, 52(2):  567-577.  DOI: 10.6052/0459-1879-19-171
Asbtract ( 310 )   HTML ( 3)   PDF (4600KB) ( 53 )
Figures and Tables | References | Related Articles | Metrics

Gas hydrate is extensively found in the voids of sediments under deep sea. It is considered one of the promising sources of future clean energy, and has attracted global attentions of research and development. Submarine slope instability could be triggered during gas production from oceanic hydrate reservoirs due to reduction in soil strength and build-up of pore pressure. This paper investigates the impact of hydrate extraction on submarine slope stability within the framework of the limit equilibrium slope stability analysis considering the thermo-hydro-chemical (THC) coupled processes of hydrate dissociation. A THC coupled simulator TOUGH + HYDRATE was employed to simulate the process of hydrate production under depressurization and thermal stimulation with horizontal wells. The simulation captured the propagation of hydrate dissociation front and the evolution of pore pressure during hydrate production. The safety factor of the slope stability during and after production was then acquired from the limit equilibrium slope analysis by means of SLOPE/W. A parametric study was conducted to investigate the effect of well locations and production methods on the slope stability. The results indicate that, for a steep slope with a sandy reservoir covered by a tight overburden layer, the slope is stabilized during single-well depressurization production due to increase in the effective stress caused by depressurization. Compared to wells installed elsewhere, a depressurization well installed at the mid-height of the slope provides the highest safety factor during production. As the pore pressure recovers back to hydrostatic conditions and the cohesion of the host soils decreases, the slope becomes less stable after production. For scenarios using double-well thermal stimulation production, the critical slip surface passes through dissociated zones, and the safety factor of the slope significantly drops due to high pressure build-up developed around the injection well. Submarine landslides could be triggered by hydrate production with thermal stimulation.

ANALYSIS OF WATER-SEABED-STRUCTURE DYNAMIC INTERACTION EXCITED BY PLANE WAVES 1)
Chen Shaolin, Sun Jie, Ke Xiaofei
2020, 52(2):  578-590.  DOI: 10.6052/0459-1879-19-354
Asbtract ( 262 )   HTML ( 1)   PDF (18558KB) ( 49 )
Figures and Tables | References | Related Articles | Metrics

The seismic response analysis of marine engineering structures is important to ensure the seismic safety of these structures. Marine engineering structures are always built in a complex environment where there are both water and soil. So both fluid-structure interaction and soil-structure interaction need to be considered in seismic response analysis of marine engineering structures. This paper is based on the unified calculation framework for wave motion in water-saturated seabed-bedrock system, and uses the Davidenkov model and modified Masing rule to describe the nonlinear characteristics of saturated seabed. The dynamic responses of marine sites and marine engineering structures when the pulse SV wave incident vertically are analyzed. First, the nonlinear responses of marine sites under the input of linear free field and nonlinear free field are analyzed comparatively. The results show that the nonlinear response of the marine site under the input of the linear free field is unreasonable. So it means the consistent constitutive model should be used for the free field analysis and finite element analysis of marine site. Then, the characteristics of responses of the marine site and water-saturated seabed-bedrock system when the seabed is linear or nonlinear are analyzed comparatively. Compared with the case that the seabed is linear, the influence of nonlinearity on response of seabed is mainly controlled by the following two factors: on the one hand, the modulus of saturated seabed decreases because of nonlinearity, so the wave impedance ratio between saturated seabed and bedrock decreases, and then, the reflection coefficient and transmission coefficient from bedrock to saturated seabed increase, resulting in the increase of response; on the other hand, the nonlinearity of saturated seabed causes the increase of damping, resulting in the decrease of response. For the example in this paper, the effect of damping dominate the result when seabed is nonlinear.

WAVE FUNCTIONS METHOD FOR CALCULATING THE DYNAMIC RESPONSE OF A TUNNEL IN UNSATURATED SOIL 1)
Guo Huiji, Di Honggui, Zhou Shunhua, He Chao, Zhang Xiaohui
2020, 52(2):  591-602.  DOI: 10.6052/0459-1879-19-280
Asbtract ( 266 )   HTML ( 3)   PDF (777KB) ( 96 )
Figures and Tables | References | Related Articles | Metrics

This paper presents a wave functions method to calculate the three dimensional dynamic interaction between a circular tunnel and the surrounding unsaturated foundation soil. The tunnel lining is conceptualized as an infinite Flügge cylindrical shell, and the unsaturated foundation soil is conceptualized as a three-phase medium consisting of fluid, solid and gas. The motion equations of infinite Flügge cylindrical shell are solved by the separation of variable method, and the Helmholtz decomposition theorem is used to solve the unsaturated soil governing equations. Based on the boundary conditions of displacement, stress and pore fluid pressure at the unsaturated soil-tunnel interface and the ground surface, the dynamic interaction between the unsaturated soil and circular tunnel is solved when the tunnel invert subjected to a unit harmonic load. The transformation properties of plane wave functions and cylindrical wave functions are used to apply the boundary conditions expressed in both the rectangular and cylindrical coordinate systems before the coupling. The proposed method is verified by comparing with the results obtained by the existing 2.5 dimensional coupled FE-BE method for single phase elastic medium, 2.5 dimensional coupled FE-BE method for two phase saturated porous medium and Pip in Pip semi analytical method for three-phase unsaturated medium. Finally, by using the proposed method, the dynamic response of the unsaturated soil-tunnel system under different soil saturations is analyzed by a case study. The results show that saturation has a great effect on the dynamic response amplitude of soil displacement and excess pore water pressure. Because the three-phase medium simulating unsaturated soil can be reduced to two-phase medium or single-phase medium, this method can also be used to calculate and analyze the dynamic response of a circular tunnel buried in two-phase saturated soil or in single-phase elastic soil after the degeneration of the calculation parameters of unsaturated soil.