Solving the World’s Toughest CFD Problems

# Immersed Boundary Method # Immersed Boundary Condition

The accurate prediction of forces and energy losses is crucial in the success of modeling many engineering problems such as discharge from an orifice plate, flow past an obstacle and flow in a sudden contraction pipe. Since the fractional cell area and volume method FAVOR™ was introduced more than 30 years ago, the standard solver in FLOW-3D has used a simple method to approximate momentum fluxes near walls (Hirt and Sicilian 1985). When computing the momentum advection term near walls or free surfaces, any spatial velocity derivative that requires velocity values located within solid or void regions is set to zero to eliminate the appearance of numerical boundary layers. From a physical point of view, this method applies a free-slip (non-penetration) boundary condition to the advection at walls, suppressing the artificial boundary layer (Hirt 1993). The loss of flux in the momentum equation is compensated by the pressure. In certain situations, the portion of pressure to compensate for the flux loss grows in time, and may result in a numerical instability called “secular instability,” which is expressed as a monotonic growth of velocities. In order to prevent this instability from developing, an empirical technique has been used to “correct” the flux at the location where the instability may emerge. However, this approach does not resolve the flux loss from the source, and could occasionally introduce a nonphysical behavior of the solution such as pressure oscillations.

A technique for approximating the advection term based on the ghost-cell immersed boundary method (Mittal et al. 2008) has been developed for FLOW-3D v12.0. The immersed boundary method technique provides more accurate pressure and force predictions near walls. The ghost-cell immersed boundary method has emerged in recent years as an enhanced boundary treatment in conventional Cartesian grid approximations in problems involving complex geometries. Since the method is only a means to handle boundaries, it can be easily added as a model into the existing FLOW-3D solver with relatively few changes to the underlying solver structure, and is compatible with other physical models in FLOW-3D. It uses a weighted-average probe technique in conjunction with various interpolation methods to handle different geometry configurations. The new model works for flows with 3D mesh blocks or hybrid 3D/shallow water mesh blocks, but not for shallow water meshes alone.

## Immersed Boundary Method Results

A simple example of the new Immersed Boundary Method model is the discharge of water from a 1m-diameter circular orifice. The water container is 10m long by 10m wide, and the initial elevation of the water surface to the center of the orifice is 6m. The drop of the elevation, h, and the volume flow rate at the orifice, Q, follow the quadratic and linear curves, respectively, as shown in the animation. Based on the equation $latex {{C}_{d}}=\frac{{4Q}}{{\pi {{d}^{2}}\sqrt{{2gh}}}}$, the average coefficient of discharge Cd for the simulation is 0.660, which is about 8% larger than the asymptotic value of 0.611 (Swamee and Swamee, 2010).

Another example is the calculation of the total resistance force on a NAVY vessel model hull. In this case, the hull length is 5.72m, and the draft is 0.248m. Based on the hull length and mean flow velocity of 2.10m/s, the Reynolds number in the simulation is approximately 12×106. Because the problem is symmetric, only half of the hull is modeled. The computational domain is 30m long by 8m wide and 5.5m deep. The experimentally obtained average total resistance force for half of the hull is 22.62N (Larsson et al. 2003). The average force was 22.43N, only 0.8% lower for the immersed boundary solver. In addition, the immersed boundary solver finished in about 40 hours, which is 8 hours faster than the standard solver.

### References

Hirt, C., & Sicilian, J. (1985). A porosity technique for the definition of obstacles in rectangular cell meshes. International Conference on Numerical Ship Hydrodynamics, 4th. Washington, D.C.

Hirt, C. (1993). Volume-fraction techniques: powerful tools for wind engineering. Journal of Wind Engineering and Industrial Aerodynamics, 46 & 47, 327-338.

Mittal, R., Dong, H., Bozkurttas, M., Najjar, F., Vargas, A., & von Loebbecke, A. (2008). A versatile sharp interface immersed boundary method for incompressible flows with complex boundaries. Journal of computational physics, 227(10), 4825-4852.

