Computing hydrodynamic interactions in conned doubly-periodic geometries in linear time Aref Hashemi1Ra ul P. Pel aez1 2Sachin Natesh1 3Brennan

2025-04-27 0 0 4.03MB 36 页 10玖币
侵权投诉
Computing hydrodynamic interactions in confined
doubly-periodic geometries in linear time
Aref Hashemi,1, ∗Ra´ul P. Pel´aez,1, 2 Sachin Natesh,1, 3 Brennan
Sprinkle,1, 4, †Ondrej Maxian,1Zecheng Gan,1, 5 and Aleksandar Donev1, ‡
1Courant Institute, New York University, New York, NY, US
2Departamento F´ısica Te´orica de la Materia Condensada,
Universidad Aut´onoma de Madrid, Madrid, Spain
3Department of Applied Mathematics, University of Colorado Boulder, Boulder, CO, US
4Department of Applied Mathematics and Statistics, Colorado School of Mines, Golden, CO, US
5Department of Mathematics, The Hong Kong University of Science and Technology, Hong Kong SAR, China
We develop a linearly-scaling variant of the Force Coupling Method [K. Yeo and M. R.
Maxey, J. Fluid Mech. 649, 205–231 (2010)] for computing hydrodynamic interactions among
particles confined to a doubly-periodic geometry with either a single bottom wall or two
walls (slit channel) in the aperiodic direction. Our spectrally-accurate Stokes solver uses the
Fast Fourier Transform (FFT) in the periodic xy plane and Chebyshev polynomials in the
aperiodic zdirection normal to the wall(s). We decompose the problem into two problems.
The first is a doubly-periodic subproblem in the presence of particles (source terms) with
free-space boundary conditions in the zdirection, which we solve by borrowing ideas from a
recent method for rapid evaluation of electrostatic interactions in doubly-periodic geometries
[O. Maxian, R. P. Pel´aez, L. Greengard and A. Donev, J. Chem. Phys. 154, 204107 (2021)].
The second is a correction subproblem to impose the boundary conditions on the wall(s).
Instead of the traditional Gaussian kernel, we use the exponential of a semicircle kernel to
model the source terms (body force) due to the presence of particles, and provide optimum
values for the kernel parameters that ensure a given hydrodynamic radius with at least two
digits of accuracy and rotational and translational invariance. The computation time of our
solver, which is implemented in graphical processing units, scales linearly with the number
of particles, and allows computations with about a million particles in less than a second for
a sedimented layer of colloidal microrollers. We find that in a slit channel, a driven dense
suspension of microrollers maintains the same two-layer structure as above a single wall, but
moves at a substantially lower collective speed due to increased confinement.
I. INTRODUCTION
The development of more efficient, accurate, and scalable methods for suspensions of rigid and
flexible particles in Stokes flow remains a key challenge in soft condensed matter physics and
chemical engineering. In addition to the long-ranged nature of hydrodynamic interactions, the
inclusion of Brownian motion and the presence of confining boundaries pose particular difficulties.
While confining boundaries partially screen the hydrodynamic interactions, they continue to decay
algebraically rather than exponentially [1], and must be captured to resolve particle dynamics.
A key component of all computational methods for Stokes flow is the efficient evaluation of
the action of the singular or regularized Green’s function for Stokes flow among a large number of
particles. Specifically, given forces (and sometimes also torques) on many source points, the goal
is to compute the resulting (linear and sometimes also angular) velocities on many target points,
usually the same as the source points. While boundary integral methods typically use the singular
∗aref@cims.nyu.edu
†bsprinkl@mines.edu
‡corresponding author; donev@courant.nyu.edu
arXiv:2210.01837v4 [physics.comp-ph] 3 Apr 2023
2
Green’s function, most computational methods used for large-scale suspensions are based on a
regularized Green’s function. There are four popular regularizations: the Rotne-Prager-Yamakawa
(RPY) [2,3] tensor, regularized Stokeslets [4], immersed boundary (IB) kernels [5], and the Force
Coupling Method (FCM) kernel [6–10]. The RPY, IB, and FCM kernels regularize (smooth) the
singularity at both the source and target; this is crucial to maintain the symmetric positive definite
(SPD) nature of the hydrodynamic mobility matrix, as is necessary for Brownian Dynamics (BD)
methods.
Here, we develop a fast method for Stokes suspensions in doubly-periodic (DP) geometries with
one or two confining walls in the aperiodic direction. Our method is closely related to the Force
Coupling Method of Yeo and Maxey [8], but with several important differences that increase both
the flexibility and the efficiency of the approach. In brief, we employ a non-Gaussian envelope
function [11] that allows accurate computations with far fewer grid cells per particles, as in im-
mersed boundary methods (IBM) [12], and also develop a novel spectral Stokes solver based on
recent work by some of us on fast methods for electrostatics in doubly-periodic geometries [13].
The hydrodynamic interaction between two distinct spherical particles or blobs of radius aat
positions r(1) and r(2) can be captured by a 3 ×3 mobility tensor that gives the velocity of one
of the particles for a given force acting on the other. This hydrodynamic mobility tensor can be
approximated in the far field as
Or(1),r(2)=Zδar(1) −r0Gr0,r00δar(2) −r00dr0dr00,(1)
where Gis the Green’s function for Stokes flow with the specified boundary conditions. Here
δa(r) is a regularized “delta function” or envelope function [6] that is typically radially-symmetric,
δa(r)≡δa(r), when the particles are sufficiently far from boundaries. It is important to point out
that the blobs do not have to represent actual physical spherical colloids; one can build physical
particles as a collection of blobs, including non-spherical rigid (passive or active) particles [14–18],
which we have termed as rigid multiblobs, or flexible particles such as semiflexible fibers [19,20].
If the blobs represent actual spherical colloids, additional near-field lubrication corrections can be
added to improve upon the far-field accuracy [10,15,21–23]; this requires including torques in
addition to forces, as we do in Sec. II.
Regardless of the context, a key task is to evaluate for Nblobs their far-field velocity Ufrom
the applied forces Fthrough the action of the mobility matrix M,U=MF,
∀i:U(i)=
N
X
j=1
Or(i),r(j)F(j),(2)
in time linear in the number of blobs N. Note that the form of (1) guarantees that the mobility
matrix is SPD by construction since Gis an SPD kernel, and the regularization is applied both at
the source and the target.
The specific choice of δaas a delta function on the surface of a sphere of radius a,δa(r) =
4πa2−1δ(r−a), results in the widely-used Rotne-Prager-Yamakawa (RPY) tensor [24–26].
When the kernels do not overlap the boundaries, one can transform (1) into a differential form by
employing the Fax´en differential operator,
Or(1),r(2)≈I+a2
6∇2
r0I+a2
6∇2
r00 G(r0,r00)
r0=r(1)
r00=r(2) .(3)
This analytical simplification allows for explicit evaluation of the RPY kernel, not just in an
unbounded domain but also in a half-space above a no-slip wall [24], because for a single no-slip
3
boundary there is a relatively simple image construction for Gdue to Blake [27]. This image
construction has also enabled the development of fast methods for evaluating (2) for a single wall,
based on either (flexible periodicity) Fast Multipole Method [28,29] or the Fast Fourier Transform
(FFT) [30,31]. For triply-periodic (TP) domains, the Positively Split Ewald (PSE) method [26]
for evaluating (2) (and also generating Brownian velocities) for the RPY kernel provides the basis
for Fast Stokesian Dynamics [15]. However, it remains a challenge to construct a similarly-efficient
method for confined suspensions since PSE uses Fourier representations in a key way for each
component of the method (Stokes solver, generating Brownian increments, and Ewald splitting).
In principle, a greater flexibility with respect to boundary conditions can be achieved by re-
placing analytical Green’s functions with a grid-based Stokes solver [32]. This requires replacing
the singular surface delta function form of δa(r) with a smooth function that can be resolved on a
grid. In the FCM, a Gaussian envelope function δais used, which allows for analytical calculation
of Oin an unbounded domain [6,7]; in numerical methods the Gaussian is truncated. By con-
trast, in the IBM, δais a discrete grid function specifically constructed to maximize grid invariance
[5,12,17], and the double convolution in (1) is discrete. This makes all IBM results grid- and
solver-dependent, without a direct continuum limit.
We develop a method that combines favorable features of FCM and IBM. Namely, we maintain
the continuum representation (1) from FCM, however, we do not use a Gaussian envelope but
rather use the “exponential of a semi-circle” (ES) kernel proposed by Barnett for its efficiency
in the context of non-uniform FFTs [11]. Not only does this kernel allow for greater flexibility in
tuning the hydrodynamic radius of the particles [33], but also allows us to use many fewer grid cells
per particle than for a Gaussian kernel, comparable to the IBM, while still solving the continuum
equations to several digits of accuracy. Specifically, for about two-three digits of accuracy we need
four grid cells per particle (in each dimension) if only translational velocity is required, and five
or six if rotational velocities are also required. The grid independence of the results makes them
transferable, and allows us to separate the Stokes solver from the physics (i.e., the problem has a
solution without a discretization). We find that Ocomputed with the ES and Gaussian kernels
are indistinguishable in practice. We employ an image construction following Yeo and Maxey [8]
(also used in the IBM [17,34]) to generalize (1) to the relevant-in-practice case when some of the
blob kernels overlap the boundaries; this is considerably harder to do for the RPY kernel, and has
required ad hoc fixes in past work by some of us [35].
For TP domains one can easily solve the Stokes equations spectrally in Fourier space, but this
is considerably harder when boundaries are present. In Sec. III we use some of the ideas applied to
the Poisson equation in [13] to develop a Stokes solver for DP geometries that are unbounded or
confined by one or two boundaries in the aperiodic direction. Specifically, we focus on bottom wall
(BW) DP geometries, with a single wall at z= 0, and slit channel (SC) DP geometry, confined by
two walls at z= 0 and z=H, but the Stokes solver can handle more flexible boundary conditions
in the aperiodic direction. Our novel fluid solver makes it possible to handle geometries that are
partially unbounded in one direction, unlike in existing IBM or FCM implementations based on
more traditional grid solvers like finite differences or finite elements [36]. Our fluid solver has the
additional advantage that its implementation requires only calls to the three-dimensional (3D) FFT
with an oversampling factor of 2 in the aperiodic direction, combined with trivially parallelizable
one-dimensional pentadiagonal boundary value solvers in the zdirection. We use this to implement
the method on Graphical Processing Units, achieving linear scaling up to as many as one million
particles, with scaling constants much better than existing fast methods [29] for the types of
problems we study here. It is important to note that Srinivasan and Tornberg have recently
developed a sophisticated method based on 3D FFTs for the singular Stokes Green’s function with
a bottom wall geometry (using images) and flexible periodicity [30,31]. Our approach is different
and specialized to doubly-periodic geometries, is conceptually simpler, and affords flexibility in the
4
boundary conditions in the unbounded direction. It is beyond the scope of this work to give a
thorough comparison of the different approaches.
In Sec. III we also study how to choose the parameters of the ES kernel to achieve a desired
accuracy and effective blob hydrodynamic radius with the coarsest possible solver grid. In Sec. IV
we perform a number of validation tests examining the self and pair mobility of particles in BW and
SC geometries, and show that the method produces results in agreement with existing theoretical
or numerical predictions. In Sec. Vwe demonstrate that it is possible to generate stochastic
(Brownian) particle displacements with covariance proportional to Musing a Lanczos method
[37,38] in a modest number of iterations independent of the number of particles, for both BW and
SC geometries. This allows us to replace the core hydrodynamic routines used in the Stokesian
dynamics method developed in [22] and the rigid multiblob BD methods developed in [18] with new
linear-scaling implementations that are substantially more efficient for sufficiently large number of
particles.
These new computational developments enable us to perform larger-scale studies of the dynamics
of confined microroller suspensions than previously feasible, in Sec. V. After validating that the
results presented by some of us in [22] are free of finite-size effects, we explore what a top wall does
to a driven dense suspension of colloidal microrollers. We find that while the suspension maintains
a two-layer structure as above a single bottom wall, in a slit channel the collective velocity is
substantially reduced due to the increased confinement, as expected.
II. MODEL FORMULATION
We develop a method to solve the Stokes equations,
η∇2u−∇p=−f,(4)
∇·u= 0,(5)
in an xy doubly-periodic domain of size x y∈[−Lx, Lx]×[−Ly, Ly] [39], with z∈[0,∞); we
refer to this as a bottom wall (BW) geometry. We will return later to the case of a slit channel
(SC) geometry with two walls at z= 0 and z=H, for which z∈[0, H]. Here u=u v w|is
the velocity field of the fluid, pthe pressure, ηthe viscosity, and the body force f=f g h|
represents particles. To close the problem, we impose tangential slip boundary conditions (BCs)
at z= 0, i.e.,
u|z=0 =uwall =uwall vwall 0|,(6)
where uwall(x, y) is a smooth function giving a prescribed slip velocity (uwall =0for no-slip
boundaries) along the wall (e.g., electrophoretic slip), and we require boundedness of uas z→ ∞.
For all of the results and tests presented in this study, we use uwall =0.
The source term fis restricted to the domain limits (i.e. f(z /∈[0, H]) = 0) and consists of
regularized monopoles and dipoles with envelope functions ∆Mand ∆Dcentered around particles
or blobs, i.e.,
f(x) =
N
X
j=1 F(j)∆M(x−r(j)) + 1
2∇×(τ(j)∆D(x−r(j))),(7)
where x=x y z|and r(j)=x(j)y(j)z(j)|is the jth particle’s location with j= 1, . . . , N,
and F(j)and τ(j)are the force and torque of the jth particle. Here, ∆Mand ∆Dare compactly-
supported, smooth, regularized delta functions called kernels in the immersed boundary method
5
(IBM), or envelopes in the Fast Coupling Method (FCM). In the classical FCM, the envelopes are
(truncated) Gaussians,
∆M/D(x) = 1
q8π3gM/D6exp −kxk2
2g2
M/D !,(8)
where gM=Rh/√πand gD=Rh/(6√π)1
3, with Rhas the effective hydrodynamic radius of a
particle/blob.
If a particle is closer to the wall than a distance zim so that its kernel overlaps the wall, we add
the negative of the particle’s envelope centered at its image point about the wall. This ensures
that the force/torque decays to zero as a particle approaches the wall. As a result, if z(j)< zim,
we replace the envelopes ∆Mand ∆Din (7) with
∆W
M/D x−r(j)= ∆M/D x−r(j)−∆M/D x−r(j)
im ,(9)
where r(j)
im =r(j)−2ˆ
ez(ˆ
ez·r(j)) is the particle’s point of reflection about the bottom wall, and
ˆ
ez=001|. In the FCM, Yeo and Maxey [8] put a positive sign in (9) for ∆D, to ensure that
the angular velocity of a particle in simple shear flow is unaffected by the proximity to the wall.
We put a negative sign to ensure that all components of the particle mobility go to zero at the wall.
This is a modeling choice, and it is important to realize that (9) is not an exact image construction
for a no-slip wall [30]. For a free-slip boundary such as a gas-liquid interface, it is straightforward
to make an exact image construction; see, for example, Appendix A in [40].
Our goal is to compute the motion of each of the particles for dynamical simulations. Suppose
we’ve solved the Stokes equations for the fluid velocity u; then the linear and angular velocities
of the particles, denoted as U(j)and Ω(j), can be obtained through the following volume integrals
over the envelopes,
U(j)=Zx
u(x)∆M(x−r(j))dx,(10)
Ω(j)=1
2Zx
(∇×u(x)) ∆D(x−r(j))dx,(11)
with the understanding that ∆M/D get replaced by ∆W
M/D for particles close to a wall.
In our own variant of FCM, we will not use a Gaussian kernel; rather, we will employ a
compactly-supported Exponential of a Semicircle (ES) kernel [11] optimized for numerical effi-
ciency in spectral methods, as we discuss in detail in Sec. III E. Other kernels can be used as well,
as long as their Fourier transform decays sufficiently fast in the wavenumber.
III. STOKES PROBLEM
We develop a doubly-periodic+correction Stokes solver based on the observation that, the orig-
inal problem (4)–(5) can be separated into two subproblems, namely the doubly-periodic “DP”
problem and the “correction” problem. For the DP subproblem, we keep the source terms fbut
there is no wall present (i.e., open BCs in z). The correction subproblem, has the walls but has no
forcing, so it can be solved analytically using a plane-wave expansion. The total solution is then
the sum of the solutions to the two subproblems.
Below, we detail our solution method for the bottom wall (BW) geometry with a no-slip wall.
Appendix Asummarizes the method for the slit channel (SC) geometry with somewhat more
摘要：

