Ruchao Shi and Xiaowang Sun

（1. School of Science, Jiangsu Ocean University, Lianyungang 222005, China；2. School of Mechanical Engineering, Nanjing University of Science and Technology, Nanjing 210094, China）

Abstract: Elastic wave refraction at the air-solid interface and wave propagations in the vicinity of the air-solid interface are numerically studied. The modified ghost fluid method (MGFM) and isobaric fix methods are combined to solve the fluid and solid statuses at the air-solid interface and construct a continuous boundary condition for the air-solid interface. The states in the ghost domain are evaluated by the MGFM-algorithm. The solid governing equations are solved with second order spatial discretization. Numerical tests verify the correctness of the presented continuous boundary condition and show that the combined method is convergent in the vicinity of the air-solid interface. The 3D numerical results by the combined method are close to those of the Arbitrary-Lagrangian-Eulerian (ALE) method. The combined method is robust when applied for multi-dimensional problems. A compression stress wave impacting on the air-solid interface result in a compression wave in air. A tension stress wave impacting on the air-solid interface result in an expansion wave in air.

Key words: fluid-structure interaction；finite difference method；level set method；continuity condition；stress wave propagation

Elastic stress wave refraction at the air-solid interface appears in many research areas, such as in material mechanics[1?2], measuring techniques[3?4], military engineering[5?7], tunnel engineering[8], maritime engineering[9], and so on. The difficulty in numerical simulations of stress wave refraction at the air-solid interface is in solving the interfacial status between air and solid.Presently, several methods are available for this kind of problem. One of them, the Arbitrary-Lagrangian-Eulerian (ALE) method[10], is based on the finite element method and body-fitted meshes and has the widest range of applications. The second one is the smoothed particle hydrodynamics method (SPH). Recently, it has been shown that this method is powerful when applied to fluid-structure interactions, although great quantities of computational resources are required[11?12].Another successful method is the ghost material method (GMM)[13]. Recent improvements and increasing applications demonstrate that this kind of method is becoming more and more mature[14].Another emerging method is the immersed boundary method[15]. The latest review on this method can be found in Ref.[16]. There are other coupled methods that work successfully for fluidstructure interactions. For instance, the Lattic Boltzmann method can be combined with the immersed boundary method[17?19], the finite element method and the SPH method can be combined[20?21], the SPH method can be coupled into the ALE method[22], etc.

Among the numerical methods mentioned above, the GMM has been verified to be robust,simple, and convenient for multi-material problems. The GMM comprises several numerical algorithms: the original ghost fluid method (OGFM-algorithm)[23], the modified ghost fluid method (MGFM-algorithm)[24], and the real ghost fluid method (rGFM-algorithm)[25]. The MGFM-algorithm can avoid most numerical oscillations at material interface. For numerical simulations of multi-material flow, including gas-water flow and fluid-structure interaction, the numerical errors of entropy at the material interface often accumulate very fast. Cumulative entropy errors lead to numerical oscillations of density or temperature at the material interface. These numerical oscillations of density or temperature are figuratively called the “overheating effect”. To overcome the “overheating effect”, the isobaric fix method was developed. The original idea of the isobaric fixing method was to use the entropy of the real fluid near the interface to repair entropy errors at the interface[26]. In practical simulations,the use of isobaric fixing is not changeless. The isobaric fixing method has been verified and validated to be efficient and accurate for the gas-water interface and solid walls.

In this work, we combine the MGFM-algorithm and the isobaric fixing method to simulate elastic stress wave refraction at the air-solid interface. In Section 1, governing equations for air and solids are elaborated. In Section 2, the isobaric fix method is coupled into the MGFMalgorithm. Then, specific implementation steps are given in Section 3. Next, numerical results and comparisons with different methods are presented. Finally, short conclusions are summarized in the last section.

1 Governing Equations

1.1 Euler equations

1.2 Elastodynamics equations

For 3D problems, the governing equation for solids is expressed as[28]

1.3 Equations of state

γ-law is used to solve air,

wherepis pressure,ρis density,eis internal energy per unit mass, andγis the specific heat ratio of gas.

1.4 Level set method

2 Combining the MGFM-Algorithm and the Isobaric Fixing Method

2.1 Interfacial continuity condition

2.2 Approximate Riemann problem solver

2.3 Isobaric fix method

2.4 Evaluation algorithm