Swamee, P., & Swamee , N., (2010). Discharge equation of a circular sharp-crested orifice. Journal of Hydraulic Research, 48(1), 106-107.

# Immersed Boundary Method

Zongxian Liang, Developer and Adwaith Gupta, CFD Engineer at Flow Science, contributed this blog.

The accurate prediction of forces and energy losses is crucial in the success of modeling many engineering problems such as discharge from an orifice plate, flow past an obstacle and flow in a sudden contraction pipe. The upcoming release of FLOW3D v12.0, features a new numerical option, the Immersed Boundary Method, that accurately predicts flows for such problems.

Immersed Boundary Method, or precisely, ghost-cell based Immersed Boundary Method, improves the accuracy of numerical flux calculations at the solid-fluid interfaces. In this blog, Zongxian Liang, developer at Flow Science provides validation examples of a sudden contraction pipe and a ship hull where the Immersed Boundary Method improves the accuracy of head loss and resistance force. Brief mathematical details of the Immersed Boundary Method and the ghost cell approach are discussed at the end of the blog.

## Sudden Contraction Pipe

One of the fluid problems that demonstrate the accuracy of the immersed boundary solver is estimating the head loss of water in a contraction pipe, as shown in Figure 1. The pipe experiences a sudden contraction from a large section of 3m in diameter to a small section of 1m. Flow rate at the inlet of the large pipe is 4m3/sec. Velocity heads are measured 3.5m upstream and downstream with respect to the contraction location. Based on the different assumptions used in the literature [1,2], the theoretical value of the head loss should be in the range of 0.494m to 0.711m. In this simulation, one Cartesian mesh block is used for the whole geometry and a pressure boundary is specified with zero stagnation pressure at the exit. The two-equation k-ω model is enabled with maximum turbulent mixing length set to 0.07m. Gravity is active in the negative x direction. The finish time is set to 15s at which the flow is considered to be steady and fully developed.

A two-variable parametric study was conducted to investigate the effect of cell size and mesh alignment with respect to the sudden contraction location. The cell sizes of the mesh are 0.1m and 0.05m. Four mesh alignments were tested, as listed in Table 1: “Aligned” denotes the baseline case where the mesh plane is aligned with the sudden contraction location, and “Z-shifted: X%” indicates that the mesh plane is shifted an X% of the cell size away from the contraction location in the z-direction.

The predictions given by the immersed boundary solver all fall in the theoretical range of 0.494 to 0.711. Especially for the cases with a smaller cell size of 0.05m, the head losses are within 4.6% of the median of the theoretical values, 0.603.

 Grid Alignment Aligned Z-shifted: 25% Z-shifted: 50% Z-shifted: 75% Cell size (m) 0.1 0.706 0.710 0.647 0.666 0.05 0.607 0.620 0.575 0.589

Table 1: Head loss (m) predicted by the immersed boundary solvers

## Resistance Force on a Model of Ship Hull

The accurate prediction of forces in external flow dynamics problems is typically of interest in design phases. An example of an external flow problem with a free surface is the calculation of the total resistance force on a NAVY vessel model hull. In this case, the hull length is 5.72m, and the draft is 0.248m. Based on the hull length and mean flow velocity of 2.10m/s, the Reynolds number in the simulation is approximately 12×106. Because the problem is symmetric, only half of the hull is modeled. The computational domain is 30m long by 8m wide and 5.5m deep, and has three nested mesh blocks with the smallest cell size equal to 0.02m. Velocity and outflow conditions are used for the inlet and outlet boundaries, respectively. RNG turbulence model with dynamic turbulent length scale calculation is used to model the turbulent flow, and the second-order, monotonicity preserving scheme for the momentum advection.

