Mathematical Problems in Engineering

Volume 2015 (2015), Article ID 438086, 17 pages

http://dx.doi.org/10.1155/2015/438086

## A Simple Immersed Boundary Method for Compressible Flow Simulation around a Stationary and Moving Sphere

^{1}Tokai University, 4-1-1 Kitakaname, Hiratsuka, Kanagawa 259-1292, Japan^{2}ISAS/JAXA, 3-1-1 Yoshinodai, Chuo-ku, Sagamihara, Kanagawa 252-5210, Japan

Received 19 June 2015; Accepted 6 September 2015

Academic Editor: Yihua Lan

Copyright © 2015 Yusuke Mizuno et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

This study is devoted to investigating a flow around a stationary or moving sphere by using direct numerical simulation with immersed boundary method (IBM) for the three-dimensional compressible Navier-Stokes equations. A hybrid scheme developed to solve both shocks and turbulent flows is employed to solve the flow around a sphere in the equally spaced Cartesian mesh. Drag coefficients of the spheres are compared with reliable values obtained from highly accurate boundary-fitted coordinate (BFC) flow solver to clarify the applicability of the present method. As a result, good agreement was obtained between the present results and those from the BFC flow solver. Moreover, the effectiveness of the hybrid scheme was demonstrated to capture the wake structure of a sphere. Both advantages and disadvantages of the simple IBM were investigated in detail.

#### 1. Introduction

Acoustic wave from rocket nozzle on the rocket launch may critically damage satellites inside the payloads of the rocket since the wave sometimes becomes quite strong. The prediction of the crucial acoustic field around rocket launching sites has been conducted by the past enormous studies and testing data [1]. However, more accurate and reliable prediction is required to construct a new rocket launch site more safely. Various numerical studies have been conducted to realize qualitative prediction [2–4] until now. In real large-scale liquid launch vehicles, generated acoustic waves are suppressed by water injection. For high accurate prediction of the acoustic wave, the interaction between the water droplets and the flows should be clarified from the physical aspect.

Fukuda et al. showed that acoustic waves are primarily attenuated by the interactions between particles and turbulence [5] from theoretical analysis. These particles that interfered with turbulence in their study are alumina particles exhausted from the rocket nozzle or water droplets from water injection at the rocket launch. However, the mechanism of the phenomenon for the attenuation has not been well understood.

We develop the Cartesian flow solver with an immersed boundary method (IBM) [6] to investigate the interactions between particles and shocks by the direct numerical simulation (DNS). The objective of this study is to develop the flow solver to investigate the supersonic flow as rocket plumes around micro-meter-sized particles like alumina particles or water droplets. Accordingly, the flow becomes high-Mach- and low-Reynolds-number flow. To resolve real interaction between flows and particle in this study, an Euler-Euler approach is employed in the equally spaced Cartesian mesh with sharp interface IBM techniques. A hybrid scheme is also utilized to deal with flows of both shocks and vortical wake around a fixed or moving sphere. A drag coefficient is compared with numerical results obtained from the high accurate boundary-fitted coordinate (BFC) grid solver [7] and theoretical values [8]. We extended two-dimensional numerical methods [6] to three-dimensional one to solve flows around multiple particles passing through shock waves. The aim of this study is to investigate gas-particle flow involving real interaction between flows and particles by large-scale flow simulation. One of the most significant points for this purpose is to keep the algorithm simple for optimization and fast computation. Both simplicity and applicability of the code should be preserved in this study.

#### 2. Numerical Method

##### 2.1. Governing Equations

The governing equations of the present flow simulations are three-dimensional compressible Navier-Stokes equations:where , , and are inviscid fluxes in the , , and directions, , , and are viscous fluxes, and contains conservative variables. The pressure is related to the total energy per unit mass by an equation of state:Here, viscous stress tensor components are given as

All variables are normalized by the freestream values of the density, the sound speed, and the reference length. The above equations are discretized on an equally spaced Cartesian mesh with cell-centered arrangement. To eliminate additional numerical dissipation everywhere, except in the vicinities of shock waves and potential flows, the inviscid terms are computed by a hybrid scheme that combines the third-order monotone upstream-centered scheme for conservation laws- (MUSCL-) Roe scheme [9] and the second-order pseudo skew-symmetric scheme [10].

In our previous studies, we proposed the hybrid scheme estimated by the Ducros sensor [11] , pseudo skew-symmetric scheme , and MUSCL-Roe scheme as follows:Here, the subscript denotes cell interface to evaluate the numerical flux along with -axis in this study. The sensor in (4) for switching the scheme is estimated by

Here, = 10^{−15} is a small value that prevents division by zero, and = 0.6 is the switching threshold. The vector is a velocity vector . The numerical fluxes for and in (4) are calculated by the formulas of (6) and (7):

Here a value is set to be 0.05 to suppress the spurious numerical oscillationwhere is the flux Jacobian and the subscript Roe denotes the Roe-averaged quantity of the left and right states. Herewhere and are the right and left eigenmatrices of , respectively, and is a diagonal matrix.

In the previous study [6], the symmetric central difference parts of (7) that are first and second terms of right hand side were replaced by the ones of the pseudo skew-symmetric form without losing the formal order of accuracy of the equation. However, in this study, the skew-symmetric form is not applied to (7) from the preliminary investigation. Detailed information for the investigation is described in Appendix A.

##### 2.2. Immersed Boundary Method

A simple IBM based on the level set function is used to represent object boundaries in this study. The level set function is determined as a signed distance of minimal distance from the object surface. In the case of multiple objects, multiple level set functions are calculated based on simple minimum distance approach [12]. Whole cells are classified into three categories of fluid cells (FC), ghost cells (GC), and object cells (OC) by using the level set function expressed as (9). The ghost cells behave as guard cells between the fluid cells and object cells and are assigned for two layers under the present definition:

Figure 1 shows distributions of the classified cells around a cut sphere in the present IBM. The blue, brown, and green regions are assigned to the fluid, ghost, and object cells. The blue, brown, and green hemispheres are surfaces of the image points (), zero level set (), and the threshold of the object cells (). Here, the circle symbols of GC, IP, and IB represent a cell centroid of ghost cell, an image point, and an immersed boundary point. Developed numerical method in the previous study [6] in 2D is extended to 3D directly. Flow quantities of ghost cells are determined by using those of the image points. The flow variables of the image points at the edge of a “probe” are defined by trilinear interpolation from surrounding fluid cells. The probe is a vector from GC through IB to IP in the direction normal to the zero level set surface. The probe length expressed as in (10) is the fixed value of 1.75 times of the mesh size to avoid recursive reference in this study. The flow variables of the ghost cell are determined by assumed distributions on the probe. In this study, three velocity components of the ghost cell are calculated by the linear extrapolation while the density and the pressure are determined by Neumann condition along with the probe expressed as