EI、Scopus 收录

## 2021 Vol. 53, No. 5

Display Method:
The orbital dynamics in the three-body system is a classical problem in the field of astrodynamics. It has rich theoretical and engineering significance, and have played an important role in the process of space activities extending from near-earth space to deep space. This paper reviews and summarizes the progress of the study of orbital dynamics in the three-body system. Combined with the development trend of deep space exploration in the future, the hotspots and challenges in the research of three-body orbital dynamics are prospected. First, the research background and significance of the three-body problem are surveyed, and the development of the dynamics model of the three-body system is briefly reviewed. Secondly, characteristics of the local motion near the equilibrium point in the three-body problem are systematically summarized. The analytical and numerical methods for periodic orbits are introduced. The latest development of quasi-periodic motion is discussed. Meanwhile, the characteristics and research progress of global periodic motions in the three-body system including resonance orbits, cycler trajectories, and free return orbits are summarized. Next, the research progress of the low-energy transfer and capture trajectory design in the three-body system is analyzed from two aspects of invariant manifold theory and weak stability boundary theory. Finally, the applications of orbital dynamics in the three-body system in formation flight and navigation constellation design are summarized. Several orbital dynamics and control problems in the design of landing trajectories for full lunar-surface coverage, the low thrust orbit optimization of the three-body system, and the utilization of non-linear equilibrium points in the three-body system are addressed for future study.
2021, 53(5): 1223-1245. doi: 10.6052/0459-1879-20-367
The interaction of shock wave with cylindrical interface is fundamental in study of the Richtmyer-Meshkov (RM) instability. Although the RM instability of two-dimensional (2D) cylindrical interfaces under a single shock wave has been extensively studied previously, the interaction of reflected shock (short for reshock) with cylindrical interfaces, especially three-dimensional (3D) cylindrical interfaces has not been investigated thoroughly, with relevant development rules and underlying mechanisms unclear. When the shock wave interacts with the evolving interface after reshock, new baroclinic vorticity appears on the interface and this will have a major influence on the evolution of the interface. In this work, the HOWD (high order WENO and double-flux) solver developed in our group is used to numerically study the reshock effect on the evolution of 2D and 3D concave cylindrical N$_{2}$/SF$_{6}$ (inner/outer phases) interfaces with incident planar shock strength of $Ma=1.29$. This work will focus on the evolution of 2D and 3D concave cylindrical interfaces after the reshock under different reflected distances, which is defined as the distance between the end wall and the center of the gas cylinder. Series of data have been extracted both before and after the reshock, including the schlieren and vorticity images of the evolving gas cylinder and the quantitative data of the geometric position of the feature points on the gas cylinder. The geometrical characteristics of the distorted interface and the generation and distribution of baroclinic vorticity in different stages are analyzed. The results indicate that for different reflected distances, the shapes of the evolving interface and the reshock at the interaction instance affect the generation and distribution of the baroclinic vorticity, resulting in distinct evolution characteristics of the RM instability. For the 3D concave cylindrical interface, the baroclinic vorticity distributed in 3D space at different heights can induce complicated 3D structures of the evolving interface.
2021, 53(5): 1246-1256. doi: 10.6052/0459-1879-21-042
For the reason that the scale-to-scale energy transfer characteristics of wall bounded turbulence exhibit strong anisotropy in three dimensional physical space, understanding the spatial distribution of energy scale-to-scale transfer is of great importance for the further construction of high-fidelity anisotropic large eddy simulation models. Direct numerical simulation were performed at Reynolds number $Re_{\tau } =180$ by using lattice Boltzmann method (LBM). In comparison with the open source channel flow turbulence database, both the mean velocity and turbulent variables, such as Reynolds stress and velocity fluctuations, agree well with the results of Kim et al. (1987) and Moser et al. (2015). The ability and validity of LBM on simulating turbulent channel flow were verified. The filter-space technique is used to obtain the scale-to-scale transport of kinetic energy, and the statistical results show that the backward scatter and forward scatter contributions to the SGS dissipation were comparable. Combined with the novel structure identification method, the distribution and geometry properties of the energy flux structure is quantitatively investigated by the clustering methodology. The results show that small scale structures account for the majority in number, and the volume probability density distribution presents a -4/5 power law. The structures can further be divided into wall attached structures and wall detached structures according to their minimal distance to the wall. In spite of the relatively small number fraction, the attached structure accounts for most of volume fraction, which shows that most of the attached structures are large scale structures. Further statistics of the attached structures exhibit that the size of the structure has a certain power law relationship, indicating that the scale-to-scale energy transport structure also has the self similarity of the attached eddy hypothesis which proposed by Townsend. The forward and backward energy flux structure pairs tend to be arranged side by side along the span direction.
2021, 53(5): 1257-1267. doi: 10.6052/0459-1879-20-432
In the present paper, the tip-leakage vortex (TLV) cavitation around a NACA0009 hydrofoil is numerically investigated with large eddy simulation method combined with a new Euler-Lagrangian cavitation model, which takes into account the effect of nuclei on tip-leakage cavitating flow. A satisfying agreement between the numerical and experimental results is obtained. Based on the numerical data, the evolution behavior of TLV cavitation with different gap sizes and its influence on the strength and radius of TLV, the nuclei distribution in the vortex core and the distribution of tangential velocity are discussed in detail. The results show that once cavitation occurs, the intensity of TLV is mainly affected by the evolution behavior of the sheet cavitation, while the tip-leakage cavitating flow has little effect on the strength of TLV. In addition, the results suggest that a smaller gap size will result in a more unstable sheet cavity, and the strength of TLV in turn presents corresponding quasi-periodic fluctuation as influenced by the evolution of the sheet cavitation. As the tip clearance progressively increases, the strength of the sheet cavitation gradually decreases with the reduction of instability, and the intensity of TLV gradually returns to the level without cavitation, together with its fluctuation level. Cavitation has a more significant influence on the nuclei distribution in the vortex core, and its degree of influence depends on the spatial stability of TLV and the intensity of TLV cavitation after the occurrence of cavitation. Moreover, the numerical results obtained by our simulations also indicate that cavitation has a significant influence on the radius of TLV. After cavitation occurs, the radius of TLV increases to a certain extent, and a rigid body rotation'' tangential velocity distribution is formed on the periphery of cavitation region, which is mainly caused by the expansion process induced by cavitation growth and the viscous effect of flow.
2021, 53(5): 1268-1287. doi: 10.6052/0459-1879-20-415
Both the experiment method and numerical simulation method are applied in this paper to investigate the energy transformation mechanism of a gas bubble collapseing in the free field. The bubble radius, velocity and acceleration of the bubble evolution process are obtained according to the experimental results via the schlieren method. These parameters are substituted~into the bubble potential energy and kinetic energy equation to explain the energy changing. By using the CFD simulation method, a three-dimensional model with reformulated mass conservation equation and momentum equation considering the weakly compressibility, is introduced to discuss the bubble collapse process. The pressure and velocity distribution around the bubble are extracted from the simulation results in order to analyze the energy transformation mechanism. The results show that (i) the relation between the potential energy and bubble radius maintains the positive correlation, with the increasing of the potential energy, the kinetic energy decreases significantly. The value of potential energy is maximum at the end of the bubble expanding, and the kinetic energy of free field returns to zero at the same time. (ii) During the shrinking stage, a high-pressure area appears around the bubble, and gradients of velocity and pressure are higher than the other area in the free field. The high-pressure area is shrinking gradually when the bubble is collapsing. (iii) The potential energy transforms into kinetic energy during the whole process of the bubble evolution, and the kinetic energy is form of wave energy. A shock wave is captured when the bubble collapse, and the most of bubble energy transform into wave energy of the shock wave at this time.
2021, 53(5): 1288-1301. doi: 10.6052/0459-1879-21-006
As inspired by the flexible characteristics of the feathers of birds, the influence of an anisotropic compliant wall on the spatial evolution of the T-S wave in the subsonic boundary layer flow is studied numerically in this article. First of all, the numerical results are in good agreement with that predicted by the linear stability theory and also the numerical results obtained by others. This demonstrates the reliability of the cell-centered finite difference method with high-order accuracy, which is adopted in this article. After that, a part of the rigid wall is replaced with an anisotropic compliant wall. The obtained results show that the compliant wall is able to reduce or even remove the unstable growth region of the T-S wave, i.e., suppress the amplitude growth of the T-S wave. Thus, the compliant wall has the potential to delay laminar-turbulent transition. The deformation of the compliant wall not only has components corresponding to the T-S wave waveform, but also the wall vibration with longer wavelength and same frequency compared to the T-S wave, which is caused by the leading edge of the compliant wall. The parameter study shows that with increased damping, the wall vibration by the leading edge is weakened. An increase in the rigidity, tension or elastic coefficient of the compliant wall will result in the wall stiffer and thus lead to smaller deformation. The greater the incline angle of the supporting lever arm, the stiffer the compliant wall. The increase of any of the above parameters will lead to a reduction in the suppression effect of the compliant wall on the T-S wave. When the flow direction is reversed, the suppression effect is also weakened. This research is aimed to reveal some of the mysteries of bird's efficient flying, and provides new ideas for passive flow control to reduce drag.
2021, 53(5): 1302-1312. doi: 10.6052/0459-1879-20-460
Thermocapillary migration of a droplet placed on a non-uniformly heated solid surface appears in a variety of practical applications, such as microfluidics, inkjet printing, et al. The flow stability analysis is crucial for the precise control of droplet migration. In the present work, the convective instability in thermocapillary migration of a wall-attached viscoelastic droplet is examined by linear stability analysis. The relation between the critical Marangoni number ($Ma_{\rm c}$) and the elastic number is obtained at different Prandtl numbers ($Pr$). The flow fields and energy mechanisms of preferred modes are analyzed. The results show that more kinds of preferred modes are excited by the elasticity. The preferred modes at small $Pr$ are the oblique and streamwise waves, while those at moderate and high $Pr$ are oblique waves and spanwise stationary modes. The strong elasticity significantly reduces the $Ma_{\rm c}$, while the weak elasticity slightly enhances the flow stability. $Ma_{\rm c}$ increases with $Pr$ at moderate $Pr$. For the oblique wave, the amplitude of perturbation temperature may appear in the middle region of flow field, while the amplitudes of other two modes only exist on the free surface. The distribution of streamlines is almost symmetric at high $Pr$. The energy analysis shows that the work done by the basic flow changes from positive to negative when the elastic number increases. The work done by the perturbation stress may either dissipate or provide energy at small $Pr$, while the work done by the basic flow is negligible at high $Pr$. For the downstream streamwise wave, the perturbation velocity and the work done by perturbation stress fluctuate several times in the vertical direction. Comparing the droplet migration with thermocapillary liquid layers, it can be found that due to the differences of basic flow and boundary conditions, there are quite different between their preferred modes and perturbation flow fields.
2021, 53(5): 1313-1323. doi: 10.6052/0459-1879-20-443
In this paper, the effect of sinusoidal alternating electric field on evaporation characteristics of three-dimensional suspended water droplet in supercritical nitrogen was studied by molecular dynamics method. And the effects of amplitude and frequency of alternating electric field on evaporation life and instantaneous evaporation rate of single droplet were mainly considered. The evaporation of droplet with 8000 water molecules in the environment of 27 000 nitrogen molecules was simulated and analyzed. Firstly, the physical parameters of water such as thermal conductivity at different temperatures and pressures were simulated using molecular dynamics method, and compared with the theoretical or experimental data. Then, the effect of uniform electric field on single droplet evaporation under subcritical conditions was studied, which was in good agreement with the results in literature, and verified the correctness of the molecular model and evaporation model. Lastly, the evaporation of single water droplet under the action of alternating electric field with different amplitudes and frequencies was simulated. The results showed that compared with the case of no electric field or uniform electric field, the alternating electric field could promote the evaporation of water droplet more significantly. When the electric field frequency was constant, with the increase of the electric field amplitude, the evaporation life of water droplet decreased gradually. Besides, the instantaneous evaporation rate, droplet temperature, and the arrangement structure of water molecules would produce oscillation characteristics with the frequency twice that of the applied electric field at a certain amplitude, and the larger the electric field amplitude, the greater the oscillation amplitude. However, when the electric field amplitude was constant, the evaporation life of droplet did not change monotonically with the increase of frequency, but reached a maximum and a minimum respectively at $f=5$GHz. The reasons for this phenomenon were explained from two aspects of water droplet energy and water molecular arrangement structure.
2021, 53(5): 1324-1333. doi: 10.6052/0459-1879-20-410
Interfaces play a key role in the load transfer of particle reinforced composites, and the interface properties have an important impact on the overall mechanical behavior of composites. However, due to the complex internal structure of composites, it is difficult to determine the interfacial strength and fracture toughness, especially for the respective prediction of normal and tangential interfacial strength. In this paper, a new method is proposed to predict the interfacial mechanical properties of particle-reinforced composites based on zirconia particle (ZrO) reinforced polydimethylsiloxane (PDMS) composite materials. Firstly, the uniaxial tensile stress-strain curves of pure PDMS matrix material and single particle filled PDMS sample are obtained, and the uniaxial tensile constitutive relationship of the PDMS matrix material is achieved. Secondly, the finite element model (FEM) consistent with single particle filled sample is established, and a specific cohesive zone model is chosen to describe the mechanical behavior of interface. Comparison of the experimental and numerical results of tensile mechanical response of samples at different stages can yield the normal strength, tangential strength and fracture toughness of the interface, respectively. The achieved interfacial mechanical parameters are further adopted in FEM calculations of samples filled with different sizes and numbers of particles. The rationality of the proposed method of predicting interfacial properties is verified by comparing the corresponding experimental and numerical results. Such a method is simple, easy to operate and has high accuracy, which should be helpful to quantitatively predict the mechanical property of particle reinforced composites, and also has certain reference significance to quantitatively predict the interfacial properties of fiber reinforced composites.
2021, 53(5): 1334-1344. doi: 10.6052/0459-1879-21-076
High temperature superconducting tapes has been widely investigated in the field of superconducting technology due to its outstanding advantages of high current carrying, low loss and high cost performance. However, the superconducting tape is a multi-layer structure consisting of the Hastelloy alloy substrate, the buffer layer, the YBCO superconducting layer and the protective layer. Its multi-layer structure characteristics including the difference in the properties of each material layer and the complex external field environment will cause local debonding. It is very easy to degrade the transmission performance and mechanical properties. In this paper, the mechanical response of locally detachment YBCO high-temperature superconducting tapes under the excitation of an external vertical magnetic field is studied. Based on the superconducting critical state Bean model and the elasticity plane strain method, the interface between the normal stress of the superconducting film and the substrate is given. The governing equations related to the shear stress, numerically calculated the normal stress in the superconducting film and the shear stress at the substrate interface with the external magnetic field. The results show that the normal stress in the superconducting film and the shear stress at the substrate film interface increase rapidly near the detachment region. At the same time, the maximum shear stress appears at the edge of the structure. The properties of the substrate material, especially the Young's modulus, has a significant effect on the stress in the structure. In the structure of soft substrate material, there is a large normal stress in the superconducting film, and when the Young's modulus of the substrate is larger, a larger shear stress appears at the substrate film interface. These factors will cause the degradation of mechanical and electrical properties of superconducting coating structure. The work may provide some theoretical guidance for fabricating the superconducting tapes and and suppression of layer debonding.
2021, 53(5): 1345-1354. doi: 10.6052/0459-1879-21-043
Studies show that the macroscopic failure of structure results from the microscopic damage within materials, such as cleavage and slip plane decohesion. Therefore, it is helpful to understand the deformation and failure process of materials by studying the damage evolution at micro-scale. Based on the crystal plasticity theory, the microscopic damage in material is studied by analyzing the stress and deformation of slip system, and the microscopic damage model is proposed to consider the resolved normal stress on crystal slip plane. This study provides a new approach for the analysis of cleavage fracture of crystalline materials. First, the gradient tensor of damage deformation is introduced in addition to the crystal elastic-plastic deformation configuration. The constitutive equation with damage energy dissipation is established from the deformation kinematics analysis, and the plastic flow equation and the damage evolution equation are derived. Second, the numerical method is established including the updating algorithm of stress and state variables and the derivation of Jacobian matrix. After that, the single crystal copper with $[100]$ orientation is studied as an example. Through comparing the results obtained by finite element computation and by experimental test, the 11 material microscopic parameters are calibrated using the particle swarm optimization algorithm. Finally, the proposed microscopic damage model is applied to the simulation of RVE under uniaxial tension. The curve of stress versus strain considering the damage effect is obtained, and the development of plastic flow and damage evolution are analyzed. The results show that the proposed model is able to compute the damage accumulation of materials and reasonably reflect the microscopic damage mechanism of crystalline materials.
2021, 53(5): 1355-1366. doi: 10.6052/0459-1879-20-454
During curing of concrete, hydration and thermal transfer inevitably result in expansion and shrinkage and hence, large tensile stresses in early-age concrete structures. As the mechanical properties of young concrete are still very low, structures are vulnerable in the construction stage to defects induced by crack nucleation, propagation and evolution, severely threatening the integrity, durability and safety of concrete structures and infrastructures like nuclear containment vessels, bridges and tunnel linings, hydraulic and off-shore structures. In order to predict the fracture property of early-age concrete and quantify its adverse effects on structural performances, it is pressing to investigate the computational modeling of early-age cracking in concrete structures under the chemo-thermo-mechanically coupled environment. To the above end, in this work we propose a multi-physically coupled phase-field cohesive zone model within our previously established framework of the unified phase-field theory. The interactions between the crack phase-field with the hydration reaction and thermal transfer are accounted for, and the dependence of the characteristics of crack phase-field evolution, e.g., the strength-based nucleation criterion, the energy-based propagation criterion and the variational principle based crack path chooser, etc., on the hydration degree and/or temperature, are quantified. Moreover, the numerical implementation of the proposed model in the context of the multi-field finite element method is also addressed. Representative numerical examples indicate that, with the couplings among hydration, thermal transfer, mechanical deformations and cracking as well as the competition between thermal expansion and autogenous shrinkage both properly accounted for, the proposed multiphysical phase-field cohesive zone model is able to reproduce the overall cracking process and fracture property quantitatively. Remarkably, the numerical predictions are affected by neither the phase-field length scale nor the mesh discretization, ensuing its promising prospective in fracture control of early-age concrete structures.
2021, 53(5): 1367-1382. doi: 10.6052/0459-1879-21-020
The ice-going ship with podded propulsor can avoid maneuverability problems, such as changing direction, turning round and so on, when navigating in the ice-covered waters. Therefore, the podded propulsor is widely used in ice-going ship. In the milling process of ice-podded propulsor, podded propulsor is subjected to extreme ice load, which brings serious harm to the structural strength of podded propulsor and the safety of ice-going ship. In order to study the change regularities of ice load in the process of podded propulsor with different maneuvering states and ice milling, first of all, the theoretical basis of the peridynamic method to study the fracture problem of objects is introduced in detail in the paper, and the feasibility of the method to simulate ice material is analyzed. Based on the coupling of the peridynamic method and the panel method, the material failure criterion for ice breaking simulation and the ice load calculation method suitable are derived. Secondly, a method to judgment the contact between podded propulsor with different maneuvering states and ice in the milling process is proposed. The numerical calculation model of ice-podded propulsor milling state is established, and the numerical simulation of ice-podded propulsor milling dynamic process is realized. Finally, the changes of ice breaking, ice load of propeller and single blade, overall torque of podded unit during podded propulsor with straight ahead, azimuthing condition, dynamic azimuthing condition and ice milling are analyzed. The results show that the contact judgment method proposed in this paper can real simulate the milling process of the ice-podded propulsor, and obtain the ice breaking phenomenon and ice load variation characteristics in the process. It can provide guidance for the development of numerical prediction technology for ice load on marine structure in ice area, the optimization design and operation of ice class podded propulsor.
2021, 53(5): 1383-1401. doi: 10.6052/0459-1879-21-001
The study of the interaction between beams and saturated soils is of great significance in both mechanics and engineering fields. In this paper, the fractional Merchant model is adopted to solve the rheological consolidation of saturated soils, which can simulate the time-depending characteristics of the soils more accurately than the common integer order viscoelastic models. Based on the solution of consolidation for multilayered cross-anisotropic viscoelastic saturated soils, the finite element method (FEM) and the boundary element method (BEM) are coupled to investigate the interaction between beams and fractional viscoelastic saturated soils. The beam is discretized into a number of elements according to the Timoshenko beam theory, and then the global stiffness matrix equation of the beam is obtained. The precise integration solution of the viscoelastic soils is considered as the kernel function of the boundary integral, and the flexibility matrix equation of soils is established by the BEM. Finally, by coupling the FEM and the BEM, the solution of the interaction between multilayered fractional viscoelastic saturated soils and the Timoshenko beam is derived by introducing the displacement coordination condition and equilibrium condition for forces between them. The soil model adopted in this study is degenerated into the Kelvin model, and the results obtained are compared with those in the existing literature, which shows a good consistency. On this basis, the effects of the fractional order and stratification of soils on the interaction between beams and viscoelastic soils are discussed. Numerical results show that: the consolidation velocity of viscoelastic saturated soils with higher fractional order is obviously faster; for layered soils, the reinforcement of topsoil can effectively control the ground settlement and reduce the differential settlement. In practical engineering, the effects of rheology of saturated soils and soil stratification should be well considered to analyze the interaction between beams and soils more accurately.
2021, 53(5): 1402-1411. doi: 10.6052/0459-1879-20-447
Most mechanical vibrations are detrimental that not only generate noise but also reduce the service life and operating performance of the equipment. As two common components, grounded stiffness and inerter can change the natural frequency of the system, which has good effect in the field of vibration control. However, most of the current research only focuses on the impact of a single component on the system, and the vibration absorber is gradually difficult to meet the growth of performance demand for vibration control. Based on the typical Voigt-type dynamic vibration absorber, a novel dynamic vibration absorber model with inerter and grounded stiffness is presented. The optimal parameters of the presented model are studied in detail, and the analytical solution of the optimal design formula is derived. First of all, the motion differential equation of the two degree-of-freedom system is established through Newton's second law, and from the system analytical solution it is found that the system has three fixed points unrelated to the damping ratio. The optimal frequency ratio of the dynamic vibration absorber is obtained based on the fixed-point theory. When screening the optimal grounded stiffness ratio, it is found that the inappropriate inerter coefficient will cause the system to generate instability. Then the best working range of the inerter is derived, and finally the optimal grounded stiffness ratio and approximate optimal damping ratio are also obtained. The working condition when the inerter coefficient is not within the best range is discussed, and the suggestions in practical application are given. The correctness of the analytical solution is verified by numerical simulation. Compared with other dynamic vibration absorbers under harmonic and random excitations, it is verified that the presented DVA can greatly reduce the amplitude of the primary system, widen the vibration reduction frequency band, and provide a theoretical basis for the design of new type of DVAs.
2021, 53(5): 1412-1422. doi: 10.6052/0459-1879-21-058
The dynamical systems with the coupling of different scales observed widly in engineering problems often behave in the bursting oscilltions, characterized by the alternations between large-amplitude oscilltations and small-amplitude oscillations, the generation mechanism of which has been one of the hot topics in nonlinear science at home and abroad. The traditional geometric pertubation method can be employed to explore the mechanism of the oscillations only in the systems with two scales in time domain, which can not be directly used to investigate the interaction between different scales in frequency domain. Meanwhile, most of the results are obtained in the vector fields with codimension-1 fold or Hopf bifurcations. Here we focus on the complicated behaviors in the vector field with codimension-3 fold-folfd-Hopf bifurcation when two scales in frequency domain exist. Based on the normal form as well as its universal unfolding with the nonlinear terms up to the third order, all the possible bifurcations are derived, which divide the two unfolding parameter plane into several regions with different dynamics. By introducing a slow-varying periodic excitation instead of one of the unfolding parameters, two types of routes for the tarjectory visiting those regions can be observed, which may result in four classes of bursting oscillations, i.e., periodic Hopf/LPC, Hopf/LPC/Hopf/LPC, fold/LPC/Hopf/Homoclinic and fold/LPC bursting attractors with the variation of the exciting term. It is found that there may exist delay between the locations of the theoretical bifurcation points and the real bifurcation points on the trajectory. The delay may increase with the exciting amplitude, since the inertia of the movement along the equilibrium states may increase. Especiallty, when the exciting amplitude increases to an extent, the trajectory may pass acorss the corresponding regions before the related bifurcation occurs, which leads to the qualitative change of the oscillations. It is shown that, the slow-fast effect with local bifurcations can be investigated by using the periodic perturbation on the unfolding parameters in the normal form of the vector field, which can therefore to present all possible types of bursting patterns as well as the mechanism.
2021, 53(5): 1423-1438. doi: 10.6052/0459-1879-21-024
2021, 53(5): 1439-1448. doi: 10.6052/0459-1879-21-009
Muscle injury and other diseases often occurs in high-intensity physical workers, so the research on the deformation characteristics and the stress distribution of skeletal muscles are of increasing importance. It is important to obtain the correct constitutive parameters for the study of mechanical behavior of biological soft tissues, and the determination of the constitutive parameters is essentially an inverse process, which possesses challenges. In this paper, two inverse methods based on machine learning are proposed to determine the constitutive parameters, which are k-nearest neighbor (KNN) model and support vector machine regression (SVR) model combined with nonlinear finite element simulation. Firstly, based on the principle of nonlinear mechanics, a finite element model is established to simulate the nonlinear deformation of skeletal muscles under compression, and the corresponding deformation characteristics and stress distribution. At the same time, the dataset of nonlinear relationship between nominal stress and principal stretch of skeletal muscles is established by using the finite element model. Then KNN model and SVR model are used to build the machine learning intelligent algorithms for the inversion of constitutive parameters of skeletal muscle tissues, and the corresponding datasets are trained. Combined with the experimental data of uniaxial compression experiment, the constitutive parameters of skeletal muscles are predicted. Finally, intensive studies also have been carried out to compare the performance of KNN model with SVR model to identify the hyperelastic material parameters of skeletal muscles. And the validity of two inversion methods were verified numerically by introducing the correlation coefficient $(R)$ and the decision coefficient ($R^{2})$. The results show that KNN model and SVR model combined with finite element method are effective and accurate method to identify the hyperelastic material parameters of skeletal muscles. This method can also be further extended for the predictions of constitutive parameters of other types of nonlinear soft materials.
2021, 53(5): 1449-1456. doi: 10.6052/0459-1879-21-038
Roll wave is a kind of unstable gravity driven free surface fluctuation in inclined open channel, which could be classified into two types: periodic roll wave with relatively stable waveform and wave speeds, and irregular roll wave with constantly changing waveform and wave speeds (natural roll wave). The development of irregular roll wave is a complicated process, which is studied however far from enough as compared with the regular roll waves. In this study, a numerical model based on the two-dimensional Reynolds averaged Navier-Stokes equation and the renormalization group $k$-$\varepsilon$ turbulence model was adopted. A large number of numerical simulations and statistical analyses were conducted to provide more dynamical and turbulence information about coalescence processes in irregular roll waves. The evolution processes of both absorption and overtake modes of coalescence were systematically studied to obtain the waveform, wave speeds, velocity profile and turbulence viscosity etc. Results imply that coalescence plays an important role on natural roll wave development, where under certain conditions, roll wave growth changes from regular natural growth mode to irregular mode dominated by coalescence, during which a following wave chases, climbs, and merges with the leading wave successively, along with adjustments of internal flow fields, to finally produce a new combined roll wave with larger wave length and wave height. Furthermore, it is found for the first time that when three roll waves were located closely enough, double coalescences would take place. The two tailing roll waves first coalesce to form a new following wave, which then interacts with the leading wave to produce another stronger roll wave propagating downstream.
2021, 53(5): 1457-1470. doi: 10.6052/0459-1879-20-459
For additively manufactured structure, the poor forming precision and surface roughness may cause surface layer heterogeneity, which leads to uncertain material properties and/or uncertain structure geometry. In order~to obtain a structure with less sensitive to the uncertainty, a rubost topology optimization method accounting for the uncertain surface thickness of structures is proposed, in which two key problems need to be solved to study the surface layer thickness uncertainty caused by the heterogeneity of the structure surface layer. One is accurately identifying the structure surface layer. The other is to carry out propagation analysis and stochastic response estimation of uncertainty. First of all, an erosion-based surface layer identification method is adopted to establishing the equivalent model of surface layer heterogeneity through smooth filtering based on Helmholtz partial differential equation(PDE) as well as discrete mapping based on Heaviside filtering and tanh function, which is called a two-step filtering/projection process. Secondly, while the thickness of the heterogeneous surface layer is regarded as a random variable subject to Gaussian distribution, the uncertain propagation is analyzed and the system stochastic response is predicted based on the perturbation finite element method. Taking the weighted sum of the mean value and standard deviation of structural compliance as the optimization objective, a robust topology optimization model considering the uncertainty of surface layer thickness is established, and the sensitivities of the objective function with respect to design variables are derived. Finally, several numerical examples are given to demonstrate the effectiveness of the proposed method. The numerical results show that the structural configuration with stronger uncertainty resistance can be obtained by considering the influence of surface thickness uncertainty on the structural performance during the design process, which effectively improves the robustness of the structural performance. Therefore, for additive manufacturing structures, it is of great significance to consider the influence of surface layer thickness uncertainty on structural performance in topology optimization design.
2021, 53(5): 1471-1479. doi: 10.6052/0459-1879-20-419
Up to now there have been a dazzling number of artificial boundary conditions (ABCs) in the field of numerical simulation of wave propagation. In order to choose the most appropriate ABCs and assess their performance in complicated wave problems, the related systems of theory and formula that can be used to classify or merge these ABCs still need to be improved. In this work we develop the theory of extrapolation-type artificial boundary condition which can merge a series of classical ABCs that have the common feature that the motion of each artificial boundary node is extrapolated from the motions of a set of adjacent nodes at several previous time steps. These ABCs include Liao's multi-transmitting formula (MTF), paraxial-approximation absorbing boundary conditions, Higdon boundary conditions, the auxiliary-variable-based ABCs of Givoli-Neta, Hagstrom-Warburton or AWWE, et al. Due to the fact that the existing boundary formulas usually have somewhat imperfections, thus we propose two new basic formulas for the extrapolation-type ABCs. One formula is an optimized MTF which incorporates a set of adjustable artificial wave velocities as the control parameters. The other formula is a unified Higdon boundary formula which is defined in a local coordinate system centered at the boundary node and uses various artificial wave velocities as control parameters. The two basic boundary formulas are the most simple and effective extrapolation-type ABCs. Other local ABCs of this type can mostly be transformed from the two basic boundary formulas, or have connections with them via some kind of equivalent intermediate formulas. Numerical experiments are conducted to validate the effectiveness of the proposed theory and boundary formulas. As compared to traditional ABC employing a single artificial wave velocity, the superiority of using ABCs with adjustable artificial wave velocities is preliminarily revealed in this work. It can be expected that the superiority will be more remarkable in simulating complicated wave problems that have several distinct physical wave velocities, such as elastic waves in soft soils with large ratio of longitudinal and transversal wave velocities, dispersive waves in ocean acoustics or meteorology and so forth.
2021, 53(5): 1480-1495. doi: 10.6052/0459-1879-20-408
The existing expansion theory of cylindrical cavity has been able to provide theoretical basis for such as the wellbore stability evaluation in petroleum engineering and side pressure and cone penetration experiment analysis. But it is rarely applied in practical engineering problems such as unsaturated foundation pressure grouting and composite foundation treatment. Based on the theory of elastoplasticity and the principle of unsaturated soil mechanics, this paper adopts the unified strength theory to analyze the problem of cylindrical cavity expansion in unsaturated soil. Firstly, the soil around the cylindrical cavity is divided into the elastic zone and plastic zone and consider following the small strain theory in the elastic zone and the large strain theory in the plastic zone, and considering the influence of intermediate principal stress and inter-grain suction on the strength of unsaturated soil. Secondly, applying the unified strength criterion expressed by effective stress, based on basic equations such as constitutive relations, geometric equations, and momentum balance equations, combined with the corresponding boundary conditions. Finally, the analytical expressions of the stress field, strain field, displacement field and limit reaming pressure in the surrounding elastic-plastic region when the cylindrical cavity expands under different drainage conditions are obtained. Through numerical examples and parameter analysis, while degenerating verification with the existing expansion theory of cylindrical cavity in saturated and unsaturated soils, the influence laws of suction, dilatancy parameters, intermediate principal stress effect parameters and initial radial effective stress on the stress field, strain field and displacement field in the elasto-plastic region are analyzed. So as to verify the correctness and effectiveness of the theory in this article, in order to provide a reasonable theoretical basis for the later practical engineering problems.
2021, 53(5): 1496-1509. doi: 10.6052/0459-1879-20-439

This article for the first time publishes a journal review opinion given by the editor-in-chief of Chinese Journal of Theoretical and Applied Mechanics (CJTAM), Professor Qian Xuesen in 1961. It carefully interprets Professor Qian's guiding principle for running the journal and reviewing academic papers, as well as his high standards for journal publication. This article also introduces the historical background of Professor Qian’s review comment and its effect and influence. The author hopes that by reviewing this historical past we can revisit Professor Qian’s scientific spirit for running scientific journal and for promoting theoretical mechanics. The author is motivated to promote the legacy of CJTAM and developing unique journal features in the practice of running scientific and technological journals in the new era.

2021, 53(5): 1510-1514. doi: 10.6052/0459-1879-21-117