PaperlayerView exact revision on arXiv

Source-linked AI summary

Efficient and accurate simulations of deformable particles immersed in a fluid using a combined immersed boundary lattice Boltzmann finite element method

Timm Krüger, Fathollah Varnik, Dierk Raabe

arXiv:1004.2416v1cond-mat.soft

TL;DR

The paper asks how to make coupled deformable-capsule simulations accurate and efficient at small resolutions. It evaluates an immersed-boundary lattice Boltzmann finite element method against analytic small-deformation results and finds hydrodynamic discretization more influential than membrane-mesh details. The analysis identifies effective capsule properties and resolution choices for reducing cost while maintaining accuracy.

  • Problem

    Accurate simulations of multiple deformable particles require efficient methods that capture dominant physics at small spatial resolutions.

  • Method

    The study couples lattice Boltzmann fluid dynamics, immersed-boundary coupling, and finite elements with a Skalak membrane model, testing mesh, spatial-resolution, and interpolation-stencil effects.

  • Results

    Hydrodynamic resolution and the IBM interpolation stencil strongly affect accuracy, whereas membrane tessellation and resolution have comparatively minor influence.

  • Takeaways & Limitations

    Accuracy is best improved first by increasing fluid-lattice resolution, while computational cost can be reduced by lowering mesh resolution at fixed hydrodynamic resolution.

Abstract

from arXiv · show

The deformation of an initially spherical capsule, freely suspended in simple shear flow, can be computed analytically in the limit of small deformations [D. Barthes-Biesel, J. M. Rallison, The Time-Dependent Deformation of a Capsule Freely Suspended in a Linear Shear Flow, J. Fluid Mech. 113 (1981) 251-267]. Those analytic approximations are used to study the influence of the mesh tessellation method, the spatial resolution, and the discrete delta function of the immersed boundary method on the numerical results obtained by a coupled immersed boundary lattice Boltzmann finite element method. For the description of the capsule membrane, a finite element method and the Skalak constitutive model [R. Skalak et al., Strain Energy Function of Red Blood Cell Membranes, Biophys. J. 13 (1973) 245-264] have been employed. Our primary goal is the investigation of the presented model for small resolutions to provide a sound basis for efficient but accurate simulations of multiple deformable particles immersed in a fluid. We come to the conclusion that details of the membrane mesh, as tessellation method and resolution, play only a minor role. The hydrodynamic resolution, i.e., the width of the discrete delta function, can significantly influence the accuracy of the simulations. The discretization of the delta function introduces an artificial length scale, which effectively changes the radius and the deformability of the capsule. We discuss possibilities of reducing the computing time of simulations of deformable objects immersed in a fluid while maintaining high accuracy.

I. INTRODUCTION

The paper targets efficient, accurate simulations of deformable capsules at small resolutions, motivated by multiscale blood-flow applications. It studies how coupling, membrane discretization, and mesh construction affect simulation behavior.

  • Simulating blood microcirculation requires resolving hydrodynamics, membrane dynamics, cell interactions, vessel-wall interactions, and scale-bridging behavior.
  • Analytic small-deformation solutions for capsules in shear flow provide reference results for evaluating numerical approaches.
  • The explicit IBM coupling is computationally efficient, but may imperfectly enforce no-slip and permit capsule-volume drift.For soft capsules and short simulations, the authors report negligible volume drift and no instabilities.
  • The paper does not address IBM stability issues or improved no-slip implementations.
  • The authors use an unstructured membrane mesh, noting that structured meshes contain coordinate singularities while unstructured meshes require gradient approximation.
  • The study examines the combined IBLBFEM and membrane tessellation for accuracy and numerical efficiency at small resolutions.

II. NUMERICAL METHODS

The simulation algorithm combines lattice Boltzmann fluid dynamics, immersed-boundary fluid–membrane coupling, and finite elements for membrane deformation.

  • The fluid solver uses the D3Q19 Bhatnagar-Gross-Krook lattice Boltzmann method.
  • Finite elements evaluate membrane strains derived from the membrane constitutive model.
  • The immersed boundary method captures interactions between the fluid and the membrane.
  • The particle mesh is included as a separate simulation property because it influences the simulations.

A. Lattice Boltzmann method

