## Abstract

A modified boundary condition for the distribution function in the Lattice Boltzmann method at the interface between solid and fluid that takes into account a finite mass ratio between two phases and an inelastic reflection is proposed. The new boundary condition is built into the immersed boundary method to compute the interaction between graphene and laminar flow. Numerical simulations are carried out at a Reynolds number of 40, and drag and lift acting on the graphene and its deformation are examined by changing the mass ratio and the coefficient of restitution. It is found that the amplitude of the oscillating motion of the graphene is enhanced when compared with the case of infinite mass ratio with perfect collision.

## 1. Introduction

In the numerical computation of a viscous fluid motion, it is very common that the fluid domain is bounded by a solid wall and that the viscous fluid adheres to the wall so that the fluid velocity is equal to the wall velocity, ** v**(

*x*_{B},

*t*)=

*U*_{B}, where

*x*_{B}and

*U*_{B}are the position and velocity vectors of the boundary, respectively. In the Lattice Boltzmann method, this macroscopic boundary condition is translated into the distribution function of an ensemble of fluid particles and is achieved by imposing a bounce-back condition or redistribution under the local thermal equilibrium. Ladd [1] proposed that when the solid is moving with the velocity

*U*_{B}, the condition for the distribution function at the solid surface is given by 1.1 where (

^{∼}) means the post-collision state at the boundary,

*c*_{i}is the lattice velocity vector that points to the

*i*th direction,

*j*is an index satisfying

*c*_{j}=−

*c*_{i}and

*E*

_{i}is given by

*E*

_{0}=1/3,

*E*

_{1}=⋯=

*E*

_{6}=1/18,

*E*

_{7}=⋯=

*E*

_{18}=1/36, with

*c*

_{s}being the speed of sound,

*c*

^{2}

_{s}=

*c*

^{2}/3 for the D3Q19 model [2]. After taking the first moment of

*f*

_{i}with the lattice velocity

*c*_{i}in equation (1.1) and considering the relationship

*c*_{i}+

*c*_{j}= 0, we obtain the well-known post-collision state of the momentum 1.2 where

*ρ*

_{f}is the fluid density and subscripts ‘

*α*’, ‘

*β*’, ⋯ refer to the spatial directions (

*x*,

*y*,

*z*), respectively. Equation (1.2) is the same as the momentum exchange equation between one particle with a velocity

**and the other particle with velocity**

*u*

*U*_{B}under the assumption that the mass of the particle with velocity

*U*_{B}is very large and that the particles always undergo perfect elastic collisions. The model has successfully been applied to many problems of the flow with rigid objects.

When the solid moves and/or is deformed as a result of momentum exchange between fluid and solid, such as a feather in the air, the mass ratio of the two phases should affect the dynamics. Moreover, microscopically, fluid molecules interact with the solid molecules through the potential, and they may transfer the energy from one phase to the other, which appears macroscopically as inelastic collisions. This view has increasingly become important in microfluidics, in which the characteristic length scale of a system has become smaller and smaller owing to the dramatic developments in fine-scale technology. For example, in the study of motion of graphenes at micron scale in a laminar flow [3,4], it was argued that the graphene is expressed in terms of coarse-grained particles that are generated by an appropriate statistical average under a given set of constraints and that they interact with the fluid. Although the obtained results are reasonable, for the precise prediction of the deformation of the graphene, there are demands to know the effects of the finite mass ratio between the two phases and the imperfectness of the collision.

In this paper, we propose a modified boundary condition for the distribution function in the Lattice Boltzmann method at the interface, which includes the effects of the finite mass ratio between the fluid and solid and the imperfectness of the collision. We apply it by using the immersed boundary method [5] to the problem of the motion of graphene in a laminar fluid.

## 2. Interface modelling

For the purpose of developing a new boundary condition for the distribution function in the Lattice Boltzmann method at the interface, it is sufficient to assume that a solid phase is described by a set of hypothetical particles and that the interface is expressed by those particles at the solid surface, although we consider the interaction between coarse-grained particles for graphene and air in §3. We also use the immersed boundary method to compute the effects of the solid phase on the fluid motion.