In Fig. 2, the MGFM evaluation algorithm is shown. In Fig.2a, ghost cells are distributed for the solid. Interfacial stressσIsand velocityusIare solved by Eqs.(13)(19)(20), and then extended into the ghost domain by the constant extrapolation method (Eq.(11)). Likewise, in Fig. 2b,ghost grid nodes are distributed for air. The interfacial status, pressure, density, and velocityare solved and then extended into the ghost domain by Eq.(11).

Fig. 2 Schematics of evaluation algorithm.

If the problem is multi-dimensional, interfacial velocities (and) are replaced by interfacial normal velocities (and). For the solid,the tangential velocityand the shear stress τ at the cell center just next to the interface are extended into the ghost domain of the solid. For air, the tangential velocityof the grid node of the air near the interface is extended into the ghost domain of air.

3 Numerical Implementation

① Define the air-solid interface using the level set function.

② Seek pairs of grid nodes,and, to construct pairs of Riemann problems.

③Solve the interfacial status for the next time step using the approximate Riemann problem solver and the isobaric fixing method.

④Evaluate the grid nodes of air and cells of the solid using the MGFM evaluation algorithm.

⑤Solve the flow field, the status of the solid,and the location of the air-solid interface for the next time step.

4 1D Applications

Case 1. As shown in Fig.3, a stress wave in the solid propagates towards the air-solid interface, resulting in a reflected wave in the solid and a transmitted wave in air (see Fig.4). The computational domain is set as [–1,1]. The boundary conditions at the two ends are defined as infinite.801 nodal points of air and 800 cells of solid are distributed. The air-solid interface is located atxI=0.

Fig. 3 Incident wave on the side of solid propagating towards interface

The convergent errors for air and the solid are tested and presented in Tab.1. Supposing the numerical algorithm isn-th order, for two different step sizes Δx1and Δx2, we have Error1=O() and Error2=O(). Dividing Error1by Error2, one can solve for the numerical ordern. It can be found that the MGFM-algorithm coupled with the isobaric fixing method has first order numerical accuracy.

Fig. 4 Reflected wave on the side of solid and transmitted wave on the side of air

Fig. 5 Numerical results on the reflected wave in solid and comparisons with exact solutions at t=6.928×10-4(Case 1)

Case 2. We set the incident wave in Fig.3 as a tension wave. The boundary conditions at the two ends are defined as infinite. The incident wave is located initially at the air-solid interface.The initial status on the side of the solid are set as=?5.0×10?1and=1 309.155. The initial conditions on the side of air are=0.0,=1 000.0 , and=1.80×10?1. The physical parameter of the gas isγ=1.4. The physical parameters of the solid are the elastic modulusYs=2.151×106, the Poisson ratioνs=0.283, and the densityρs=7.7. The statusesandcan be solved theoretically as==993.169 and==4.320×10?1. The CFL numberμis set to be 0.9. From Fig.7, there are good concurrences of results between the presented method and theoretical analysis. In Fig.7, the density of the solid is not plotted because it remains constant (ρs=7.7) and is much higher than the density of the air. According to the results of stress and pressure, it can be found that a tension wave in the solid impacting on the air-solid interface results in an expansion wave in air.

Fig. 6 Numerical results on the transmitted wave in air and comparisons with exact solutions at t=1.108×10-4 (Case 1)

Tab. 1 Numerical accuracy of MGFM-algorithm coupled with isobaric fix method at t=7.993×10-3

Fig. 7 Comparison of results between the presented method and theoretical analysis at t=1.125×10-3 (Case 2)

Fig. 8 Wave penetrates two air-solid interfaces (Case 3)

5 Multi-Dimensional Problem

Fig. 9 Numerical results of wave propagation for comparisons against theoretical results at t=5.892×10–3(Case 3)

Fig. 10 Numerical results for comparisons against theoretical results at t=7.576×10–3 (Case 3)

Fig. 11 Schematic of ship cabin (Case 4)

Fig. 12 Stress distribution in ship hull

Fig. 13 Propagation of transmitted compression wave on the side of air

Fig. 14 Comparisons of dynamic responses between the combined method and ALE method of ANSYS software

6 Conclusion

In this work, the MGFM-algorithm is combined with the isobaric fixing method. Tests on the numerical errors demonstrate that the combination of these two methods is convergent. If a compression stress wave refracts at the air-solid interface, the transmitted wave on the side of air is also a compression wave. If a tension stress wave refracts at the air-solid interface, the transmitted wave on the side of air is an expansion wave. The combined method is applied for a 3D stress wave penetrating a ship hull. The 3D numerical results on the reflected wave are compared between air-solid interface and free surface and also between the combined method and the ALE method. The combined method is robust when applied for multi-dimensional problems.