The total resistance force is 22.62N for half of the hull model, based on a total drag coefficient of 0.0423 in experiments . The force from the simulation is the sum of the shear and pressure forces in the x direction, averaged from 35 to 50s. The immersed boundary solver gives 22.43N, 0.8% lower than the experimental result.

## Conclusion

We see from these two validations that the Immersed Boundary Method approach estimates the head loss and the force correctly within the theoretical range provided in the literature. Immersed Boundary Method improves the accuracy of estimating fluxes at a very fundamental level and is expected to improve simulation results for many, if not most, applications. To know more about the ghost-cell based Immersed Boundary Method development, continue reading. Otherwise, stay tuned for our next blog post!

## Ghost-Cell Based Immersed Boundary Method

In FLOW-3D, the free-slip boundary condition is applied to the advection of velocities to eliminate numerical boundary layers caused by fractional cell areas and volumes near solid boundaries. In order to estimate a rational flux for the control volume, the immersed boundary solver computes a fluid velocity in the solid that implicitly satisfies the boundary condition, as illustrated in Figure 2. The fluid cell in the solid is called ghost cell and this method is generally referred to as the ghost-cell approach. Figure 2. Flux at the left face of control volume (enclosed by blue dash-dot lines) is computed using the velocity of the ghost cell u_(i-1) inside the solid.

To enforce the boundary condition, an image-point (open diamond denoted by IP) of the ghost cell is created in the fluid region by extending a line segment from the ghost cell (red diamond, denoted by GC) to the wall in the normal direction and that intersects the wall at the midpoint (open circle denoted by BI for boundary-intercept point) between GC and IP. For a non-penetrating boundary condition at BI and assuming the tangential velocities at GC and IP are the same as the wall surface velocity, the velocity of the ghost cell is estimated by the following equation:

$latex \displaystyle {{\vec{u}}_{{GC}}}={{\vec{u}}_{{IP}}}-2\left[ {\left( {{{{\vec{u}}}_{{IP}}}-{{{\vec{u}}}_{{BI}}}} \right)\cdot \vec{n}} \right]\vec{n}$

where $latex \displaystyle {{\vec{u}}_{{GC}}}$, $latex \displaystyle {{\vec{u}}_{{IP}}}$  and  $latex \displaystyle {{{{\vec{u}}}_{{BI}}}}$ are the fluid velocity at the ghost cell, the image-point and the boundary-intercept point, respectively, and   is the unit normal vector at the boundary. We use trilinear interpolation to approximate the value of velocity  using fluid velocities in the neighboring cells:

(2)    $latex \displaystyle {{\vec{u}}_{{IP}}}=\sum\limits_{{i=1}}^{8}{{\left( {{{\alpha }_{i}}{{u}_{i}}} \right)}}\hat{i}+\left( {{{\beta }_{i}}{{v}_{i}}} \right)\hat{j}+\left( {{{\gamma }_{i}}{{w}_{i}}} \right)\hat{k}$

where u1v1 and w1 are the velocity of interpolation nodes surrounding the image-point and α1β1 and γ1 are interpolation coefficients. The details of computing the interpolation coefficients can be found in Ref. 4.

The interpolant invokes the velocity value of other ghost cells, and a coupled system of the ghost cells is thus generated. We use a Jacobi based iteration method in conjunction with convergence acceleration techniques to obtain a fast solution of the resulting system.

### References

1. White, F. M., Fluid Mechanics (McGraw-Hill Book Company, 2003).
2. Saleh, J., Fluid Flow Handbook (McGraw-Hill Professional, 2002).
3. Larsson, L., Stern, F. & Bertram, V., Benchmarking of computational fluid dynamics for ship flows: the Gothenburg 2000 workshop. Journal of Ship Research 47 (1), 63–81 (2003).
4. Mittal, R. et al., A versatile sharp interface immersed boundary method for incompressible flows with complex boundaries. Journal of computational physics 227 (10), 4825-4852 (2008).