The lattice Boltzmann method represents fluid motion through populations streaming and colliding on a regular lattice. Here it uses a D3Q19 model, body-force coupling, moving-wall shear flow, and low-Mach, low-Reynolds assumptions.

  • LBM evolves mesoscopic populations f_i that stream and collide on a regular lattice instead of directly discretizing pressure and velocity.
  • The BGK relaxation parameter τ controls kinematic viscosity and can also change the simulation time step.
  • The D3Q19 model uses 19 discrete velocities whose lattice weights ensure fluid isotropy.
  • A body-force density is incorporated into the lattice Boltzmann equation for coupling with immersed membranes.
  • Macroscopic density, velocity, and deviatoric shear stress are recovered locally from population moments.
  • Moving-wall bounce-back generates simple shear flow between walls at z = ±H/2, with shear rate γ̇ = 2u_w/H.
  • The simulations require sufficiently small wall speeds for the small-Mach-number limit and Re ≪1 for the small-deformation capsule theory.

B. Membrane model and force computation

The membrane model computes elastic forces from triangular finite elements using strain invariants and the Skalak constitutive law. The simulations include shear energy but omit bending, volume, and surface energies.

  • Membrane forces are derived from strain-energy contributions, with total energy decomposed into strain, bending, surface, and volume terms.
  • The Skalak model describes isotropic membrane elasticity through principal stretches and strain invariants.
  • The shear modulus k_s and area-dilation modulus k_α control membrane resistance to shear and dilation.
  • The membrane is discretized into flat triangular face elements that remain flat even under large deformations.
  • Linear shape functions map node displacements to a spatially constant displacement gradient within each triangular face.
  • Face deformation is determined from current and undeformed edge lengths and included angle, yielding the principal stretches λ_1 and λ_2.
  • The simulations retain only shear energy, setting W = W^S while omitting bending, volume, and surface energies.
  • Membrane forces acting on the fluid are computed from the membrane energy using the principle of virtual work.

C. Immersed boundary method

The immersed boundary method couples membrane forces and motion to a fixed fluid lattice through force spreading and velocity interpolation. Each time step alternates finite-element force calculation, lattice-Boltzmann fluid update, and membrane advection.

  • The immersed boundary method couples an arbitrary Lagrangian membrane surface to the Eulerian fluid lattice as a front-tracking method.
  • Membrane forces are converted into an Eulerian body-force density through the discrete spreading operation.The force density is then supplied to the lattice Boltzmann equation.
  • The discrete delta kernel has finite support, and its width is a free immersed-boundary parameter subject to conservation requirements.The method uses a separable kernel δ(r) = φ(x)φ(y)φ(z), with common interpolation functions of different support widths.
  • Membrane-node velocities are interpolated from the updated lattice velocities, assuming the membrane moves with the ambient fluid velocity under the no-slip condition.The new velocities use new lattice velocities and old membrane-node positions before explicit advection.
  • The discretized explicit method does not guarantee exact volume conservation, but volume-drift correction is unnecessary for the paper’s short simulations.
  • Each time step computes membrane forces, spreads them to the lattice, updates the fluid, interpolates membrane velocities, and advances membrane positions.The scheme sets ∆t = 1 and repeats after optionally writing the membrane and fluid states to disk.

D. Membrane tessellation

The paper compares three approaches for generating spherical membrane meshes and examines their geometric regularity and construction constraints. Successive subdivision of a highly symmetric base mesh produces the most homogeneous and isotropic meshes, but restricts available face counts.

  • Mesh-generation approaches: The study considers implicit-surface tessellation with CGAL, finite-element mesh generation with Gmsh, and successive subdivision from a coarse symmetric mesh.
  • CGAL: CGAL can tessellate arbitrary implicit surfaces while controlling face count and node spacing, but its meshes have reduced homogeneity and isotropy.
  • Gmsh: Gmsh tessellates constructed geometric objects such as spheres, although the mesh does not completely preserve the initial sphere’s high isotropy.
  • Successive subdivision: An icosahedral subdivision starts from 20 equilateral triangular faces, 12 nodes, and 30 equal-length edges, then creates midpoint nodes and projects them onto the sphere.For triangular closed surfaces, node and face counts satisfy 2N_n = N_f + 4.
  • Subdivision limitation: Subdivision produces strong geometric regularity but cannot generate an arbitrary number of faces from a fixed base polyhedron.The restriction can be relaxed by choosing a lower-symmetry starting body with a different face count.
  • Mesh comparison: At approximately 1280 faces, the subdivided mesh has the smallest scatter in face area, edge length, and angle distributions among the compared meshes.The icosahedron-based mesh is used by default unless otherwise stated.