It is appropriate to briefly explain the graphene expressed in terms of the coarse-grained particles. The graphene is a single layer of carbon atoms connected by an sp^{2}-bond, the characteristic length of a unit cell is of the order of 10^{−10} m and the characteristic time of atomic vibration is 10^{−14} s [6]. It is impossible to directly and simultaneously compute the dynamics of carbon atoms of the graphene and the fluid because there exists a huge gap in the characteristic length and time scales of these atoms and fluid motion. In order to bridge this gap, Ogata and co-workers have developed a coarse-grained particle method in which the Hamiltonian for the carbon atoms with crystalline structure is statistically projected onto the one for a small number of consolidated atoms. The statistical projection is made by applying the phonon approximation to the Hamiltonian and by taking the ensemble average with a given set of constraints in the thermal equilibrium. The coarse-graining process is made for a relatively small number of particles at one time and followed by recursive operations on the resultant system until the system size reaches the desired size, the length scale at which the fluid description is meaningful. Starting from the original carbon atoms of the graphene, we have successfully obtained the coarse-grained Hamiltonian for the graphene. Indeed, it is found that the elastic moduli of the coarse-grained graphene agree well with those computed directly by molecular dynamics (MD), except the term of shear–strain correlation *C*_{xyxy}, which is about 7 per cent smaller than that using MD. Therefore, it should be stressed that the property of the graphene represented in terms of the coarse-grained particles has a direct connection to the actual graphene and is properly built into the equation of motion of the coarse-grained particles. For more details, the reader may refer to Inoue *et al.* [3,4] and Kobayashi *et al.* [7].

The collision of fluid and solid is assumed to be described as an interparticle collision (figure 1). Mass and velocity of the fluid particle as a working object are given by sums of those of the fluid molecules within a rectangular box, of which the height is the same as the fluid lattice and the bottom face is as large as the nearest solid particle cell (figure 1)
2.1
where *m*_{f} is the mass of the fluid molecule (mass number^{−1}) and *N*_{f} is the so-called number density of the fluid molecule (number volume^{−1}). Similarly, and denote the mass and velocity of the particle at the solid surface, respectively.

In the collision, momentum is conserved, 2.2

Change of the relative velocity of the particles in the post-collision state is macroscopically expressed in terms of ‘a coefficient of restitution’ *e* (0≤*e*≤1) between fluid and solid particles, which is defined by
2.3
Substituting equation (2.3) into equation (2.2) and dividing both sides by , the post-collision state of the fluid velocity can be obtained as follows:
2.4
where is the mass ratio.

Equation (2.4) is a general form of the momentum exchange at the interface. In the case of a perfectly elastic collision, *e*→1, equation (2.4) can be rewritten as
2.5
and when the fluid mass is light enough to be neglected compared with that of the solid, meaning that *χ*→0, equation (2.4) is
2.6
If the above two conditions are satisfied simultaneously, equation (2.4) becomes
2.7
and equation (1.2) is recovered when both sides are multiplied by *ρ*_{f}.

A simple way to construct an improved model for the distribution function, which is consistent with equation (2.4), is to set
2.8
It is easily shown that the first moment of equation (2.8) yields
2.9
In the immersed boundary method, the condition equation (2.8) is applied at the particle position *X*_{B} on the solid surface and *f*_{i}(*X*_{B}) is interpolated from those on the Lattice Boltzmann method grid points. The body force mimicking the existence of the interface is computed from the momentum exchange between the particle and fluid and redistributed on the Lattice Boltzmann method grid points nearby. The reactive part of the body force acts on the particles on the surface [3].

## 3. Results