Computinghydrodynamicinteractionsincon neddoubly-periodicgeometriesinlineartimeArefHashemi,1,RaulP.Pelaez,1,2SachinNatesh,1,3BrennanSprinkle,1,4,yOndrejMaxian,1ZechengGan,1,5andAleksandarDonev1,z1CourantInstitute,NewYorkUniversity,NewYork,NY,US2DepartamentoFsicaTeoricadelaMateriaCondensada,Uni...

展开>> 收起<<
Computing hydrodynamic interactions in conned doubly-periodic geometries in linear time Aref Hashemi1Ra ul P. Pel aez1 2Sachin Natesh1 3Brennan.pdf

共36页,预览5页

还剩页未读， 继续阅读

声明：本站为文档C2C交易模式，即用户上传的文档直接被用户下载，本站只是中间服务平台，本站所有文档下载所得的收益归上传人(含作者)所有。玖贝云文库仅提供信息存储空间，仅对用户上传内容的表现方式做保护处理，对上载内容本身不做任何修改或编辑。若文档所含内容侵犯了您的版权或隐私，请立即通知玖贝云文库，我们立即给予删除！

相关推荐

更多
立即下载
分类：图书资源 价格：10玖币 属性：36 页 大小：4.03MB 格式：PDF 时间：2025-04-27

开通VIP享超值会员特权

  • 多端同步记录
  • 高速下载文档
  • 免费文档工具
  • 分享文档赚钱
  • 每日登录抽奖
  • 优质衍生服务

作者详情

相关内容

更多

热门标签

高考理综 人际关系 力的合成 配电装置 动力学连接体 全宋诗作者索引 公务员考试
/ 36
评分 收藏
分享
立即下载
客服
关注