III. THEORY

The theory defines capsule deformation, inclination, and tank-treading measures for initially spherical capsules in simple shear flow, using small-deformation analytical results and geometric characterization methods.

  • Equal internal and external fluid density and viscosity simplify the capsule-flow computations without losing too much generality.
  • The steady capsule is characterized by deformation D, inclination angle θ, and angular velocity ω during tank-treading motion.
  • The capsule inclination angle is measured between the largest semiaxis and the flow direction, while membrane angular velocity is below γ̇/2 for deformable capsules.
  • For small deformations in Stokes flow, the Taylor deformation parameter D quantifies capsule shape change, with D = 0 for a sphere.
  • The inertia-tensor and linear-fit methods produce virtually identical D and θ, although the inertia tensor slightly underestimates semiaxes because of mesh discretization and flat elements.

IV. SIMULATIONS AND RESULTS

The simulations fix the Skalak model and baseline lattice-Boltzmann parameters, then examine how resolution and numerical parameters affect capsule deformation and IBLBFEM convergence.

  • The simulations use the Skalak model with ks = kα, dimensionless lattice spacing, and default lattice-Boltzmann relaxation parameter τ = 1.
  • The study varies simulation-box size, Reynolds number, reduced shear rate, relaxation parameter, membrane tessellation, and immersed-boundary interpolation stencils.

A. General simulation parameters

The study identifies parameter ranges that approximate unbounded, low-Reynolds-number linear shear flow and examines deviations from the analytic deformation response. Finite size and inertia are limited under selected conditions, whereas relaxation, nonlinear response, and IBM discretization can affect accuracy.

  • H/r = 10 is sufficient to model unbounded simple shear flow, with plateau differences below 0.5% compared with H/r = 12.Finite-box effects are negligible for H/r ≥10.
  • Re ≤0.02 keeps inertia effects unimportant for the steady-state deformation response.Differences between Re = 0.01 and 0.02 appear only slightly in the transient evolution.
  • G = 0.01 preserves the linear small-deformation approximation, whereas deviations become significant for G ≥0.02.The nonlinear membrane response becomes evident at higher reduced shear rates.
  • τ ≈1 is a practical compromise because larger relaxation parameters reduce time steps but can produce drastic numerical artifacts.Results are relatively independent of τ for τ ≤1 and diverge when τ becomes substantially larger than unity.
  • H/r = 10, Re = 0.02, G = 0.01, and τ = 1 provide excellent approximations to unbounded Couette flow in the linear elastic limit.Under these settings, the theoretical Skalak-membrane value is D/G = 6.25, but simulated deformation parameters are larger.
  • The excess deformation is attributed most plausibly to spatial discretization or discrete-delta interpolation and spreading, which can effectively soften or enlarge the capsule.The IBM introduces an artificial length scale that may alter the capsule response.

B. Mesh discretization

Membrane mesh resolution and tessellation have little influence on the capsule’s deformation state at small reduced shear rates, although some angular details can be misrepresented.

  • Mesh resolution changes from 320 to 5120 faces while keeping the fluid resolution fixed, enabling direct assessment of membrane discretization.The tested average edge lengths are l̄/∆x = 1.50, 0.75, and 0.38 for coarse, intermediate, and fine meshes.
  • Mesh-node spacing need not generally be smaller than ∆x/2, allowing membrane resolution to be reduced without significantly decreasing simulation accuracy.Because IBM interpolation and spreading scale with node count, this reduction can decrease computing time.
  • Different tessellation methods produce no difference in deformation-parameter evolution, although Gmsh and CGAL meshes misrepresent inclination at small deformations.Their angular time evolution also shows strong deviations in the reported tests.
  • Mesh tessellation and resolution do not significantly change the deformation state, and even a small mesh is accurate at small G.The icosahedron-based mesh was therefore used for the remaining simulations.

C. Interpolation and spreading

