Lagrangian manners

For the immersed interface X, we use a collection of discrete points sk = k∆s, k = 0, 1, · · · , M , with the interface mesh width ∆s so the Lagrangian markers of the interface are represented by X_k = X(s_k). Other interfacial quantities, such as force density F and velocity U are defined on position Xk so the notations Fk = F(sk) and U_k = U(s_k) are used. One can see the illustration of IB framework of regular uniform grid and set of discrete Lagrangian markers in Figure 3.3.

Figure 3.3: The basic framework of IB method, spreading discrete Lagrangian markers in a computational domain with uniform mesh.

Especially pay attention that we often have to compute derivatives of interfacial quantities. For instance, the vesicle membrane forces which consist of tension and bending force involve evaluating derivatives of position Xk. Therefore calculation of derivative along the interface plays an important role to the dynamics of the interface.

Here we provide two ways to compute derivative of a function along interface. Without loss of generality, for any function defined on the interface ψ(s), we approximate the partial derivative ^∂ψ_∂s by

• Finite difference method: The value of ^∂ψ_∂s can be approached by standard second-order central difference scheme as

Dfψ(s) = ψ(s + ∆s/2) − ψ(s − ∆s/2)

∆s .

• Spectral method: Since ψ is assumed to be periodic along interface, we can adopt the spectral Fourier discretization to achieve higher-order of accuracy.

The interface is represented by the discrete Fourier series expansion as

ψ(s) =

M/2−1X

k=−M/2

cψ_ke^iks,

where cψk is the corresponding Fourier coefficients for ψ(s). This Fourier coef-ficients can be performed very efficiently by using the Fast Fourier Transform (FFT). Then the Fourier coefficients of P -th derivative of ψ(s) can be obtained by taking an inverse FFT to (ik)^Pcψ_k.

Chapter 4

Fractional step immersed boundary method

In this chapter, we develop a fractional step method based on IB formulation for Stokes flow with an inextensible interface (vesicle without bending stiffness). Solving the fluid variables such as the velocity and pressure, the present problem involves finding extra unknown elastic tension such that the surface divergence of the velocity is zero along the interface. Once the velocity is found, the interface is moving accord-ing to the local fluid velocity. Since the fluid is incompressible and the interface is inextensible, both the area enclosed by the interface and the total length of the in-terface should be conserved. Notice that, the dynamics of vesicles are determined by their boundary rigidity, inextensibility, and the hydrodynamical forces. Our present model for inextensible interface can be regarded as a simplified vesicle model without bending effect.

The rest of the chapter are organized as follows. In the next section, we describe our governing equations for the Stokes flow with an inextensible interface based on IB formulation. We also show the skew-adjoint property between the spreading oper-ator acting on the tension and the surface divergence of the velocity. In Section 4.2, the symmetry of the resultant matrix equation is provided first, and then a numerical algorithm based on the fractional step method is developed. Finally, the numerical re-sults including the convergence tests and the tank-treading motion for an inextensible interface under a simple shear flow are shown in detail in Section 4.4.

4.1 Governing equations

We begin by stating the mathematical formulation of the Stokes flow with an inex-tensible interface. Consider there is a moving, immersed, inexinex-tensible interface Γ(t) in the fixed fluid domain Ω. We assume that the fluids inside and outside of the interface are the same so the governing equations in IB formulation can be written

as follows.

Eqs. (4.1) and (4.2) are the familiar incompressible Stokes equations with a singular force term arising from the interface. Eq. (4.3) represents the inextensibility constraint of the interface which is equivalent to the zero surface divergence of the velocity along the interface. Here, the interfacial velocity U is simply the interpolation of the fluid velocity at the interface which is defined as in traditional IB formulation. Eq. (4.4) simply represents that the interface moves along with the local fluid velocity (the interfacial velocity). The interaction between the fluid and the interface is linked by the two-dimensional Dirac delta function δ(x) = δ(x)δ(y). Unlike the traditional IB formulation in which the elastic tension σ(s, t) is either known or a function of immersed boundary configurations, here, the tension is a part of solution needed to be determined. In this model, we consider the Stokes flow; however, the numeri-cal method can be extended to Navier-Stokes flow straightforwardly by treating the nonlinear advection terms explicitly during the time evolution.

The difficulties in solving the above interfacial problem are as follows. Firstly, since the fluid is incompressible and the interface is inextensible, both the area enclosed by the interface and total length of the interface should be conserved simultaneously.

Furthermore, the local inextensibility constraint (4.3) is more stringent than the con-servation of the total interfacial length since the latter one is a global constraint.

Secondly, the elastic tension must be treated as an unknown function which is needed to be solved with the fluid variables simultaneously. In previous literature, most of related work is based on boundary integral methods, see for example, [71, 66, 59] and the references therein. However, boundary integral methods generally assume infi-nite domains, and cannot be generalized to full Navier-Stokes equations since there is no corresponding Green function. Until recently, Kim and Lai [28] have applied a penalty IB method to simulate the dynamics of inextensible vesicles. By introduc-ing two different kind of Lagrangian markers, the authors are able to decouple the fluid and vesicle dynamics so that the computation can be performed more efficiently.

One potential problem of this approach is that the time step depends on the penalty number and must be chosen smaller as the penalty number becomes larger. In [35], a new finite difference scheme based on immersed interface method (IIM) has been developed for solving the present problem in Navier-Stokes flow. The authors treat the unknown elastic tension as an augmented variable so that the augmented IIM can

be applied. In this chapter, we discretize the equations (4.1)-(4.3) directly without decoupling and use a fractional step method to solve the resultant linear system of equations. The numerical algorithm will be given in next section.

Before to continue, we first prove that spreading operator acting on the function σ and the surface divergence operator of the velocity are skew-adjoint with each other.

To proceed, let us define the spreading operator S of σ and the surface divergence operator ∇_s of U as follows.

S(σ) = We also define the inner product of functions on Ω and Γ in the following.

hu, viΩ=

where l in Eq. (4.7) is the arc-length parameter. Then we have hS(σ), uiΩ =

ds (intergartion by parts and the closed interface)

=

By comparison, it leads that the spreading operator and the surface divergence oper-ator are skew-adjoint.

The reason for showing the skew-adjointness of those two operators are two fold.

Firstly, since the surface divergence of velocity is zero in Eq. (4.8), from the above derivation, we have hS(σ), ui_Ω = 0. That is, the present elastic tension does not do work to the fluid which is not surprising since it is the Lagrange multiplier for the inextensible constraint. However, if we add the bending force along the interface as the one in vesicle problems [66, 28], then the bending force does do work to the fluid.

Secondly, the skew-adjointness is also satisfied in discrete sense (see in next section) so that the resultant matrix equation is symmetric; therefore, some efficient iterative solvers such as conjugate gradient (CG) method can be easily applied.