We consider the motion of a graphene of height *H* and width *W*, which is fixed perpendicularly at the bottom plate in a laminar flow. Incoming flow is assumed to be steady Blasius-type flow for which the fluid velocity is zero at the bottom and gradually increases and becomes uniform at *y*>0.1*L*, where *L* is the length of a computational domain in the wall-normal direction. The computational domain for the fluid phase is a rectangular box whose size is 2*L*×*L*×1.5*L* and is divided into 120×60×90 cubes in the streamwise (*x*), wall-normal (*y*) and spanwise (*z*) directions, respectively. Details of the numerical set-up is described in Inoue *et al*. [3]. The Reynolds number is defined as , where is the maximum speed at the entrance and *ν* is the kinematic viscosity. The Reynolds number is 40, which corresponds to , where *c* is the characteristic speed of the fluid phase. In order to obtain a well-posed initial flow field, the Lattice Boltzmann equation is integrated for 6.7 *μ*s (which corresponds to the time that the mean flow passes through the whole domain about 7.5 times) to attain a steady state of the flow during which the graphene is sustained to be rigid. After this time, the graphene is released to freely move under the action of the flow field. Four cases of the mass ratio, *χ*=0,0.45,0.5,1, and two cases of the coefficient of restitution *e*=0.5,1 are tested. *χ*=0.45 is the theoretical value of the mass ratio between the air and graphene in this simulation. Since *χ* and *e* are chosen independently, there are eight combinations of *χ* and *e*.

Figure 2 shows the time series of the streamwise hydrodynamic force (drag) acting on the graphene for several mass ratios *χ* with perfect reflection (*e*=1) and inelastic reflection (*e*=0.5). Figure 3 shows the drag and wall-normal hydrodynamic forces (lift), and their oscillation amplitudes at the final time of the oscillation (approximate equilibrium state). When *χ* increases, the average value and the oscillation amplitude at the early times of the drag acting on the graphene decrease, while the approach to the asymptotic equilibrium state becomes slightly slower. When *e* is decreased, the same behaviour of the average value and the oscillation amplitude of the drag is observed. In addition, the reduction of the equilibrium value of the lift force acting on the graphene is observed for an increase of *χ* and/or for a decrease of *e*.

These observation can be understood as follows. Suppose that there is no fluid. Then, the initially perturbed graphene keeps oscillating because the dynamics inside the graphene is described by the coarse-grained Hamiltonian for which the total energy is conserved [3]. Now, switch on the fluid motion. From equation (2.9), the change of fluid momentum owing to the collision is computed as
3.1
This equation shows that the momentum change decreases with increase in *χ* and with decrease in *e*, provided that *u*_{f}−*u*_{s} is unchanged, which implies that the graphene experiences less resistive force from the fluid. Therefore, the oscillation amplitude at early times and the equilibrium value at later times become smaller with an increase in *χ* and/or with a decrease in *e*. The oscillation period is almost unchanged, irrespective of the change in *χ* and *e*. This is because the period is determined by the graphene’s property alone.

Figure 4 shows the time series of the streamwise displacement of the tip of the graphene for several *χ* for figure 4*a* (*e*=1) and figure 4*b* (*e*=0.5), respectively, and figure 5 shows the final position of the tip in the streamwise and wall-normal directions (symbols) and its oscillation amplitude. Note that the amplitude of oscillation in the streamwise direction is 10 times as large as the wall-normal one. Initially, the graphene bends downstream, then oscillates with decreasing amplitude in time. The displacement becomes the largest when *χ*=0, while the amplitude of oscillation becomes larger when *χ*=1. Decrease of *e* causes a reduction in the displacement and an increase in the amplitude of oscillation.

## 4. Conclusions

We have developed the modified interface model that takes into account the mass ratio and coefficient of restitution between a fluid and solid. The model is applied to the numerical computation of motion of the graphene in the laminar flow by using the immersed boundary method. It was observed that the streamwise hydrodynamic force acting on the graphene decreases when the mass ratio *χ* increases and the coefficient of restitution *e* decreases. On the other hand, the oscillation amplitude of the force increases and the oscillation period is unchanged. The same trend was also observed for the displacement of the graphene in the downstream direction. These observations were explained in terms of the momentum exchange formula equation (3.1), the smaller the drag force becomes owing to the increase of *χ*, the more freely the graphene is allowed to oscillate. It is important for more accurate computation of the fluid–solid interaction to study these effects theoretically and numerically with the help of MD. Comparison of the present model with another numerical model using MD simulation is now under way, and the results will be reported somewhere.

## Acknowledgements

This work was supported by the Japan Science and Technology Agency. The authors thank Prof. Ogata and Dr Kobayashi for their fruitful discussions, and the Information Technology Center at Nagoya University for computational support.

## Footnotes

One contribution of 25 to a Theme Issue ‘Discrete simulation of fluid dynamics: applications’.

- This journal is © 2011 The Royal Society