The interpolation stencil sets a finite numerical membrane width that can alter capsule deformation, while mesh tessellation and mesh ratio have comparatively limited effects over validated ranges. Increasing hydrodynamic resolution improves agreement with analytic predictions, but stencil choice requires balancing accuracy against fluctuations and permeability.

  • Interpolation-stencil effects: The interpolation stencils introduce a finite membrane width, so computed solutions converge to analytic predictions as LI/LM approaches zero.LI represents numerical membrane thickness, while LM represents the membrane radius.
  • Mesh resolution: Mesh resolution has little effect on D and θ, whereas angular velocity is extremely sensitive to icosahedron-based mesh resolution.The 5120-element mesh can be replaced by 1280- or 320-element meshes with virtually no change in deformation or inclination.
  • Interpolation-stencil effects: The 2-point stencil fails for ¯l/∆x > 2 because fluid penetrates the membrane, whereas φ3 and φ4 show no similar behavior at ¯l/∆x = 2.The larger interpolation range of φ3 and φ4 probably helps prevent fluid passage.
  • Convergence: Increasing Nr brings D and θ toward analytic values, with coupled-system convergence between orders 1 and 2.The exception is inclination-angle convergence with φ2, which is attributed to mesh degradation.
  • Interpolation-stencil effects: Narrower stencils usually reduce average deviations but produce stronger fluctuations, so stencil selection balances accuracy against smoothness.The deviations are identified as an IBM artifact associated with interpolation-stencil width; φ2 at large resolution is affected by accelerated mesh degradation.
  • Effective deformability and rescaling: The finite-width effect can be partially compensated by redefining effective radius or stiffness, although this rescaling is not straightforward for interacting or nonspherical capsules.The authors keep fluid properties fixed and report considerably smaller effective inclination-angle errors after redefining G.
  • Overall assessment: The IBLBFEM accurately captures capsule physics in the infinite-resolution limit, while interpolation stencils remain the dominant source of numerical deviation at finite resolution.The influence of mesh ratio is small across the broad range 0.5 < ¯l/∆x < 1.5.
  • Mesh-ratio effects: 0.5 < ¯l/∆x < 1.5 permits reduced computational requirements without compromising capsule impermeability or accuracy in the studied cases.The authors report reasonable accuracy even at ¯l/∆x = 1.5, although smaller mesh spacing may be needed for stiff or strongly deformed particles.

D. Mesh degradation and volume drift

Mesh degradation emerges when fine membrane meshes are paired with coarse hydrodynamic resolution, while volume drift depends strongly on the interpolation stencil and resolution. In the reported simulations, volume drift is generally small but may matter in longer, lower-resolution runs.

  • Mesh degradation: Mesh degradation is strongest for φ2 at high membrane resolution, especially when hydrodynamic resolution is coarse.For φ2, degradation is most severe with small capsule radius and many faces; φ3 and φ4 preserve the mesh better.
  • Volume drift: 2 · 10^-4 to 8 · 10^-4, 5 · 10^-5 to 2 · 10^-4, and 9 · 10^-6 to 3 · 10^-5 are the observed volume deviations for φ2, φ3, and φ4, respectively.
  • Volume drift: The 4-point interpolation function controls capsule volume more effectively, and higher hydrodynamic resolution reduces volume drift.
  • Volume drift: Although current volume drift is negligible, it could become significant in long-time simulations at smaller resolution.Restoring forces from volume, surface, and bending energies may help avoid this deficiency.

V. CONCLUSIONS

The conclusions identify hydrodynamic discretization as more important than membrane discretization for accuracy and efficiency. They recommend refining the fluid grid for accuracy, reducing mesh resolution for speed, and accounting for an effective numerical membrane thickness.

  • Conclusions: τ values above 1 cause strong deviations, whereas τ ≤ 1 produces results that barely depend on its magnitude.The relaxation parameter also affects simulation time because it changes the LBM time step.
  • Conclusions: The hydrodynamic and mesh resolutions can be changed independently within a certain range without significantly changing the simulated physics.The investigated assumption that average node spacing must be below half a fluid lattice constant was not supported.
  • Conclusions: The IBM interpolation stencil φ introduces the major discretization error, while membrane discretization has a smaller effect.
  • Conclusions: To increase accuracy, first add fluid lattice nodes at fixed mesh resolution; to reduce computing time, retain hydrodynamic resolution and reduce mesh resolution.
  • Conclusions: The interpolation stencil changes effective capsule radius or stiffness through the numerical membrane thickness.This thickness decreases as the hydrodynamic resolution is refined.
  • Conclusions: The study considers only small deformations of a single capsule, so its minor mesh-resolution effect is limited to that setting.
Loading 1004.2416v1…